CODEWITHSUNDEEP
rmapsspatialr-sf
Alexa Fredston
Alexa Fredston

Reputation: 727

Draw polygons around coordinates in R

I have clusters of points from different datasets collected in different geographical regions. I'd like to plot the spatial extent of each region on a map, for which I need a single polygon per region. I can find commands to convert points to polygons by "connecting the dots", but none by drawing an outline (the geometric method isn't terribly important since this is just for plotting). Here's a reprex to generate a simplified dataframe like the one I'm using, for points in Norway, Alaska, and the continental US:

nor_lat <- c(72.117, 71.05, 71.717, 71.167, 71.817, 72.333, 71.25, 70.917, 
  70.85, 70.933)
nor_lon <- c(33.217, 30.333, 26.3, 22.333, 22.333, 22.333, 18.517, 20.133, 
             21.333, 19.917)
us_lat <- c(41.283333, 41, 40.95, 40.95, 38.683333, 40.783333, 40.733333, 
            40.516667, 38.566667, 41.266667)
us_lon <- c(-71.116667, -70.733333, -70.583333, -70.483333, -74.816667, 
            -70.5, -70.566667, -70.266667, -74.85, -71.366667)
ak_lat <- c(58.015, 57.67167, 57.33833, 56.685, 56.99333, 57.31667, 57.65, 
            57.99667, 58.32167, 58.33333)
ak_lon <- c(-158.31667, -158.36, -158.41, -159.76, -159.72333, -159.66333, 
            -159.64167, -159.605, -159.54, -160.72167)
nor_dat <- data.frame(lat=nor_lat, lon=nor_lon)
us_dat <- data.frame(lat=us_lat, lon=us_lon)
ak_dat <- data.frame(lat=ak_lat, lon=ak_lon)
nor_dat$region <- "nor"
us_dat$region <- "us"
ak_dat$region <- "ak"
dat <- rbind(nor_dat, us_dat, ak_dat)

For the object dat, I'd like to generate a simple polygon around the points in each region, using sf if possible. (I'm not too worried about figuring out how to iterate over the regions, but I can't even figure out how to draw a polygon around one set of points!)

Upvotes: 0

Views: 5303

Answers (2)

Skaqqs
Skaqqs

Reputation: 4140

library(sf)
library(dplyr)

hulls <- dat %>%
  st_as_sf(coords = c("lon", "lat")) %>%
  group_by(region) %>%
  summarize(geometry = st_union(geometry)) %>%
  st_convex_hull()
plot(hulls)

enter image description here

Upvotes: 3

lovalery
lovalery

Reputation: 4652

One approach using sf and dplyr for "geocomputation" and filtering, tmap and tmaptools for visualization.

Please find the following reprex.

Reprex

  • Convert the dataframe "dat" into a sf object
library(sf)
library(dplyr)
library(tmap)
library(tmaptools)

sf_use_s2(TRUE)

dat_sf <- st_as_sf(dat, coords = c("lon", "lat"), crs = 4326)
  • Compute a smoothed enveloping polygon for each set of points in each area
conv_nor <- st_simplify(st_buffer(st_convex_hull(st_union(st_geometry(filter(dat_sf, region == "nor")))), dist = 15000), dTolerance = 5000)
conv_us <- st_simplify(st_buffer(st_convex_hull(st_union(st_geometry(filter(dat_sf, region == "us")))), dist = 15000), dTolerance = 5000)
conv_ak <- st_simplify(st_buffer(st_convex_hull(st_union(st_geometry(filter(dat_sf, region == "ak")))), dist = 15000), dTolerance = 5000)
  • Visualization: Overview
tm_shape(dat_sf)+
  tm_dots(size = 0.2, col = "black", shape = 21, alpha = 0.3)+
  tm_shape(conv_nor)+
  tm_borders(col = "gray", lwd = 2)+
  tm_shape(conv_us)+
  tm_borders(col = "blue", lwd = 2)+
  tm_shape(conv_ak)+
  tm_borders(col = "gray", lwd = 2)+
  tm_layout(frame = FALSE)

  • Visualization: detailed views
tm_shape(filter(dat_sf, region == "nor"), bbox = bb(conv_nor, ext = 1))+
  tm_dots(size = 0.5, alpha = 0.5)+
  tm_shape(conv_nor)+
  tm_borders(col = "gray", lwd = 2)+
  tm_layout(frame = FALSE, title = "'nor' area", 
            title.size = 1.5, title.position = c(0.85, "bottom"))


tm_shape(filter(dat_sf, region == "us"), bbox = bb(conv_us, ext = 1))+
  tm_dots(size = 0.5, alpha = 0.5)+
  tm_shape(conv_us)+
  tm_borders(col = "gray", lwd = 2)+
  tm_layout(frame = FALSE, title = "'us' area", 
            title.size = 1.5, title.position = c(0., "top"))


tm_shape(filter(dat_sf, region == "ak"), bbox = bb(conv_ak, ext = 1))+
  tm_dots(size = 0.5, alpha = 0.5)+
  tm_shape(conv_ak)+
  tm_borders(col = "gray", lwd = 2)+
  tm_layout(frame = FALSE, title = "'ak' area", 
            title.size = 1.5, title.position = c(0.75, "top"))

  • Detailed views with background

Finally, if you want some contextual background, you can use the tmap_mode("view") which call leaflet - apparently you work on the marine environment ;-)

tmap_mode("view")
#> tmap mode set to interactive viewing
tm_shape(filter(dat_sf, region == "nor"), bbox = bb(conv_nor, ext = 1))+
  tm_dots(size = 0.5, alpha = 0.5)+
  tm_shape(conv_nor)+
  tm_borders(col = "gray", lwd = 2)+
  tm_layout(frame = FALSE, title = "'nor' area", 
            title.size = 1.5, title.position = c(0.85, "bottom"))


tm_shape(filter(dat_sf, region == "us"), bbox = bb(conv_us, ext = 1))+
  tm_dots(size = 0.5, alpha = 0.5)+
  tm_shape(conv_us)+
  tm_borders(col = "gray", lwd = 2)+
  tm_layout(frame = FALSE, title = "'us' area", 
            title.size = 1.5, title.position = c(0., "top"))


tm_shape(filter(dat_sf, region == "ak"), bbox = bb(conv_ak, ext = 1))+
  tm_dots(size = 0.5, alpha = 0.5)+
  tm_shape(conv_ak)+
  tm_borders(col = "gray", lwd = 2)+
  tm_layout(frame = FALSE, title = "'ak' area", 
            title.size = 1.5, title.position = c(0.75, "top"))


tmap_mode("plot")
#> tmap mode set to plotting

Created on 2021-10-20 by the reprex package (v2.0.1)

Upvotes: 1

Related Questions