Calculating distance between spatial point and spatial line

I have a simple features object with spatial points (coordinates) in it, and I have a large SpatialLinesDataframe (representing roads). I have about 9000 points, and I am wanting to calculate the shortest distance between each point and the nearest road in the SpatialLinesDataframe. As rgdal and rgeos are now deprecated, I am unsure how I can best calculate these distances efficiently (ie gDistance won't work). Any help is appreciated!

Solution

If it's OK for you to switch to sf package, you can combine or union your road network (presumably LINESTRING geometries) into a single MULTILINESTRING and use st_distance():

library(sf)
#> Linking to GEOS 3.11.2, GDAL 3.6.2, PROJ 9.2.0; sf_use_s2() is TRUE
# as an example, Roxel road network from sfnetworks package
roxel_ <- sfnetworks::roxel
roxel_
#> Simple feature collection with 851 features and 2 fields
#> Geometry type: LINESTRING
#> Dimension:     XY
#> Bounding box:  xmin: 7.522594 ymin: 51.94151 xmax: 7.546705 ymax: 51.9612
#> Geodetic CRS:  WGS 84
#> # A tibble: 851 × 3
#>    name                  type                                           geometry
#>  * <chr>                 <fct>                                  <LINESTRING [°]>
#>  1 Havixbecker Strasse   residential      (7.533722 51.95556, 7.533461 51.95576)
#>  2 Pienersallee          secondary   (7.532442 51.95422, 7.53236 51.95377, 7.53…
#>  3 Schulte-Bernd-Strasse residential (7.532709 51.95209, 7.532823 51.95239, 7.5…
#>  4 <NA>                  path        (7.540063 51.94468, 7.539696 51.94479, 7.5…
#>  5 Welsingheide          residential       (7.537673 51.9475, 7.537614 51.94562)
#>  6 <NA>                  footway     (7.543791 51.94733, 7.54369 51.94686, 7.54…
#>  7 <NA>                  footway           (7.54012 51.94478, 7.539931 51.94514)
#>  8 <NA>                  path        (7.53822 51.94546, 7.538131 51.94549, 7.53…
#>  9 <NA>                  track       (7.540063 51.94468, 7.540338 51.94468, 7.5…
#> 10 <NA>                  track       (7.5424 51.94599, 7.54205 51.94629, 7.5419…
#> # ℹ 841 more rows

# sample points 
roxel_bbox <- 
  roxel_ |>
  st_bbox() |>
  st_as_sfc()
sample_points <- st_sf(id = 1:20, geometry = st_sample(roxel_bbox, 20))

plot(st_geometry(roxel_))
plot(st_geometry(sample_points), add = TRUE, pch = 19, col = "red")


# convert linestrings to a single multilinestring
roxel_u <- st_union(roxel_)
roxel_u
#> Geometry set for 1 feature 
#> Geometry type: MULTILINESTRING
#> Dimension:     XY
#> Bounding box:  xmin: 7.522594 ymin: 51.94151 xmax: 7.546705 ymax: 51.9612
#> Geodetic CRS:  WGS 84
#> MULTILINESTRING ((7.533722 51.95556, 7.533461 5...

# add distance to a multilinestring for each point
sample_points$dist_to_road <- st_distance(sample_points, roxel_u)
sample_points
#> Simple feature collection with 20 features and 2 fields
#> Geometry type: POINT
#> Dimension:     XY
#> Bounding box:  xmin: 7.522862 ymin: 51.94221 xmax: 7.545885 ymax: 51.96073
#> Geodetic CRS:  WGS 84
#> First 10 features:
#>    id                  geometry   dist_to_road
#> 1   1   POINT (7.52293 51.9458) 191.023154 [m]
#> 2   2  POINT (7.537142 51.9531)  20.392243 [m]
#> 3   3 POINT (7.544911 51.94645)  76.427547 [m]
#> 4   4 POINT (7.523751 51.95759) 110.614399 [m]
#> 5   5 POINT (7.545885 51.96069) 173.888472 [m]
#> 6   6 POINT (7.523605 51.94263) 234.182639 [m]
#> 7   7 POINT (7.522862 51.95603) 101.840416 [m]
#> 8   8 POINT (7.525905 51.94221)  94.636355 [m]
#> 9   9  POINT (7.53068 51.95891)   3.990801 [m]
#> 10 10 POINT (7.525584 51.94721)   4.575924 [m]

^{Created on 2024-01-18 with reprex v2.0.2}