Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Combining polygons and calculating their area (i.e. number of cells) in R

Tags:

r

r-raster

I have a simple raster (created with R-package: raster). Using the function "rasterToPolygons" I get polygons of all raster cells that contain the value "1":

library(raster) 
dat = list()
dat$x = seq(1.5, by = 10, len = 10)
dat$y = seq(3.5, by = 10, len = 15)
dat$z = matrix(sample(c(0,1), size = 10*15, replace = T), 10, 15)

r=raster(dat);plot(r)

r_poly = rasterToPolygons(r, fun = function(r) {r == 1}, dissolve = F)
plot(r_poly, add = T)

I do not use "dissolve = T" to avoid that all polygons are merged into one big polygon. Instead, I wish to obtain a new SpatialPolygonsDataFrame in which all polygons that share an edge or a point are combined. Polygons that are clearly separated should be identifiable as individual polygones. Based on the new SpatialPolygonsDataFrame I would like to analyze the size of the combined polygones as follows:

b = extract(r,r_poly_new) # "r_poly_new" contains the combined polygons
str(b)                    # list of clearly separated polygons
tab = lapply(b,table)      
tab

My question is twofold: 1) How to combine polygones that share an edge or point? 2) How to get this information into a format which allows analyzing the areas of the combined polygones? Thanks very much for your feedback.

like image 533
Alex Z Avatar asked Dec 18 '13 13:12

Alex Z


People also ask

How do you count points in a polygon in R?

If you transform the points and France polygons to sf objects, you can use st_intersects() to count the number of points in each polygon.

How do you calculate the area of a raster in R?

The „area“ function of R's „raster package“ results in a vector of the cell sizes in the raster object (in km2), which differ from north to south. In order to calculate the entire area's size can be obtained by adding up all cell sizes or multiplying the median cell size by the number of cells (which I did here).

What is Geomerge?

geomerge is a framework for geospatial data integration that merges raster, spatial polygon, and (dynamic) spatial points data into a spatial (panel) data frame at any geographical resolution. Details. The geomerge function conducts a series of spatial joins for Geographic Information Systems (GIS) data.

What is a SpatialPolygonsDataFrame?

Spatial polygons are a set of spatially explicit shapes/polygons that represent a geographic location. Spatial polygons are composed of vertices which are a set of a series of x and y coordinates or Spatial points. Spatial polygons can be combined with data frames to create what's called a SpatialPolygonsDataFrame .


2 Answers

You could first use raster::clump() to identify clusters of connected raster cells and then apply rasterToPolygons() to "polygonize" those cells. (Do note, though, that each clump's area can be computed directly from the RasterLayer without converting it to a SpatialPolygonsDataFrame, as shown below):

library(rgeos) ## For the function gArea

## Clump and polygonize
Rclus <- clump(r)
SPclus <- rasterToPolygons(Rclus, dissolve=TRUE)

## Check that this works
plot(SPclus, col = seq_along(SPclus))

## Get cluster areas from RasterLayer object
transform(data.frame(freq(Rclus)), 
          area = count*prod(res(Rclus)))

## Get cluster areas from SpatialPolygons object
transform(data.frame(SPclus), 
          area = gArea(SPclus, byid=TRUE))

enter image description here

like image 53
Josh O'Brien Avatar answered Nov 04 '22 12:11

Josh O'Brien


The rgeos package has many polygon manipulation tools. gUnion will union together touching polygons:

require(rgeos)
uni <- gUnion( r_poly , r_poly )
plot( uni , col = 2 )

enter image description here

like image 23
Simon O'Hanlon Avatar answered Nov 04 '22 12:11

Simon O'Hanlon