On Mon, 15 Aug 2011, Michelle Convey wrote:
Hi there,
I am working on a project with the spgwr package. I am able to map the
coefficients from my GWR analysis and am able to get summary data on the
exponentials on the same data, but am not sure how to map them. I'm probably
missing something quite obvious as I am fairly new to this, but any help would
be greatly appreciated. Please see a snippet of the code below:
library spgwr
library maptools
Not needed.
library GISTools
There is no such package on CRAN. All examples must be fully reproducible.
## GWR
# Get the bandwidth:
a.bw <- gwr.sel(aa~bb, data = data.pt, adapt = TRUE)
a.bw
# Create model:
gwr.a = gwr(aa~bb, data=data.pt, adapt=a.bw, fit.points=grid1)
# Co-efficients summary:
m <- length(gwr.a$lm$coefficients)
CM <- t(apply(as(gwr.a$SDF, "data.frame")[, (1 + (1:m)), drop = FALSE],
2, summary))[, c(1:3, 5, 6)]
# This puts it in a matrix:
if (is.null(dim(CM)))
CM <- t(as.matrix(CM))
Please do mot mess with things you don't need to; gwr.a$SDF is a
Spatial*DataFrame, so simply do summary() on the relevant columns.
Find the columns by inspecting names(gwr.a$SDF) and say:
summary(gwr.a$SDF[,c("var1", "var2", ...)])
For inspection, summary(gwr.a$SDF) is sufficient. You have copied code
from within print.gwr() without understanding the underlying objects.
To get the exponents, use the fact that gwr.a$SDF is a
Spatial*DataFrame, so simply say:
gwr.a$SDF$evar1 <- exp(gwr.a$SDF$var1)
...
and repeat summary() with the new variables, and plot() or whatever - I
suggest spplot() to map coefficients.
Also be aware that GWR is a rather doubtful procedure, as its authors have
not replied to published documented criticism, and that this package is
intended to assist research in finding out why it is unreliable. Any
coefficient variability may indicate model misspecification (but may occur
when not misspecification is present), which you should correct, rather
than use GWR.
Roger
# Exponentials:
CM <- cbind(exp(CM), exp(coefficients(gwr.a$lm)))
# This prints out a summary matrix:
colnames(CM) <- c(colnames(CM)[1:5], "Global")
printCoefmat(CM)
## MAKE MAP - This only maps the coeffiecients - I would like to be able to map
the exponentials
# Plot the results:
gwr.a = SpatialPixelsDataFrame(gwr.a$SDF, data.frame(gwr.a$SDF))
# Changes colour and map quantiles:
gwr.shades = auto.shading(gwr.res_a$bb, n=4,
cols=brewer.pal(4,'Blues'),cutter=quantile.cuts)
level.plot(gwr.res_a, index = "LTI", shades = gwr.shades)
Any help would be gratefully received.
Many thanks,
M
_______________________________________________
R-sig-Geo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
--
Roger Bivand
Department of Economics, NHH Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: [email protected]
_______________________________________________
R-sig-Geo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
