Package {pkmapr}

Title:	Pakistan Spatial Data Toolkit
Version:	1.2.1
Depends:	R (≥ 4.1.0)
Description:	Provides a tidy interface to Pakistan's official administrative boundary data from the United Nations Office for the Coordination of Humanitarian Affairs (OCHA). Downloads and caches spatial data at country, province, district, and tehsil levels as 'sf' objects compatible with the 'tidyverse' and geospatial ecosystem. Includes utilities for geographic dictionary lookup, coordinate reference system selection, spatial measurement, and neighbour structure construction for use with 'spdep', 'ggplot2', 'leaflet', and related packages.
License:	GPL-3
Encoding:	UTF-8
RoxygenNote:	7.3.3
Imports:	sf, dplyr, rlang, cli, spdep, jsonlite
Suggests:	ggplot2, leaflet, rmapshaper, httptest2, testthat (≥ 3.0.0), tmap, knitr, readr, rmarkdown, spelling
VignetteBuilder:	knitr
URL:	https://abdullahumer1101.github.io/pkmapr/, https://github.com/abdullahumer1101/pkmapr/
BugReports:	https://github.com/abdullahumer1101/pkmapr/issues
NeedsCompilation:	no
Packaged:	2026-05-08 10:25:36 UTC; abdul
Author:	Abdullah Umer [aut, cre]
Maintainer:	Abdullah Umer <abdullahumer1101@gmail.com>
Repository:	CRAN
Date/Publication:	2026-05-13 07:30:11 UTC

pkmapr: Pakistan Spatial Data Toolkit

Description

Provides a tidy interface to Pakistan's official administrative boundary data from the United Nations Office for the Coordination of Humanitarian Affairs (OCHA). Downloads and caches spatial data at country, province, district, and tehsil levels as sf objects compatible with the tidyverse and geospatial ecosystem.

Value

No return value, called for package documentation and namespace management.

Key functions

Data loading

get_districts, get_tehsils, get_provinces, get_country: Download administrative boundaries at different levels.

Lookup

pk_dictionary: Look up official names and PBS geocodes.

Joining

pk_join: Join external data with administrative units, with match checking.

Spatial weights

pk_neighbors: Build spatial weights for spdep and spatialreg.

Visualisation

pk_map, pk_map_interactive: Quick choropleth maps with ggplot2 or leaflet.

Package options

pkmapr uses the following global options:

pkmapr_use_cache

Logical. If TRUE, caches downloaded data in a persistent user directory. Default FALSE (uses temporary directory).

pkmapr_refresh

Logical. Re-download data even if cached copy exists. Default FALSE.

Author(s)

Maintainer: Abdullah Umer abdullahumer1101@gmail.com (ORCID)

References

United Nations Office for the Coordination of Humanitarian Affairs (OCHA). Pakistan Administrative Boundaries. https://data.humdata.org/

See Also

Useful links:

• Package website

• Source code

• Bug reports

Get national boundary for Pakistan.

Description

Get national boundary for Pakistan.

Usage

get_country(crs = 4326)

Arguments

Value

Returns an sf object (class "sf" and "data.frame") with:

The output represents the complete national boundary of Pakistan.

Examples

pakistan <- get_country()
plot(sf::st_geometry(pakistan))
print(pakistan$area_km2)

Get Pakistan district boundaries

Description

Get Pakistan district boundaries

Usage

get_districts(province = NULL, crs = 4326)

Arguments

Value

Returns an sf object (class "sf" and "data.frame") with:

When province is specified, the output contains only districts within that province. The output represents administrative boundaries at the district level.

Examples

# All districts
all_districts <- get_districts()

# Filter to Punjab province (case-insensitive)
punjab_districts <- get_districts(province = "Punjab")
punjab_districts <- get_districts(province = "punjab")  # Same result

plot(sf::st_geometry(punjab_districts))

Get Pakistan province boundaries

Description

Get Pakistan province boundaries

Usage

get_provinces(crs = 4326)

