[ESRI-L] outline polygons of point clumps

On Fri, 12 May 2006, Xiaohua Dai wrote:

My idea is to use your procedure to obtain food-tree patches of
rhinos, then relate the patches to rhino movement pattern. For
example, rhinos may move more slowly in food patches than outside the
patches. I forgot one thing, could you also tell me how to write these
hulls into a polygon shape file? write.polylistShape {maptools} or
convert.to.shapefile {shapefiles}?

set.seed(1)
xy <- matrix(runif(500, 0, 10), ncol=2)
xy_clusts <- hclust(dist(xy), method="complete")
cl_10 <- cutree(xy_clusts, 10)
which_cl_10 <- tapply(1:nrow(xy), cl_10, function(i) xy[i,])
chulls_cl_10 <- lapply(which_cl_10, function(x) x[chull(x),])
# so far as before - now:

n <- length(chulls_cl_10)
library(sp)
polygons <- lapply(1:n, function(i) {
 chulls_cl_10[[i]] <- rbind(chulls_cl_10[[i]], chulls_cl_10[[i]][1,])
# the convex hulls do not join first and last points, so we copy here
 Polygons(list(Polygon(coords=chulls_cl_10[[i]])), ID=i) })
out <- SpatialPolygonsDataFrame(SpatialPolygons(polygons),
 data=data.frame(ID=1:n))
plot(out)
# note standard-violating intersecting polygons!
tempfile <- tempfile()
library(maptools)
writePolyShape(out, tempfile)
in_again <- readShapePoly(tempfile)
plot(in_again, border="blue", add=TRUE)

The trick is to construct a SpatialPolygons object defined in the sp
package, and write it and a single data field out to a shapefile. The
shapefile is strictly wrong, because the polygons can intersect, but at
least it's something to start with.

Roger

Thanks
Xiaohua

On 5/12/06, Roger Bivand <Roger.Bivand at nhh.no> wrote:

On Fri, 12 May 2006, Xiaohua Dai wrote:

Hi Roger,

Thank you for your kindly help. It is exactly what I want. I will try
to understand your programming.  I think the program should be useful
in many fields as you indicated.

I did not expect the problem can be solved so easily. Again, R is very
powerful with so many good functions. I am quite interested in lapply
and tapply. I will find some references/websites on their
applications.

Note that using a different method= argument to hclust, or a different
metric to dist(), you will get different clusters, as well as getting
different numbers for different cutree() values. If you have specific
range distances (from an off-list message, the points are rhino), you
could also use other methods to create your own clusters where the first
nearest neighbour distance is less than your empirically observed range -
there should be a literature on this among ecologists, which others on the
list may be able to help with. Is there anything in the adehabitat
package, or others?

Roger

Also thanks to remind me of the post guide.

Best wishes
Xiaohua

On 5/12/06, Roger Bivand <Roger.Bivand at nhh.no> wrote:

On Fri, 12 May 2006, Xiaohua Dai wrote:

Sorry, I did not make my question clear. Since I have a point theme
with many points, some of them may clump together. the problems here
are:

A more specific follow-up after my reply directly to R-help.

I'm assuming you have a running R, and that convex hulls are what you
want.

Take a random data set:

set.seed(1)
xy <- matrix(runif(500, 0, 10), ncol=2)
xy_clusts <- hclust(dist(xy), method="complete")
# complete linkage hierarchical clustering
plot(xy_clusts)
# shows the clustering tree
cl_10 <- cutree(xy_clusts, 10)
cl_20 <- cutree(xy_clusts, 20)
cl_30 <- cutree(xy_clusts, 30)
# cut the tree - the objects contain the memberships
which_cl_10 <- tapply(1:nrow(xy), cl_10, function(i) xy[i,])
chulls_cl_10 <- lapply(which_cl_10, function(x) x[chull(x),])
# construct convex hull polygons for each cluster
plot(xy)
res <- lapply(chulls_cl_10, polygon)
# and repeat for cl_20 and cl_30
which_cl_20 <- tapply(1:nrow(xy), cl_20, function(i) xy[i,])
chulls_cl_20 <- lapply(which_cl_20, function(x) x[chull(x),])
plot(xy)
res <- lapply(chulls_cl_20, polygon)
which_cl_30 <- tapply(1:nrow(xy), cl_30, function(i) xy[i,])
chulls_cl_30 <- lapply(which_cl_30, function(x) x[chull(x),])
plot(xy)
res <- lapply(chulls_cl_30, polygon)

If you need the list of convex hulls out as a shapefile, we can do that,
if you need a raster, I'd suggest using kernel density with different
bandwidths for a start, and NA out the densities below a chosen threshold.

Hope this helps,

Roger Bivand

PS. Use package maptools, function readShapePoints() to read your
shapefile, and coordinates() of the input object to extract the
coordinates. If need be, project from geographical coordinates using
transform methods in package rgdal.

1.  how to find clumps in a point theme?
2.  the convex-hull extension I found only deal with all the points in
a theme at each time?  how to make each convex hull around each point
clump automatically?

Thanks.

Xiaohua



On 5/12/06, Bob Booth <bbooth at esri.com> wrote:

Xiaohua,

That would be one way to do it. There are others.
Try searching ArcScripts for "convex hull"
http://arcscripts.esri.com/

For example:
http://arcscripts.esri.com/details.asp?dbid=14535
or
http://arcscripts.esri.com/details.asp?dbid=12084

Bob


-----Original Message-----
From: ESRI-L [mailto:ESRI-L at esri.com] On Behalf Of Xiaohua Dai
Sent: Friday, May 12, 2006 2:33 AM
To: ESRI-L at esri.com
Subject: [ESRI-L] outline polygons of point clumps

Dear all,

How to generate one outline polygon for each point clump? Are there
any present functions in ArcView, ArcGIS, R or some freewares? I just
had a quick look at the package adehabitat and did not find the
function.

To my knowledge, I could do it as follows: 1) make a grid map of my
study area with cell values = 0; 2) assign 1 to the cells containing
at least one point; 3) convert the 1-value cells into polygons.

Am I right?

Thanks
Xiaohua