Spatial Analysis with pkmapr

library(pkmapr)
library(dplyr)

This vignette covers the spatial analysis utilities in pkmapr, including centroids, buffers, distance calculations, point-in-polygon assignment, boundary dissolving, and coordinate reference system selection.

Centroids

pk_centroid() converts polygon geometries to points, useful for labelling and distance calculations:

districts <- get_districts()
centroids <- pk_centroid(districts)

pk_map(districts) +
  ggplot2::geom_sf(data = centroids, color = "red", size = 0.5)

Buffers

pk_buffer() creates buffer zones around administrative units. Distances are specified in kilometres:

lahore <- get_districts(province = "Punjab") |>
  filter(district_name == "Lahore")

lahore_buffer <- pk_buffer(lahore, dist_km = 10)

pk_map(lahore_buffer, title = "10 km buffer around Lahore") +
  ggplot2::geom_sf(data = lahore, fill = "red", alpha = 0.5)

Distance calculations

pk_distance() computes centroid-to-centroid distances between two sf objects, returning a distance matrix:

# Distance matrix between all provinces
provinces <- get_provinces()
dist_matrix <- pk_distance(provinces, provinces)

# Distance from each province centroid to Karachi
karachi <- get_districts(province = "Sindh") |>
  filter(district_name == "Karachi")

distances <- pk_distance(provinces, karachi)

Point-in-polygon

pk_points_in() assigns GPS point locations to the administrative unit they fall within:

facilities <- data.frame(
  name = c("Hospital A", "Clinic B"),
  lon  = c(74.3, 74.5),
  lat  = c(31.5, 31.6)
) |>
  sf::st_as_sf(coords = c("lon", "lat"), crs = 4326)

facilities_with_district <- pk_points_in(facilities, districts)

facilities_with_district |>
  sf::st_drop_geometry() |>
  select(name, district_name)

Dissolve boundaries

pk_union() aggregates finer administrative units to a coarser level by dissolving shared boundaries:

tehsils <- get_tehsils()
districts_from_tehsils <- pk_union(tehsils, by = "district_name")

# Areas should be approximately equal
pk_area(get_districts())
pk_area(districts_from_tehsils)

Coordinate reference systems

pkmapr returns data in WGS84 (EPSG:4326) by default, which measures in degrees. For metric operations such as area, buffer, and distance calculations, a projected CRS is required. Use pk_crs_suggest() to get a recommended CRS for your data’s spatial extent, then reproject with pk_project():

# Get a recommended projected CRS
pk_crs_suggest(get_districts(province = "Punjab"))

# Reproject for metric operations
districts_utm <- pk_project(districts, crs = 32642)