Arguments

Value

Returns an sf object (class "sf" and "data.frame") with:

The output represents the administrative boundaries of Pakistan's provinces and territories.

Examples

provinces <- get_provinces()
plot(sf::st_geometry(provinces))
head(provinces)

Get Pakistan tehsil boundaries

Description

Get Pakistan tehsil boundaries

Usage

get_tehsils(district = NULL, province = NULL, crs = 4326)

Arguments

Value

Returns an sf object (class "sf" and "data.frame") with:

The output represents the finest available administrative boundaries in pkmapr, suitable for high-resolution spatial analysis, local-level mapping, and joining with tehsil-level census or survey data.

Examples

# All tehsils
all_tehsils <- get_tehsils()

# Filter to Sindh province
sindh_tehsils <- get_tehsils(province = "Sindh")
sindh_tehsils <- get_tehsils(province = "sindh")  # Case-insensitive

# Filter to Lahore district
lahore_tehsils <- get_tehsils(district = "Lahore")
lahore_tehsils <- get_tehsils(district = "lahore")  # Case-insensitive

plot(sf::st_geometry(lahore_tehsils))

Recalculate area in km² for an sf object

Description

Computes accurate areas by reprojecting to UTM Zone 42N before measurement. Returns the input object with area_km2 added or updated.

Usage

pk_area(x)

Arguments

Value

Returns the input sf object (class "sf" and "data.frame") with the

Examples

districts <- get_districts()
districts <- pk_area(districts)
head(districts$area_km2)

Leaflet basemap centred on Pakistan

Description

Returns a leaflet map pre-zoomed to Pakistan's extent as a starting point for building custom interactive maps.

Usage

pk_basemap(provider = "CartoDB.Positron")

Arguments

Value

Returns a leaflet object (class "leaflet" and "htmlwidget") with the following configuration:

The output is a fully customizable leaflet basemap ready for adding layers (e.g., administrative boundaries, point data, markers). Use this as a starting point for building custom interactive maps.

Examples

# Default basemap with CartoDB.Positron tiles
pk_basemap()

# Alternative tile provider
pk_basemap(provider = "OpenStreetMap")

Get a bounding box for a named administrative unit

Description

Retrieves the spatial extent (bounding box) of a specific administrative unit by name and level. Useful for setting map views or cropping other spatial data.

Usage

pk_bbox(name, level = c("province", "district", "tehsil"))

Arguments

Value

Returns a bbox object (class "bbox") with named elements:

The output is suitable for use with ggplot2::coord_sf(xlim = c(xmin, xmax), ylim = c(ymin, ymax)) or leaflet::fitBounds(lng1 = xmin, lat1 = ymin, lng2 = xmax, lat2 = ymax). It represents the rectangular extent that exactly contains the requested administrative unit.

Note

If you see an error like ⁠object 'xxx' not found⁠ when using this function, the issue is likely in your data preparation, not pk_bbox(). Test the function directly: pk_bbox("Punjab", level = "province"). If that works, check that your data has the expected column names.

Examples

# Get bounding box for Lahore district
bb <- pk_bbox("Lahore", level = "district")
print(bb)

# Use with ggplot2

  library(ggplot2)
  districts <- get_districts()
  bb_punjab <- pk_bbox("Punjab", level = "province")
  ggplot() +
    geom_sf(data = districts) +
    coord_sf(xlim = c(bb_punjab["xmin"], bb_punjab["xmax"]),
             ylim = c(bb_punjab["ymin"], bb_punjab["ymax"]))

Create buffers around sf geometries in km

Description

Reprojects to UTM Zone 42N internally so buffer distances are in kilometres, then returns buffers in the original CRS.

Usage

pk_buffer(x, dist_km)

Arguments

Value

Returns an sf object with the same attributes as x but with geometries transformed to buffers of radius dist_km kilometres.

The output CRS is identical to the input CRS (reprojected back from UTM Zone 42N after buffering). This ensures buffers are circular in projected space with accurate kilometre distances.

