Calculating driving distances from shapefiles in R

I am working with the R programming language.

I am interested in learning about how to calculate the driving distance (e.g. based on road networks) between two sets of coordinates.

For example:

• CN Tower: 290 Bremner Blvd, Toronto, ON M5V 3L9 (-79.61203, 43.68312)
• Toronto Airport: 6301 Silver Dart Dr, Mississauga, ON L5P 1B2 (-79.61203, 43.64256)

To solve this problem, I tried to download the shapefile for the Canadian Road Network and subset it for the Province of Ontario:

``````library(sf)
library(rgdal)
library(sfnetworks)
# Set the URL for the shapefile
url <- "https://www12.statcan.gc.ca/census-recensement/2011/geo/RNF-FRR/files-fichiers/lrnf000r22a_e.zip"

temp_dir <- tempdir()
temp_file <- file.path(temp_dir, "lrnf000r22a_e.zip")

unzip(temp_file, exdir = temp_dir)

# source: https://gis.stackexchange.com/questions/456748/strategies-for-working-with-large-shapefiles/456798#456798
a = st_read(file.path(temp_dir, "lrnf000r22a_e.shp"), query="select * from lrnf000r22a_e where PRUID_R ='35'")
``````

When imported into R, the file looks something like this:

``````Simple feature collection with 570706 features and 21 fields
Geometry type: LINESTRING
Dimension:     XY
Bounding box:  xmin: 5963148 ymin: 665490.8 xmax: 7581671 ymax: 2212179
First 10 features:
OBJECTID NGD_UID       NAME TYPE  DIR AFL_VAL ATL_VAL AFR_VAL ATR_VAL CSDUID_L                           CSDNAME_L CSDTYPE_L CSDUID_R                           CSDNAME_R CSDTYPE_R PRUID_L PRNAME_L PRUID_R PRNAME_R RANK
1         4 3434819       <NA> <NA> <NA>    <NA>    <NA>    <NA>    <NA>  3526003                           Fort Erie         T  3526003                           Fort Erie         T      35  Ontario      35  Ontario    5
2         5 1358252      South LINE <NA>    <NA>    <NA>    <NA>    <NA>  3551027                Gordon/Barrie Island        MU  3551027                Gordon/Barrie Island        MU      35  Ontario      35  Ontario    5
3         8 1778927       <NA> <NA> <NA>    <NA>    <NA>    <NA>    <NA>  3512054                           Wollaston        TP  3512054                           Wollaston        TP      35  Ontario      35  Ontario    5
4        11  731932     Albion   RD <NA>    <NA>    <NA>    2010    2010  3520005                             Toronto         C  3520005                             Toronto         C      35  Ontario      35  Ontario    3
5        18 3123617  County 41   RD <NA>     640     708     635     709  3511015                     Greater Napanee         T  3511015                     Greater Napanee         T      35  Ontario      35  Ontario    3
6        20 4814160       <NA> <NA> <NA>    <NA>    <NA>    <NA>    <NA>  3553005     Greater Sudbury / Grand Sudbury        CV  3553005     Greater Sudbury / Grand Sudbury        CV      35  Ontario      35  Ontario    5
7        21 1817031       <NA> <NA> <NA>    <NA>    <NA>    <NA>    <NA>  3537028                         Amherstburg         T  3537028                         Amherstburg         T      35  Ontario      35  Ontario    5
8        24 4825761       <NA> <NA> <NA>    <NA>    <NA>    <NA>    <NA>  3554094 Timiskaming, Unorganized, West Part        NO  3554094 Timiskaming, Unorganized, West Part        NO      35  Ontario      35  Ontario    5
9        25  544891     Dunelm   DR <NA>       1       9       2      10  3526053                      St. Catharines        CY  3526053                      St. Catharines        CY      35  Ontario      35  Ontario    5
10       28 1835384 Seven Oaks   DR <NA>     730     974     731     975  3515005             Otonabee-South Monaghan        TP  3515005             Otonabee-South Monaghan        TP      35  Ontario      35  Ontario    5
CLASS                 _ogr_geometry_
1     23 LINESTRING (7269806 859183,...
2     23 LINESTRING (6921247 1133452...
3     23 LINESTRING (7320857 1089403..
``````

Then, by consulting different references (e.g. https://cran.r-project.org/web/packages/sfnetworks/vignettes/sfn01_structure.html, https://stackoverflow.com/questions/66513615/creating-a-route-from-a-list-of-points-and-overlaying-the-route-list-of-points) - I tried to calculate the distance between these two points:

``````# convert the shapefile to an sfnetwork object
net <- as_sfnetwork(a)

# define your start and end points
q1 <- st_point(c(-79.61203, 43.68312))
q2 <- st_point(c(-79.38709, 43.64256))

# set the CRS of the points to match the CRS of your shapefile
q1 <- st_sfc(q1, crs = st_crs(a))
q2 <- st_sfc(q2, crs = st_crs(a))

# find the shortest path between the two points
path <- st_network_paths(net, q1, q2)

# calculate the distance of the path in meters
distance <- sum(st_length(path))
``````

But I get the following error:

Error in UseMethod("st_geometry") : no applicable method for 'st_geometry' applied to an object of class "c('tbl_df', 'tbl', 'data.frame')"

How can I fix this problem?

Is it also possible to visualize the results?

• You don't need to read the shapefile data using `readOGR` - you don't go on to use the created R object `shapefile`. You are only using the `sf` object `a`. Commented Apr 3, 2023 at 19:46
• What package supplies the network functions? Commented Apr 3, 2023 at 19:46
• @ Spacedman: thank you for your reply! The "sfnetwork" package is using the network functions. Commented Apr 4, 2023 at 0:53
• Do you mean `sfnetworks`? There's no `sfnetwork` on CRAN. Commented Apr 4, 2023 at 8:02
• Yes, that's what I meant. Thanks! Commented Apr 4, 2023 at 13:12

The return value from `sf_network_paths` is not an `sf` spatial object and so `st_length` will fail with this error. The help for `sf_network_paths` explains:

``````Value:

An object of class ‘tbl_df’ with one row per returned path.
Depending on the setting of the ‘type’ argument, columns can be
‘node_paths’ (a list column with for each path the ordered indices
of nodes present in that path) and ‘edge_paths’ (a list column
with for each path the ordered indices of edges present in that
path). ‘'all_shortest'’ and ‘'all_simple'’ return only
‘node_paths’, while ‘'shortest'’ returns both.
``````

Code in the help page for that function probably tells you how to extract the paths from the network and the returned value, reading and understanding the working examples is always a good idea.