3

I am not a pro with R nor spatial analysis. I am looking for a way to combine the polygons inside a spatial polygons data frame based on a field in the @data slot: the equivalent of dplyr's "group_by" for spdf's.

I'm not sure if merge, join, or combine are the right words but I hope it is clear what I'm looking for.

library(sp)

#coordinates:
xy1 = cbind(c(1,2,2,1),c(1,1,2,2))
xy2 = cbind(c(2,3,3,2),c(1,1,2,2))
xy3 = cbind(c(1,2,2,1),c(2,2,3,3))
xy4 = cbind(c(2,3,3,2),c(2,2,3,3))

#polygons:
p1 = Polygon(xy1)
ps1 = Polygons(list(p1),ID = "a")
p2 = Polygon(xy2)
ps2 = Polygons(list(p2),ID = "b")
p3 = Polygon(xy3)
ps3 = Polygons(list(p3),ID = "c")
p4 = Polygon(xy4)
ps4 = Polygons(list(p4),ID = "d")

#spatial polygons:
sps_m = SpatialPolygons(list(ps1,ps2,ps3,ps4))

#dataframe:
data_m = data.frame(dt = c("Group A","Group B","Group A","Group C"),row.names = c("a","b","c","d"))

#spatial polygons dataframe:
spdf_m = SpatialPolygonsDataFrame(sps_m,data_m)

#plot spdf:
plot(spdf_m)
Claudio Paladini
  • 1,000
  • 1
  • 10
  • 20
  • 2
    It sounds like you want a dissolve, which is the `gUnaryUnion()` function. At the risk of sounding shameless, I wrote about dissolving polygons recently: https://philmikejones.me/post/2015-09-03-dissolve-polygons-in-r/ – Phil Mar 13 '18 at 15:38
  • haha yes, thank you. in the mean time I also found: https://philmikejones.wordpress.com/2015/09/03/dissolve-polygons-in-r/ – Claudio Paladini Mar 13 '18 at 16:00
  • 1
    also aggregate() from the raster library: https://www.rdocumentation.org/packages/raster/versions/2.6-7/topics/aggregate – Claudio Paladini Mar 13 '18 at 16:13
  • 1
    @Phil thi is a great post, and you updated it to do everything the `sf` way. I highly agree with you on that! – Fitzroy Hogsflesh Oct 30 '18 at 12:46

1 Answers1

5

The correct term to describe the act of joining, combining, uniting or merging (making 1 out of multiple) polygons seems to be either dissolving or aggregating.

The funtion that works for me is aggregate() from the raster package.

library(sp)
library(raster)

#coordinates:
xy1 = cbind(c(1,2,2,1),c(1,1,2,2))
xy2 = cbind(c(2,3,3,2),c(1,1,2,2))
xy3 = cbind(c(1,2,2,1),c(2,2,3,3))
xy4 = cbind(c(2,3,3,2),c(2,2,3,3))

#polygons:
p1 = Polygon(xy1)
ps1 = Polygons(list(p1),ID = "a")
p2 = Polygon(xy2)
ps2 = Polygons(list(p2),ID = "b")
p3 = Polygon(xy3)
ps3 = Polygons(list(p3),ID = "c")
p4 = Polygon(xy4)
ps4 = Polygons(list(p4),ID = "d")

#spatial polygons:
sps_m = SpatialPolygons(list(ps1,ps2,ps3,ps4))

#dataframe:
data_m = data.frame(dt = c("Group A","Group B","Group A","Group C"),row.names = c("a","b","c","d"))

#spatial polygons dataframe:
spdf_m = SpatialPolygonsDataFrame(sps_m,data_m)

groups = aggregate(spdf_m, by = "dt")

#plot spdf:
plot(spdf_m)
plot(groups)
Claudio Paladini
  • 1,000
  • 1
  • 10
  • 20