Examples

districts <- get_districts()
buffered <- pk_buffer(districts, dist_km = 10)

Extract centroids from an sf object

Description

Returns polygon centroids as a point sf object with only geometry (attributes are dropped to avoid warnings about constant attributes).

Usage

pk_centroid(x)

Arguments

Value

Returns an sf point object (class "sf") with:

Centroids represent the geometric center of each polygon, useful for point-based visualisation, distance calculations, or as nodes for spatial network analysis.

Examples

districts <- get_districts()
centres <- pk_centroid(districts)
plot(centres)

Suggest an appropriate projected CRS for a Pakistan sf object

Description

Examines the geographic extent of an sf object and recommends the most appropriate projected coordinate reference system for metric operations. Pakistan spans UTM zones 41N, 42N, and 43N; national-level analyses benefit from an equal-area projection.

Usage

pk_crs_suggest(x)

Arguments

Details

Using WGS84 (EPSG:4326) for distance, area, or buffer operations produces inaccurate results as it measures in degrees rather than metres.

Value

A named list:

epsg

Integer EPSG code for the recommended CRS.

name

Human-readable CRS name.

rationale

One-sentence explanation of the recommendation.

Examples

  pk_crs_suggest(get_country())
  pk_crs_suggest(get_districts(province = "Balochistan"))

Pakistan Administrative Boundaries Dictionary

Description

Search and explore administrative names and codes for Pakistan.

Usage

pk_dictionary(level = NULL, name = NULL, code = NULL)

Arguments

Value

Returns a data frame (class "data.frame") with the following columns:

The output represents a searchable dictionary of all administrative units in Pakistan. Use this function to explore available units, find codes, or validate input for other pkmapr functions.

Case Insensitivity

All matching in pk_dictionary() is case-insensitive. "Lahore", "lahore", and "LAHORE" all return the same results.

Examples

# All provinces
pk_dictionary(level = "provinces")

# Case-insensitive search for districts containing "lahore"
pk_dictionary(level = "districts", name = "lahore")
pk_dictionary(level = "districts", name = "LAHORE")  # Same result

# Search by code
pk_dictionary(code = "PK6")

Compute distances between two sf objects in km

Description

Reprojects internally to UTM Zone 42N for accurate metric distances.

Usage

pk_distance(x, y, by = c("centroid", "edge"))

Arguments

Value

Returns a numeric matrix (class "matrix") of distances in kilometres, with dimensions nrow(x) by nrow(y).

When by = "centroid", distances are measured between polygon centroids. When by = "edge", distances are measured between the closest points on polygon boundaries. The output represents the shortest straight-line distance between features.

Examples

provinces <- get_provinces()
d <- pk_distance(provinces, provinces)
print(d)

# Edge-to-edge distances
d_edge <- pk_distance(provinces, provinces, by = "edge")

Intersect two sf objects

Description

Returns the geometric intersection of two sf objects with CRS alignment handled automatically.

Usage

pk_intersect(x, y)

Arguments

Value

Returns an sf object containing the geometric intersection of x and y:

The output CRS matches x. Empty geometries (no intersection) are dropped. Output represents areas/features common to both input layers.

Examples

districts <- get_districts()
buffered <- pk_buffer(districts, dist_km = 10)
intersected <- pk_intersect(districts, buffered)

Join data to a pkmapr sf object with match checking

Description

Performs a left join of external data to a pkmapr spatial object, with automatic validation of matching keys. Uses code columns (e.g., district_code, tehsil_code, province_code) wherever possible to ensure reliable joins even when names change or have spelling variations.

Usage

pk_join(spatial, data, by)

Arguments

Value

Returns the spatial sf object (class "sf" and "data.frame") with all columns from data joined to the matching rows.

The output preserves:

Rows in data that do not match any spatial unit generate a warning and receive NA values for spatial attributes in the joined result. Rows in the spatial object that have no match in data retain their original attributes but receive NA for the joined data columns.

