On Tue, 6 Jul 2010, J. Sebastian Tello wrote:
Dear all,
I am currently trying to find out if there is a way to calculate a pseudo R2
value from the output of SAR error models or spatial GLS models. I have seen
this presented in publications, and the program SAM (Spatial analysis in
macroecology) produces this value... but I need to implement the calculation of
this pseudo R2 in a piece of R code I am writing. The authors of SAM say they
calculate their pseudo R2 using Nagelkerke's 1991 furmulation, which I also
have seen used in various other places; so I have tried to use it with R output
from functions errorsarlm and gls in R, but the results do not make sense.
Here is what I am doing:
For SAR error model:
my.lm<-lm(Y~X1+X2+X3, data=DATA)
coords<-as.matrix(cbind(DATA$x,DATA$y))
nb.5<-dnearneigh(coords,0,5)
nb.5.w<-nb2listw(nb.5, glist=NULL, style="W", zero.policy=TRUE)
my.sar.error<- errorsarlm(my.lm, listw=nb.5.w, na.action="na.omit")
lL.B <- logLik(my.sar.error)
lL.0 <- updata(my.sar.error, ~1)
Really? There is no update() method for sarlm objects. so do not expect
anything sensible. Do however see the Nagelkerke=TRUE argument to the
summary() methos for sarlm objects:
library(spdep)
data(oldcol)
COL.errW.eig <- errorsarlm(CRIME ~ INC + HOVAL, data=COL.OLD,
nb2listw(COL.nb, style="W"))
summary(COL.errW.eig, Nagelkerke=TRUE)
Hope this helps,
Roger
yN<-nrow(DATA)
Nagelkerke<-(1-exp((lL.B-lL.0)*(-2/yN)))/(1-exp((lL.0)*(2/yN)))
Nagelkerke[1]
For the GLS it looks almost identical:
exp.correl.structure<-corExp(form=~x+y)
my.gls<-gls(Y~X1+X2+X3, data=DATA, correlation=exp.correl.structure,
method="ML"))
lL.B <- logLik(my.gls)
lL.0 <- updata(my.gls, ~1)
yN<-nrow(DATA)
Nagelkerke<-(1-exp((lL.B-lL.0)*(-2/yN)))/(1-exp((lL.0)*(2/yN)))
Nagelkerke[1]
From running these analyses, I obtain that the pseudo R2 is of almost 0.997,
while the ols R2 is of only 0.59. I know this two numbers are not really
comparable (I do not expect them to be the same), but also the pseudo R2
produced by SAM is of around 0.55. So I am not sure what is wrong. I think the
pseudo R2 I am obtaining is for the model after accountiong for all the spatial
autocorrelation, rather than only from what can be predicted by the explanatory
variables in the model. If any one knows of a way to do this, or has any ideas
of approaches to take or resources to visit, I would appreciate the input.
Thanks in advance,
Sebastian
J. Sebasti?n Tello
Department of Biological Sciences
285 Life Sciences Building
Louisiana State University
Baton Rouge, LA, 70803
(225) 578-4284 (office and lab.)
[[alternative HTML version deleted]]
--
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.biv...@nhh.no
_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
