Skip to contents

Use case: working with spatial R packages

One may find themselves having to clean up messy coordinates as part of their project/work/etc. How does this look when fit into a workflow going all the way to visualization?

Let’s say you have the following messy coordinates that you’ve compiled from different places, leading to a variety of messy formats:

lats <- c(
  "46.4183",
  "46.4383° N",
  "46.5683° N",
  "46° 27´ 5.4\" N",
  "46° 25.56’",
  "N46°24’4.333"
)
lons <- c(
  "25.7391",
  "E25°34’6.4533",
  "25.3391° E",
  "25.8391° E",
  "25° 35.56’",
  "E25°34’4.333"
)

Parse messy coordinates

dat <- tibble::tibble(
  longitude = parse_lon(lons),
  latitude = parse_lat(lats)
)
dat
#> # A tibble: 6 × 2
#>   longitude latitude
#>       <dbl>    <dbl>
#> 1      25.7     46.4
#> 2      25.6     46.4
#> 3      25.3     46.6
#> 4      25.8     46.5
#> 5      25.6     46.4
#> 6      25.6     46.4

Combine coordinates with other data

dat$shape <- c("round", "square", "triangle", "round", "square", "square")
dat$color <- c("blue", "yellow", "green", "red", "green", "yellow")
dat
#> # A tibble: 6 × 4
#>   longitude latitude shape    color 
#>       <dbl>    <dbl> <chr>    <chr> 
#> 1      25.7     46.4 round    blue  
#> 2      25.6     46.4 square   yellow
#> 3      25.3     46.6 triangle green 
#> 4      25.8     46.5 round    red   
#> 5      25.6     46.4 square   green 
#> 6      25.6     46.4 square   yellow

Coerce to an sf object

datsf <- sf::st_as_sf(dat, coords = c("longitude", "latitude"))
datsf
#> Simple feature collection with 6 features and 2 fields
#> Geometry type: POINT
#> Dimension:     XY
#> Bounding box:  xmin: 25.3391 ymin: 46.4012 xmax: 25.8391 ymax: 46.5683
#> CRS:           NA
#> # A tibble: 6 × 3
#>   shape    color            geometry
#>   <chr>    <chr>             <POINT>
#> 1 round    blue    (25.7391 46.4183)
#> 2 square   yellow (25.56846 46.4383)
#> 3 triangle green   (25.3391 46.5683)
#> 4 round    red     (25.8391 46.4515)
#> 5 square   green   (25.59267 46.426)
#> 6 square   yellow (25.56787 46.4012)

Calculate the center of the plot view

center_lon <- mean(dat$longitude)
center_lat <- mean(dat$latitude)

Plot data using the leaflet package

if (!requireNamespace("leaflet")) install.packages("leaflet")
library("leaflet")
leaflet() %>%
  addTiles() %>%
  addMarkers(data = datsf) %>%
  setView(center_lon, center_lat, zoom = 10)

We’d like to have data only for a certain area, e.g., a political boundary or a park boundary. We can clip the data to a bounding box using sf::st_crop().

First, define the bounding box, and visualize

bbox <- c(
  xmin = 25.42813, ymin = 46.39455,
  xmax = 25.68769, ymax = 46.60346
)
leaflet() %>%
  addTiles() %>%
  addRectangles(bbox[["xmin"]], bbox[["ymin"]], bbox[["xmax"]], bbox[["ymax"]]) %>%
  setView(center_lon, center_lat, zoom = 10)

Crop the data to the bounding box

datsf_c <- st_crop(datsf, bbox)

Plot the new data

leaflet() %>%
  addTiles() %>%
  addMarkers(data = datsf_c) %>%
  setView(center_lon, center_lat, zoom = 10)