Introduction to distance-to

Overview

The distanceto package is designed to quickly sample distances from points features to other vector layers. Normally the approach for calculating distance to (something) involves generating distance surfaces using raster based approaches eg. raster::distance or gdal_proximity and subsequently point sampling these surfaces. Since raster based approaches are a costly method that frequently leads to memory issues or long and slow run times with high resolution data or large study sites, we have opted to compute these distances using vector based approaches. As a helper, there’s a decidedly low-res raster based approach for visually inspecting your region’s distance surface. But the workhorse is distance_to.

The distanceto package provides two functions:

  • distance_to
  • distance_raster

Install

# Enable the robitalec universe
options(repos = c(
    robitalec = 'https://robitalec.r-universe.dev',
    CRAN = 'https://cloud.r-project.org'))

# Install distanceto
install.packages('distanceto')

distance_to

library(distanceto)
library(sf)
#> Linking to GEOS 3.12.1, GDAL 3.8.4, PROJ 9.4.0; sf_use_s2() is TRUE

Long-lat / unprojected coordinates

# Load nc data
nc <- st_read(system.file("shape/nc.shp", package="sf"))
#> Reading layer `nc' from data source `/github/workspace/pkglib/sf/shape/nc.shp' using driver `ESRI Shapefile'
#> Simple feature collection with 100 features and 14 fields
#> Geometry type: MULTIPOLYGON
#> Dimension:     XY
#> Bounding box:  xmin: -84.32385 ymin: 33.88199 xmax: -75.45698 ymax: 36.58965
#> Geodetic CRS:  NAD27

# Set number of sampling points
npts <- 1e3

# Sample points in nc
ncpts <- st_sample(nc, npts)

# Select first 5 of nc
ncsub <- nc[1:5,]

# Measure distance from ncpts to first 5 of nc
dist <- distance_to(ncpts, ncsub, measure = 'geodesic')

# or add to ncpts
ncpts$dist <- dist

head(dist, 30)
#>  [1] 220779.24 239242.52 202437.43 162589.60 104233.70 243555.10  50995.95
#>  [8] 128209.23 234472.77  81482.86  24700.94  76465.55 110375.44  40099.31
#> [15] 135443.94 102143.42 162993.94  39967.70 116564.27  73480.90 113867.25
#> [22]  77448.39 126882.11  36999.52 242041.27  84438.19  72914.88 105166.21
#> [29]  36784.36  92904.91
hist(dist)

Projected coordinates

# Transform nc data to local projected coordinates (UTM 18N)
nc_utm <- st_transform(nc, 32618)

# Set number of sampling points
npts <- 1e2

# Sample points within nc data
nc_utm_pts <- st_sample(nc_utm, npts)

# Select one polygon within nc data
nc_utm_select <- nc_utm[1, ]

# Measure distance from seine points to seine
dist <- distance_to(nc_utm_pts, nc_utm_select)

# or add to seine points
nc_utm_pts$dist <- dist

head(dist, 30)
#>  [1] 167679.3958 227347.6486 278871.2164 216317.9648 150411.0908    934.0783
#>  [7] 328271.2004 266498.2327 368821.5136 312042.3773 430028.6610 119176.1696
#> [13] 367558.4249 114690.1906 381484.9658 336960.5982 204421.9242  38008.6047
#> [19] 467824.1236 189040.8401  35137.2760   6534.9270 242661.4757  84836.0838
#> [25] 121200.3371 237927.6456 130448.4973  17606.9104      0.0000 286140.1684
hist(dist)

distance_raster

library(raster)
#> Loading required package: sp

rdist <- distance_raster(nc_utm_select, 1e4, extent = st_bbox(nc_utm))

plot(rdist)