[R-sig-Geo] merge sliver polygons
Roger Bivand
Roger.Bivand at nhh.no
Fri Sep 5 21:25:25 CEST 2008
On Fri, 5 Sep 2008, Murray Richardson wrote:
> Thanks Roger,
>
> I have singletons with areas. I think the slivers will be quite obvious
> based on their areas relative to the rest of the polygon areas. I could even
> compute some perimeter to area index I supposed since they tend to be long a
> skinny.
>
> So, my main uncertainty is still just how to work with the neighbour list to
> merge the sliver polygons with the LARGEST non-sliver polygon they share an
> edge with. i.e. to generate the new ids to use in unionSpatialPolygons. I
> can't seem to find a good example of how to use the nb object in such a way.
OK, so the Polygons are all singletons. This isn't tried:
is_a_sliver <- AREAS < very_small
library(spdep)
nb <- poly2nb(your_SPDF)
# check length(is_a_sliver) against length(nb)
crd <- card(nb)
IDs <- sapply(slot(your_SPDF, "polygons"), function(x) slot(x, "ID"))
for (i in seq(along=nb)) {
if (is_a_sliver[i] && crd[i] > 0) {
nbi <- nb[[i]]
max_area <- nbi[which.max(AREAS[nbi])]
IDs[i] <- IDs[max_area]
}
}
then use IDs in unionSpatialPolygons()
I repeat, untried, but something like this should work.
Roger
PS. If you have the perimeter, you could try it as a sanity check if area
alone give spurious results.
>
> cheers
>
> Murray
>
>
> Roger Bivand wrote:
>> On Fri, 5 Sep 2008, Murray Richardson wrote:
>>
>> > Hi Roger
>> >
>> > I'm just getting around to trying this out. I must say, I'm not clear on
>> > how to work with the neighbour list object to accomplish this.
>> > Specifically, how would I identify the largest neighbours of the sliver
>> > candidates?
>> >
>> > Furthermore, do I first subset the polygons according to the area
>> > threshold, and then construct the neighbour list for only these
>> > candidate slivers? I'm not clear on how to identify neighbours for a
>> > subset of polygons.
>>
>> My idea (IIRC!) was to look at a distribution of the areas of the Polygon
>> objects in the dataset. A good deal will depend on whether the slivers are
>> represented as separate Polygons objects (note the _s_), or whether they
>> are Polygon objects belonging to Polygons objects. So nimbleness will be
>> needed!
>>
>> library(maptools)
>> xx <- readShapeSpatial(system.file("shapes/sids.shp",
>> package="maptools")[1],
>> IDvar="FIPSNO", proj4string=CRS("+proj=longlat +ellps=clrk66"))
>> pls <- slot(xx, "polygons")
>> AREAS <- sapply(pls, function(x) sapply(slot(x, "Polygons"), function(y)
>> slot(y, "area")))
>> summary(sapply(AREAS, length))
>> summary(unlist(AREAS))
>>
>> is not a very good example, because the "areas" are planar in non-planar
>> coordinates, but I hope you see the point. If you have a clear break in
>> the distribution, you can find out which they are, and where they are. If
>> all the Polygons objects (the ones in the "polygons" slot of your
>> SpatialPolygons object, pls) are singletons, things get much easier,
>> because the merging will be between Polygons objects, provided directly in
>> maptools.
>>
>> This should be a start - once you have singleton Polygons objects with
>> areas as attributes, we can continue.
>>
>> Hope this helps,
>>
>> Roger
>>
>> >
>> > Sorry to bother you with this again. It will be very useful if it works
>> > as I would also like to use it for region merging based on other
>> > attributes.
>> >
>> > Still waiting on that book...!
>>
>> PS. Edzer, Virgilio and I did see the only copy in Europe in Dortmund at
>> useR! three weeks ago. It turns out that it was printed in the US, and is
>> at or awaited at the publisher's warehouse for dispatch. Things seem to
>> take time in the "real" world!
>>
>> >
>> > Murray
>> >
>> >
>> >
>> > Roger Bivand wrote:
>> > > On Wed, 13 Aug 2008, Murray Richardson wrote:
>> > >
>> > > > Hello again r.sig.geo list,
>> > > > > Thanks Roger, for help on my previous question regarding
>> > > iterating > through a shapefile.
>> > > > > I'm sure once I receive my copy of "Applied Spatial Data
>> > > Analysis with > R" I will find answers to simple questions like this
>> > > on my own, but in > the meantime....
>> > > > > Is it possible to merge sliver polygons that fall below a
>> > > certain > threshold area with adjacent neighbours (e.g. perhaps using
>> > > > unionSpatialPolygons but without aggregating any polygons?). If a
>> > > > sliver shares edges with more than one polygon, it doesn't really
>> > > matter > which one it merges with, but if I had to choose a rule I
>> > > would have it > merge with the largest one.
>> > >
>> > > Not such a simple question ...
>> > >
>> > > Both the Polygon and Polygons objects in the SpatialPolygons object
>> > > have
>> > > "area" slots, with different roles. The Polygon objects have a
>> > > correct
>> > > naive area in the geometry of the coordinates taken as planar. The
>> > > Polygons objects use the "gross" area of Polygon objects belonging to
>> > > them, but "only" to provide the plot order (plot from largest to
>> > > smallest
>> > > to avoid over-painting).
>> > >
>> > > If you "trust" the area slot of the Polygons objects (beware of hole
>> > > Polygon objects), you can first find your candidate slivers by
>> > > retrieving
>> > > the areas by:
>> > >
>> > > Polygons_areas <- sapply(slot(SPobj, "polygons"),
>> > > function(x) slot(x, "area"))
>> > >
>> > > and set a cutoff. Then use poly2nb(SPobj, queen=FALSE) in spdep to
>> > > find
>> > > the neighbours (rook criterion). Next use the output object to
>> > > identify
>> > > the largest neighbours of the sliver candidates, and build a "new
>> > > Polygons" ID vector. Finally, use unionSpatialPolygons(). I'm
>> > > assuming you
>> > > wouldn't have asked if there was useful data in the slivers!
>> > >
>> > > Hope this helps,
>> > >
>> > > Roger
>> > >
>> > > > > Thanks in advance,
>> > > > > Murray Richardson
>> > > >
>>>>> _______________________________________________
>> > > > 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
More information about the R-sig-Geo
mailing list