Clean polygons to remove slivers

clean_polygons_remove_slivers.md

Trying different methods to clean polygons to remove slivers

Load libraries

library(sf)
library(dplyr)
library(lwgeom)

# using planar geometry
sf::sf_use_s2(FALSE)

Original data for reproducible example

This is how the original data looks like:

# load data
states1920 <- geobr::read_state(year = 1920)
plot(states1920)

and here's the problem. Dissolving the borders leaves some slivers

states1920$br <- 1

states1920 |>
  group_by(br) |>
  summarise() |>
  plot(col='gray90')

Failed attempts

Here are a few solutions I've tried, all of them unsuccessful.

# using make_valid and buffer
states1920 |>
  sf::st_make_valid() |>
  sf::st_buffer(dist = 0) |>
  group_by(br) |>
  summarise() |>
  plot(col='gray90')


# using st_snap_to_grid
states1920 |>
  sf::st_make_valid() |>
  st_transform(crs = 32722) |>
  lwgeom::st_snap_to_grid(size = 0.01) |>
  group_by(br) |>
  summarise() |>
  plot(col='gray90')

Update: Solution

Here are a few solutions that work, which gives an output like this one below. Thanks to a few tweeps that came up with suggestions in this thread.

using sfheaders

states1920 |>
  sf::st_make_valid() |>
  group_by(br) |>
  summarise() |>
  sfheaders::sf_remove_holes() |>
  plot(col='gray90')

using nngeo

states1920 |>
  sf::st_make_valid() |>
  group_by(br) |>
  summarise() |>
  nngeo::st_remove_holes() |>
  plot(col='gray90')

using smoothr

This is the slowest solution

area_thresh <- units::set_units(800, km^2)

states1920 |>
  sf::st_make_valid() |>
  group_by(br) |>
  summarise() |>
  smoothr::fill_holes(threshold = area_thresh) |>
  plot(col = "gray90")

using rmapshaper

It works by simplifying geometries, which I would consider and unwanted side effect in this case

states1920 |>
  sf::st_make_valid() |>
  group_by(br) |>
  summarise() |>
  rmapshaper::ms_simplify() |>
  plot(col='gray90')

Using nngeo::st_remove_holes() removes the large slivers but not sure if it messes up some of the other borders along the coastline. You can experiment changing the max_area argument

library(geobr)
library(sf)
library(dplyr)

# load data
states1920 <- geobr::read_state(year = 1920)
states1920$br <- 1

test <- states1920 |>
  # optionally switch to a more appropriate CRS
  st_transform(crs = 32722) |>
  group_by(br) |>
  summarise() |>
  st_make_valid(.) |>
  #> Try to remove slivers
  nngeo::st_remove_holes() |>
  #> Back to 4674
  st_transform(crs = 4674) |>
  #> Just to make sure
  st_make_valid()

plot(test, col='gray90')