Note

After joining, always inspect names(result) to check for column name conflicts. If your data shares column names with the spatial object (e.g., province_name, district_name), both versions will be preserved with .x and .y suffixes. Rename or select the appropriate columns before further analysis.

Examples

districts <- get_districts()
my_data <- data.frame(district_code = "PK603", value = 42)
joined <- pk_join(districts, my_data, by = "district_code")
print(names(joined))

# Example with missing match (generates warning)

  bad_data <- data.frame(district_code = c("PK603", "INVALID_CODE"), value = c(42, 99))
  joined_bad <- pk_join(districts, bad_data, by = "district_code")

Produces a ggplot2 map for rapid exploratory visualisation. Returns a ggplot object that can be extended with standard ggplot2 layers.

Description

Produces a ggplot2 map for rapid exploratory visualisation. Returns a ggplot object that can be extended with standard ggplot2 layers.

Usage

pk_map(x, fill = NULL, title = NULL, ...)

Arguments

Value

Returns a ggplot object (class "gg" and "ggplot") representing a choropleth map.

When fill = NULL, the output is an outline map with grey90 fill and white borders, useful for context or reference.

When a fill variable is provided, the output uses a viridis color scale with automatic legend, for visualizing spatial distributions of continuous variables (e.g., area, population, density).

The returned ggplot object can be extended with additional layers, themes, or scales using standard ggplot2 syntax.

Examples

# Outline map of provinces
pk_map(get_provinces())

# Choropleth map with fill variable
pk_map(get_provinces(), fill = "area_km2", title = "Province areas")

Interactive choropleth map of a pkmapr sf object

Description

Produces an interactive leaflet map. Returns a leaflet object that can be extended with standard leaflet functions.

Usage

pk_map_interactive(x, fill = NULL, popup = NULL, ...)

Arguments

Value

Returns a leaflet object (class "leaflet" and "htmlwidget") representing an interactive web map.

When fill = NULL, the output shows polygon outlines only.

When a fill variable is provided, the output renders polygons with:

The map includes a legend for the fill variable and uses the CartoDB.Positron tile provider as the basemap.

Examples

districts <- get_districts()
pk_map_interactive(districts,
                   fill = "area_km2",
                   popup = c("district_name", "area_km2"))

Construct a spatial neighbours list for Pakistan administrative units

Description

Builds a contiguity or distance-based spatial neighbours structure for direct use with spdep and spatialreg. Handles Pakistan-specific complexities including non-contiguous units and disputed boundaries.

Usage

pk_neighbors(
  x,
  style = c("queen", "rook", "knn"),
  k = NULL,
  disputed = c("exclude", "include", "flag")
)

Arguments

Value

Returns a named list (class "list") with the following:

