[R] troubles performing Moran.I test
Roger Bivand
Roger.Bivand at nhh.no
Thu Jan 8 13:58:54 CET 2009
<Barbara.Spillmann <at> agrar.uni-giessen.de> writes:
>
> dear R users,
>
> I have troubles performing Moran.I test as suggested on
> http://www.ats.ucla.edu/stat/r/faq/morans_i.htm
Using the ape package which is written for a different application area is
not necessarily a good idea - the web page you refer to also has basic
blunders, mainly calculating and using planar distances when spherical
were called for. You repeat this in your data. Distances should be measured
by Great Circle between geographical coordinates.
>
> my spatial data are longitude and lattitide of communities. The
> calculation of the inverse distance matrix according to the homepage
> (using my data)
>
...
> the data which might be spatially autocorrelated is LN(Unem05/Unem98,
> the LN of the development in unemployment rates in the communities
> between 1998 and 2005.
>
> I have 426 communities in total and I don's see what might be wrong
> with the data...However, I have some NAs in there...
> when I try to perform the test using:
>
> Moran.I(datAL$LN.Rt05.Rt98., ALdist.inv, na.rm=TRUE)
The function takes a number of short cuts that are not obvious without reading
the code, or possibly reading Gittleman & Kot - it performs a hidden row
standardisation of the weights.
If we start from ape's Moran.I() example:
set.seed(1)
tr <- rtree(30)
x <- rnorm(30)
w <- 1/cophenetic(tr)
diag(w) <- 0
Moran.I(x, w)
we can recreate the same results using:
library(spdep)
# convert w to a row standardised general weights object
lw <- mat2listw(w)
lwW <- nb2listw(lw$neighbours, glist=lw$weights, style="W")
moran.test(x, lwW, alternative="two.sided")
Note that spdep provides functions for calculating Great Circle distances,
see dnearneigh() and nbdist().
>
> I get the following error message:
> Fehler in if (obs <= ei) 2 * pv else 2 * (1 - pv) :
> Fehlender Wert, wo TRUE/FALSE nötig ist
>
Inserting an NA into x:
is.na(x[5]) <- TRUE
Moran.I(x, w, na.rm=TRUE)
and
xc <- complete.cases(x)
wc <- w[xc, xc]
lwc <- mat2listw(wc)
lwWc <- nb2listw(lwc$neighbours, glist=lwc$weights, style="W")
moran.test(x[xc], lwWc, alternative="two.sided")
agree in the coefficient but not beyond that (ei is not for the correct n
in Moran.I()). By the way, Moran's I really doesn't make sense for missing
data. There is a provision in moran.test() for subsetting, but not for the
general weights matrix you are using.
> can anyone give me a hint what is going wrong??
Since it isn't NA in your variable, it must be something else, so use
debug(Moran.I) to see which of obs and/or ei are NAs in the line you quote.
Have you had the opportunity to review the Spatial task view on CRAN?
It is not impossible that it might shed some light on your problem, possibly
more than the authority you cite, including a link to the R-sig-geo list.
Roger Bivand
>
> many thanks in advance!!
>
> Barbara
More information about the R-help
mailing list