[R-sig-Geo] Different Moran's I for same model residuals in different runs using moran.mc{spdep}

Xingli Giam xingli.giam at adelaide.edu.au
Sat Jan 31 06:05:36 CET 2009


Hello List,

I am using moran.mc{spdep} to extract values of Moran's I for a simple GLM
model residuals of binned distances.

However, everytime I run moran.mc() in a new session, Moran's I values are
different. So when I manually plot a correlogram (model residuals Moran's I vs
distance classes), the correlogram looks different every session, although I am
using the same data.

Is this normal? Or is there anything wrong with my code or neighbour weight
lists? When I try removing the glist, and just use (ME0_200.listw <- nb2listw
(nb0_200, style="W", zero.policy=T), I continue to get different Moran's I
values every time I run it.

Here's my code:
## neighbourhood bins for correlogram (0-5000km, in 200 km increments)
nb0_200 <- dnearneigh(as.matrix(dat$x_long, dat$y_lat), 0, 200, longlat=T)
nb0_200_dists <- nbdists(nb0_200, as.matrix(dat$x_long, dat$y_lat))
nb0_200_sims <- lapply(nb0_200_dists, function(x) (1/x))
ME0_200.listw <- nb2listw(nb0_200, glist=nb0_200_sims, style="W", zero.policy=T)

nb200_400 <- dnearneigh(as.matrix(dat$x_long, dat$y_lat), 200, 400, longlat=T)
nb200_400_dists <- nbdists(nb200_400, as.matrix(dat$x_long, dat$y_lat))
nb200_400_sims <- lapply(nb200_400_dists, function(x) (1/x))
ME200_400.listw <- nb2listw(nb200_400, glist=nb200_400_sims, style="W",
zero.policy=T)

...

nb4800_5000 <- dnearneigh(as.matrix(dat$x_long, dat$y_lat), 4800, 5000,
longlat=T)
nb4800_5000_dists <- nbdists(nb4800_5000, as.matrix(dat$x_long, dat$y_lat))
nb4800_5000_sims <- lapply(nb4800_5000_dists, function(x) (1/x))
ME4800_5000.listw <- nb2listw(nb4800_5000, glist=nb4800_5000_sims, style="W",
zero.policy=T)

# Check for spatial autocorrelation using permutation tests for Moran's I
(Monte Carlo sim of 10000 runs)
n1<- moran.mc(residuals(mod1), nsim=9999, listw=ME0_200.listw,
na.action=na.omit, zero.policy=T)
n2<- moran.mc(residuals(mod1), nsim=9999, listw=ME200_400.listw,
na.action=na.omit, zero.policy=T)
...

n25<- moran.mc(residuals(mod1), nsim=9999, listw=ME4800_5000.listw,
na.action=na.omit, zero.policy=T)

# plot correlogram of glm residuals
n.vec <- c(n1$statistic, n2$statistic, n3$statistic, n4$statistic, n5
$statistic, n6$statistic,
n7$statistic, n8$statistic, n9$statistic, n10$statistic, n11$statistic, n12
$statistic,
n13$statistic, n14$statistic, n15$statistic, n16$statistic, n17$statistic, n18
$statistic,
n19$statistic, n20$statistic, n21$statistic, n22$statistic, n23$statistic, n24
$statistic,
n25$statistic)

plot(n.vec, type="b", col="red", ylim=c(-0.5,0.5))

Many thanks,
Xingli




More information about the R-sig-Geo mailing list