The output is a complete spatial weights structure for use in spatial autocorrelation tests (Moran's I), local indicators of spatial association (LISA), and spatial regression models (SAR, SEM, etc.). The nb defines the neighbour graph; the listw provides the row-standardized weights matrix.

Pakistan-specific handling

Gilgit-Baltistan and Azad Kashmir might break some spatial statistics. The disputed argument controls how the Line of Control and special administrative boundaries flagged in the OCHA source data are treated, offering better control over analytical decisions.

Examples

  districts <- get_districts()
  w <- pk_neighbors(districts)

  # Calculate Moran's I using spdep
  moran_result <- spdep::moran.test(districts$area_km2, w$listw)
  print(moran_result)

  # Queen contiguity (default)
  w_queen <- pk_neighbors(districts, style = "queen")

  # K-nearest neighbours (k=5)
  w_knn <- pk_neighbors(districts, style = "knn", k = 5)

  # Flag disputed boundaries for documentation
  w_flagged <- pk_neighbors(districts, disputed = "flag")
  if (!is.null(w_flagged$boundary_note)) cat(w_flagged$boundary_note)

Assign points to administrative units (point-in-polygon)

Description

For each point in points, identifies which polygon it falls within and joins that polygon's attributes to the point record.

Usage

pk_points_in(points, polygons, return_all = TRUE)

Arguments

Value

Returns an sf object of the same class as points with:

Points that fall outside all polygons receive NA values for all polygon attributes when return_all = TRUE. When return_all = FALSE, such points are removed entirely.

The output represents point locations with their containing administrative unit attributes attached, allowing spatial aggregations and unit-based analyses.

Examples

# Get district boundaries
districts <- get_districts()

# Create sample points (or use your own sf object)
set.seed(123)
sample_points <- sf::st_sample(districts, size = 50)
sample_points_sf <- sf::st_sf(geometry = sample_points)

# Assign points to districts
points_with_districts <- pk_points_in(sample_points_sf, districts)
print(head(points_with_districts))

Project an sf object to a Pakistan-appropriate CRS

Description

Convenience wrapper around sf::st_transform() with a default of UTM Zone 42N, appropriate for most Pakistan analyses requiring accurate distance, area, or buffer operations.

Usage

pk_project(x, crs = 32642)

Arguments

Value

Returns the sf object reprojected to the specified CRS. The output has the same attributes and geometry type as the input, but coordinates are transformed to the new projection. UTM Zone 42N (EPSG:32642) preserves distances and areas accurately across most but not all of Pakistan.

Examples

districts <- get_districts()
projected <- pk_project(districts)
sf::st_crs(projected)$epsg  # Should be 32642

Search Across All Administrative Levels

Description

Search for administrative units by partial name or code across all levels (provinces, districts, tehsils).

Usage

pk_search(query, fuzzy = FALSE)

Arguments

Value

Returns a data frame (class "data.frame") with the following columns:

The output represents all administrative units that match the search query. When no matches are found, returns invisible NULL with a warning message. When fuzzy matching is enabled, the output may include approximate matches that could be useful for handling typos or spelling variations.

Case Insensitivity

All matching in pk_search() is case-insensitive by default. "Lahore", "lahore", and "LAHORE" all return the same results.

Fuzzy Matching

When fuzzy = TRUE, the function uses approximate string matching to handle typos and spelling variations. For example, "lahor" will match "Lahore". This is useful for interactive exploration but may return unexpected results for ambiguous queries.

Examples

# Case-insensitive search
pk_search("lahore")   # Returns Lahore district and tehsils
pk_search("LAHORE")   # Same result

# Fuzzy search for misspelled "lahore"
pk_search("lahor", fuzzy = TRUE)
pk_search("lahre", fuzzy = TRUE)

# Search by code
pk_search("PK6")

Dissolve sf polygons by a grouping column

Description

Merges polygons sharing the same value in a grouping column and recalculates area_km2.

Usage

pk_union(x, by)

Arguments

Value

Returns a dissolved sf object (class "sf" and "data.frame") with:

Polygons that are not spatially adjacent but share the same group value will become MULTIPOLYGON objects. Invalid geometries are repaired automatically using sf::st_make_valid().

Examples

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

Check package version and update status

Description

Compares the installed version of pkmapr with the latest release available on GitHub.

Usage

pk_version(quiet = FALSE)

Arguments

Value

Returns invisibly a list (class "list") with the following components:

Examples

pk_version()

# Silent mode for programmatic use
vers <- pk_version(quiet = TRUE)
if (!is.na(vers$latest) && vers$installed != vers$latest) {
  message("Update recommended!")
}

Deprecated functions in pkmapr

Description

These functions were deprecated. Please use the new pk_* functions instead.

Usage

pak_crs_suggest(...)

pak_dictionary(...)

pak_join(...)

pak_area(...)

pak_distance(...)

pak_centroid(...)

pak_bbox(...)

pak_project(...)

pak_neighbors(...)

pak_points_in(...)

pak_buffer(...)

pak_union(...)

pak_intersect(...)

pak_map(...)

pak_map_interactive(...)

pak_basemap(...)

Value

No return value, called for side effects (deprecation warnings).