[R-sig-Geo] Using ME() of spdep: Is there a way to make the process step by step?

[Roger Bivand]

On Mon, 15 Jun 2009, vega_l at yahoo.com wrote:

>
> Hello everyone, sorry for posting this besides the ecology related r 
> list, is the first and last time I will do it. (But in the last digest 
> of the ecology list, only my lonely message was posted, kind of 
> discouraging)

(Things take time, maybe the participants were waiting for someone else 
to reply)

>
> We are trying to remove the spatial structure of our data by employing 
> Moran Eigen Vectors with the spdep package. Our dataset is kind of 
> large, 3694 samples and 13 variables.
>
> The bad news: my Dual Core 2.0 MHz laptop with Ubuntu 9.04 (2.0 Gb RAM) 
> halted after almost 4 days of continuous processing, and all the effort 
> was lost. (Why the calculations took so much time??)
>

If you read the references, you'll see that ME does a brute-force search 
among the eigenvectors, choosing them one at a time depending on the 
degree to which they reduce autocorrelation in the residuals. Solving the 
eigenproblem for a dense 3700 by 3700 matrix takes some time, but is done 
once only. However, each loop (of a maximum of 3700 loops) is doing up to 
3700 regressions and autocorrelation tests (numbers reduce as eigenvectors 
are chosen). The method has never been claimed to scale well to larger 
numbers of observations. Almost certainly, it is trying to include many 
hundred eigenvectors. In addition, the autocorrelation test is being done 
by permutation bootstrap, adding lots of extra pain. This is the was ME 
was written in the underlying papers, so this is the way it is done. The 
alternatives are any of a large variety of GLMM fits, though you probably 
have too much data for WinBUGS to help. Other GLMM will use a spatial 
correlation structured effect. Maybe look at gamm(), or gam() with a 
spline in observation coordinates from mgcv? Why use ME() with a data set 
that is too large for the method?

> I wonder if there is a way to accelerate the calculations with ME() and 
> to save the resulting eigenvalues as they get calculated. Or if it is 
> possible to make the calculation process of ME() step by step.

You can try to fake it by running ME() under debug(), saving the eV matrix 
once it has been created to the global environment (not for the 
faint-hearted). Then just run glm() in the normal way, adding combinations 
of columns of eV, and testing the results with lm.morantest(), which will 
accept a glm object, but treat it as an lm object. Whether your choice of 
colmns will be a good choice is unknown.

If you prefer to stay with ME, make sure that verbose=TRUE (default) so 
that you can see that it is running, change nsim to say 49 (questionable), 
and tighten alpha to say 0.001. Tighten alpha first a lot (1e-8) so that 
you get an error message: base correlation larger than alpha (if you don't 
get the message, stop the function), then slacken it a little, so that 
ME() stops after selecting the first few eigenvectors that reduce 
autocorrelation most. This has the disadvantage that there will be a lot 
of residual autocorrelation.

Have you tried just using the Gabriel graph for neighbours - why do you 
believe that the inverse distances should be used (apart from their use in 
the literature)? Would ME() work faster with a less informative neighbour 
definition?

Hope this helps,

Roger

>
> The code in R is a follows:
>
> taimin.xy = taimin[c("dy", "dx")];
> coordinates(taimin.xy) <- ~dy+dx; #dy, dx, name of the respective
> coordinate columns
> coords<-coordinates(taimin.xy);
> library(spdep);
> TaiminGabrielGraph<-gabrielneigh(coords, nnmult = 12);
> tai.gab.nb<-graph2nb(TaiminGabrielGraph,sym=TRUE);
> nbtaim_distsg <- nbdists(tai.gab.nb, coords);
> nbtaim_simsg <- lapply(nbtaim_distsg, function(x) (1-((x/(4*50))^2)) );
> MEtaig.listw <- nb2listw(tai.gab.nb, glist=nbtaim_simsg, style="B");
> sevmtaig <- ME(Presabs
> ~Age+Curv+Zon2005+ZoneMeiji+Wi+Sol+Slope+Seadist+Elev+Basin,data=taimin,
> family=binomial,listw=MEtaig.listw)
>
> I hope someone can give me some help.
>
> Thank you in advance!
>
>
>      ____________________________________________________________________________________
> ¡Obtén la mejor experiencia en la web!
> Descarga gratis el nuevo Internet Explorer 8.
> http://downloads.yahoo.com/ieak8/?l=e1
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

-- 
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no