7

Here is a simple example. I have a shapefile of three US states:

Three US states

This is the code that creates the shapefile.

library(USAboundaries)
library(sf)
library(dplyr)

states <- us_states(map_date = "2000-01-01", resolution = "high", states = c("CA", "OR", "WA")) %>%
  mutate(group_var = if_else(state_abbr == "CA", 1, 2))
plot(st_geometry(states))

Now, I'm trying to combine the geometries based on the value of group_var, so California should stand alone while Oregon and Washington get lumped together into a single geometry. Unfortunately, st_combine doesn't take a grouping variable, and although the aggregate function in sf looks promising, code like this throws an error that group_var cannot be found.

x <- aggregate(states, group_var, mean)

Furthermore, aggregate requires a function as a third argument, presumably because aggregate in the stats package applies the function to the data, but in this case, there isn't any data to apply a function to. I'm just trying to combine/aggregate the shapefiles.

9

You could also do this using dplyr's group_by() and summarize() functions:

states %>%
  group_by(group_var) %>% 
  summarize(geometry = st_union(geometry))
1
  • This is great; thank you. Since I'm more familiar with dplyr, this is much more intuitive to me. – Michael A Mar 24 '19 at 3:49
3

You were close; just cast the second argument in aggregate() as a list, like this:

x <- aggregate(states, 
               by = list(states$group_var),
               FUN = mean)

x
plot(x)
2
  • This is brilliant; thank you! Out of curiosity, what is the mean function actually taking the mean over if there is no non-geometric data in the data frame? – Michael A Mar 21 '19 at 17:41
  • Running x suggests that factors get turned into NA's. – JepsonNomad Mar 21 '19 at 18:10
0

You could also simplify Chris Merkord's answer with: states %>% group_by(group_var) %>% summarize()

0
0

You are basically describing a dissolve operation. This can easily be done using the rgeos package. Note that I switched datasets because I was having difficulty installing the USAboundaries and USAboundariesData packages.

library(rgeos)
library(sp)
library(spData)
library(rgdal)

columbus <- rgdal::readOGR(system.file("shapes/columbus.shp", 
                           package="spData")[1])
   print(columbus$NSB)
   plot(Columbus)

diss <- rgeos::gUnaryUnion(columbus, id = columbus$NSB)
  plot(diss, col=c("red","green"))
    plot(columbus, add=TRUE)

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.