Re: [R-sig-Geo] overlay spatialgrid on spatial polygons data frame

[Roger Bivand]

On Thu, 20 May 2010, Agustin Lobo wrote:

I have an SpPolDF and a coarse SpatialGrid, i.e.,
srdf=SpatialPolygonsDataFrame(sr, data.frame(cbind(1:3,5:3),
row.names=c("r1","r2","r3")))
#the one used in help(overlay);
s...@data
  X1 X2
r1  1  5
r2  2  4
r3  3  3
bbox(srdf)
    min    max
x 178544 181477
y 329665 333676

test <- SpatialGrid(GridTopology(cellcentre.offset=c(178440,329600), 
cellsize=c(500,500), cells.dim=c(8,10)))
plot(srdf)
plot(test,add=T)

My goal is to "overlay" test over srdf and tranfer the data in
s...@data to a data frame associated
to test (=> test must become a SpatialGridDF) so that each cell has
the average value of the overlaid s...@data polygon(s) weighted
by area: in case 25% of a given cell in test would overlay polygon r1
and 75% polygon r2, the value of X1 for that cell would be
0.25*1 + 0.75*2

This involves a topology operation on two sets of polygons, and is not 
supported. You can approximate it by creating a finer raster and using 
that subsequently aggregate, for example using the polygons in ?overlay 
to form srdf, then:

test2 <- SpatialGrid(GridTopology(cellcentre.offset=c(178440,329600),
 cellsize=c(50,50), cells.dim=c(80,100)))
o <- overlay(test2, srdf)
o1 <- SpatialPointsDataFrame(coordinates(test2), data=data.frame(o=o))
o2 <- o1[!is.na(o1$o),]
o2$X1 <- srdf$X1[o2$o]
o2$X2 <- srdf$X2[o2$o]
u <- overlay(test, o2)
uu <- aggregate(slot(o2, "data"), list(u), mean)
X1 <- rep(as.numeric(NA), nrow(coordinates(test)))
X2 <- rep(as.numeric(NA), nrow(coordinates(test)))
X1[uu$Group.1] <- uu$X1
X2[uu$Group.1] <- uu$X2
testdf <- 
SpatialGridDataFrame(GridTopology(cellcentre.offset=c(178440,329600),
cellsize=c(500,500), cells.dim=c(8,10)), data=data.frame(X1, X2))
spl <- list("sp.lines", as(srdf, "SpatialLines"), col="grey30")
spplot(testdf, c("X1", "X2"), col.regions=bpy.colors(20), sp.layout=spl)

The finer resolution should be fine enough to capture the relative areas 
of the polygons adequately.

Roger



According to help(overlay),
"Value
a numerical array of indices of x on locations of y, or a data.frame
with (possibly aggregate) properties of x in units of y"

I've tried this:
overlay(srdf, test, fn=mean, na.rm)

    X1 X2
NA   NA NA
NA.1 NA NA

Cannot go beyond, perhaps overlay() is not suited for this: the
definition of overlay in its help page
"overlay combines points (or grids) and polygons by performing
point-in-polygon operation on all point-polygons combinations"
is a bit confusing, the first part of the sentence gives me some hope
but the second one seems to exclude what I'm trying to do.

Any help appreciated (or perhaps pointing to more doc?)


Agus

_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@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.biv...@nhh.no

_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-sig-geo