Homework 7: Spatial Autocorrelation

Starting with lagged scatterplots

mexicobirds<-readRDS("birdRichnessMexico.rds")
head(mexicobirds)
## Loading required package: sp
##               coordinates nSpecies      lat       long
## 323 (409123.8, -486270.1)      286 19.60912 -101.06139
## 242  (-216075.5, 88686.7)      250 24.77660 -107.16283
## 333 (842161.7, -564968.9)      385 18.74534  -96.94086
## 115 (-890505.6, 405332.8)      153 27.36783 -114.10095
## 453  (41190.45, 186561.5)      190 25.66340 -104.58487
## 227 (-425309.4, 466154.9)      236 28.10320 -109.36972
library(ggmap)
## Loading required package: ggplot2
mexMap <- get_googlemap(center = c(lon = -102, lat = 23), zoom = 5, maptype = 'terrain')
## Map from URL : http://maps.googleapis.com/maps/api/staticmap?center=23,-102&zoom=5&size=640x640&scale=2&maptype=terrain&sensor=false
map1 <- ggmap(mexMap) 
map1 <- map1 + geom_point(data=data.frame(mexicobirds), aes(x=long, y=lat, 
                                                      color=nSpecies),size=6)
map1 <- map1 + labs(x="Longtitude",y="Latitude",title="Bird Species Richness")
map1
## Warning: Removed 1 rows containing missing values (geom_point).

library(gstat)
hscat(nSpecies~1,data=mexicobirds,breaks = seq(0,2e6,by=2e5))

The lagged scatterplots show that there is strong spatial autocorrelation at distances shorter than 8E5 meters at which point the trend switches to very weakly negative autocorrelation. This correlated to 80km.

Variograms

birdVarCloud <- variogram(nSpecies~1, mexicobirds, cloud = TRUE)
plot(birdVarCloud,pch=20,cex=1.5,col="black",
     ylab=expression("Semivariance ("*gamma*")"),
     xlab="Distance (m)",main = "Species Richness Mexico Birds")

birdvar<-variogram(nSpecies~1, mexicobirds, cloud = FALSE)
plot(birdvar,pch=20,cex=1.5,col="black",
     ylab=expression("Semivariance ("*gamma*")"),
     xlab="Distance (m)",main = "Species Richness Mexico Birds")

birdvar2<-variogram(nSpecies~1, mexicobirds, width=1000,cloud = FALSE)
plot(birdvar2,pch=20,cex=1.5,col="black",
     ylab=expression("Semivariance ("*gamma*")"),
     xlab="Distance (m)",main = "Species Richness Mexico Birds")

The variograms show that the sill appears at roughly 800,000 meteres were the semi variance starts to peter off but then there is a point of high semivariance after than could be an outlier or possibly the result of edge effect. The width argument control how many points are displayed, so when I decreased width there were many more points and the shape of the variogram was harder to interpret, with no obvious sill.

Calculate Moran’s I

library(spdep)
## Loading required package: Matrix
library(ncf)
coords<-coordinates(mexicobirds)
nb<-tri2nb(coords,row.names=NULL)
## 
##      PLEASE NOTE:  The components "delsgs" and "summary" of the
##  object returned by deldir() are now DATA FRAMES rather than
##  matrices (as they were prior to release 0.0-18).
##  See help("deldir").
##  
##      PLEASE NOTE: The process that deldir() uses for determining
##  duplicated points has changed from that used in version
##  0.0-9 of this package (and previously). See help("deldir").
moran.test(mexicobirds$nSpecies,nb2listw(nb))
## 
##  Moran I test under randomisation
## 
## data:  mexicobirds$nSpecies  
## weights: nb2listw(nb)  
## 
## Moran I statistic standard deviate = 17.166, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic       Expectation          Variance 
##       0.703649266      -0.005025126       0.001704339
species<- spline.correlog(x=mexicobirds$long, y=mexicobirds$lat,
                         z=mexicobirds$nSpecies, resamp=100, xmax=1500, quiet=TRUE,latlon=TRUE)
plot(species,text=TRUE)

x <- species$real$predicted$x
y <- species$real$predicted$y
lCI <- species$boot$boot.summary$predicted$y["0.025",]
uCI <- species$boot$boot.summary$predicted$y["0.975",]
xx <- c(x,rev(x))
yy <- c(lCI,rev(uCI))
plot(x,y, xlim = c(0,1500), ylim = c(-1, 1), type = "n", 
     xlab = "Distance (m)", ylab = "Moran's I")
polygon(xx,yy,col="grey",border=NA)
abline(h=0,lty="dashed")
lines(x,y)

The Moran’s I spline shows that there is strong positive autocorrelation at smaller distances and small negative autocorrelation after a distance of 1000 meters

Directional Variograms

library(sp)
directional<-variogram(nSpecies~1,mexicobirds)
plot(directional)

dirt<-variogram(nSpecies~1, mexicobirds, alpha=c(0,45,90,135))
plot(dirt, main="Species Richness of Mexico Birds via Direction")

In the directional variogram, there appears to be anisotropy because depending on which direction the pairwise pairs are calculated in, the relationship between gamma an distance changes. So in the northern direction there doesn’t appear to be a sill at all, the species richness just continually increases. However, in the eastern direction the sill appears at around 80km which is similar to the overall variogram except there are two steps rather than approaching a single asymptote. In the southern direction the species richness also appears to continue increasing without approaching a sill. In the western direction, the sill appears at much shorter distances (20km) than the overall variogram. This is probably due to the number of points that exist in each direction, because Mexico is longer than it is wide.