Skip to content
Prev 25275 / 29559 Next

Moran Index n ape and spdep package

On Tue, 17 Jan 2017, Stefano Menichetti wrote:

            
library(sp)
data(meuse)
myData.dists <- as.matrix(dist(cbind(meuse$x, meuse$y)))
myData.dists.inv <- 1/myData.dists
diag(myData.dists.inv) <- 0
library("ape")
Moran.I(meuse$zinc, myData.dists.inv)

This is what you had already, and spdep can do the same, but in a much 
more flexible way, where the user has control of what is happening.

library(spdep)
nb_inf <- dnearneigh(cbind(meuse$x, meuse$y), 0, 10000)
nb_inf # dense
dists_inf <- nbdists(nb_inf, cbind(meuse$x, meuse$y))
dists_inf_inv <- lapply(dists_inf, function(x) 1/x)
lw <- nb2listw(nb_inf, glist=dists_inf_inv, style="B")
Mlw <- listw2mat(lw)
all.equal(Mlw, myData.dists.inv, check.attributes=FALSE)
Moran.I(meuse$zinc, myData.dists.inv)
Moran.I(meuse$zinc, Mlw)
moran.test(meuse$zinc, lw)

# not the same, ape::Moran.I converts to style="W" - row-standardising the 
# weights - internally and transparently

spdepI <- moran.test(meuse$zinc, nb2listw(nb_inf, glist=dists_inf_inv,
  style="W"))
apeI1 <- Moran.I(meuse$zinc, Mlw)
apeI0 <- Moran.I(meuse$zinc, myData.dists.inv)
all.equal(unname(spdepI$estimate[1]), unname(apeI0$observed))
all.equal(unname(spdepI$estimate[1]), unname(apeI1$observed))
all.equal(unname(spdepI$statistic),
  unname((apeI0$observed-apeI0$expected)/apeI0$sd))

# note that spdep::moran.test and ape::Moran.I have different defaults for 
# alternative, and that ape::Moran.I implements the randomisation variance 
# without stating that it does this

All now the same. But:

nb_tri <- tri2nb(cbind(meuse$x, meuse$y))
nb_tri # sparse
dists_tri <- nbdists(nb_tri, cbind(meuse$x, meuse$y))
dists_tri_inv <- lapply(dists_tri, function(x) 1/x)
moran.test(meuse$zinc, nb2listw(nb_tri, glist=dists_tri_inv, style="W"))

See vignette("nb") for different ways of finding neighbours. Sparse 
weights are always to be preferred over dense, because the covariance 
matrix of the observations becomes dense by inversion when spatial 
autocorrelation is present.

Moran's I can always be fooled by a wrong mean model (Moran.I and 
moran.test use the mean of the variable of interest). If there is a 
spatial trend in otherwise uncorrelated data, and the trend is not 
included in the mean model, the measure will find spurious 
autocorrelation.

Hope this clarifies,

Roger