Clustering points to polygons

Day 3 of 30DayMapChallenge
R
30DayMapChallenge
spatial
datavisualization
Author

Michaël

Published

2023-11-03

Modified

2023-11-03

A photo of a female lynx posing on a dead tree

Female lynx posing on the dead tree – CC-BY-ND by Tambako The Jaguar

Day 3 of 30DayMapChallenge: « Polygons » (previously). Using data from a previous post on lynx occurrences.

So we mapped the sightings of lynx, but there are a lot of outliers; how can we map the core habitat of the species?

We will use DBSCAN to cluster the points, then we will build a concave hull of these clusters. We arbitrary chose a minimal distance of 12 km between clusters and at least 5 points per clusters. Many coordinates are duplicated so we’ll get rid of artificial one-dimension polygons.

library(rgbif)
library(tidyverse)
library(sf)
library(leaflet)
library(dbscan)
library(concaveman)

# use existing data from GBIF
# keep valid and recent sightings in France.
lynx <- occ_download_get("0034730-231002084531237") |> 
  occ_download_import() |> 
  drop_na(decimalLongitude, decimalLatitude) |>
  filter(year >= "2015", 
         countryCode == "FR",
         occurrenceStatus == "PRESENT") |> 
  st_as_sf(coords = c("decimalLongitude", "decimalLatitude"), 
           crs = "EPSG:4326", 
           remove = FALSE)

# Create clusters
# we reproject to a Lambert projection to have units in meters
clusters <- lynx |> 
  st_transform("EPSG:2154") |> 
  st_coordinates() |> 
  dbscan(eps = 12000, minPts = 5)

# Create concave hulls
# we can chose an appropriate concavity parameter
concave_hulls <- lynx |> 
  bind_cols(tibble(cluster = clusters$cluster)) |>
  filter(cluster != 0) |> 
  group_by(cluster) |>
  summarise(nd = n_distinct(decimalLongitude, decimalLatitude),
            .groups = "drop") |> 
  filter(nd > 2) |> 
  mutate(geom = map(geometry, ~ concaveman(st_sf(st_sfc(.x)),
                                           concavity = 1.7))) |> 
  st_drop_geometry() |> 
  unnest(geom) |> 
  st_sf(crs = "EPSG:4326") |> 
  st_zm()
lynx |> 
  leaflet() |> 
  addCircleMarkers(radius = 0.1, opacity = 0.1) |> 
  addPolygons(data = concave_hulls, color = "red") |> 
  addTiles()
Figure 1: Lynx occurrence 2015-2023

We get 4 clusters: Jura, Vosges (North and South) and West of Dijon. We could estimate a rough lynx occupancy area in France:

concave_hulls |> 
  st_area() |> 
  sum() |> 
  units::set_units(km^2)
17925.6 [km^2]