Batch geocoding street intersections in R

Question

I am currently working with data that is formatted like this:

    tribble(
      ~street1, ~street2, ~county, ~state
      N BENTON WY, W TEMPLE ST, LOS ANGELES, CA,
      11TH PL, BLAINE ST, LOS ANGELES, CA,
      W 6TH ST, HOPE ST, LOS ANGELES, CA,
      S GRAND AV, W 18TH ST, LOS ANGELES, CA,
      BROADWAY, 5TH ST, LOS ANGELES, CA,
    )

This corresponds to a dataset containing around 825,000 observations with missing coordinates. These data have only the names of the nearest cross streets, county, and state information (note they not include street numbers). I need to geocode these observations and recover coordinates so that my final data will look something like this:

   tribble(
     ~street1, ~street2, ~county, ~state, ~latitude, ~longitude
     N BENTON WY, W TEMPLE ST, LOS ANGELES, CA, XX.XXXX, -YY.YYYY,
     11TH PL, BLAINE ST, LOS ANGELES, CA, XX.XXXX, -YY.YYYY,
     W 6TH ST, HOPE ST, LOS ANGELES, CA, XX.XXXX, -YY.YYYY,
     S GRAND AV, W 18TH ST, LOS ANGELES, CA, XX.XXXX, -YY.YYYY,
     BROADWAY, 5TH ST, LOS ANGELES, CA, XX.XXXX, -YY.YYYY,
   )

I have already researched a few possible solutions but haven't found a method that will work.

While the Google Maps API (ggmap package) is very good at identifying coordinates from cross streets as inputs, the cost to geocode this many observations (4.00 USD per 1000 queries according to their website) makes that option infeasible.

I've looked through the documentation of other packages such as RDSTK and tidygeocoder but they don't seem to support API queries using two street names as inputs. The Census Geocoder similarly does not have that option, allowing only single address inputs.

Using the OpenStreetMap API through the osmdata package seemed like a promising option after reading this very detailed StackOverflow answer, but attempting to replicate this code with much bigger bounding boxes has produced runtime errors every time.

See for example the following code using Los Angeles county, following the format of user hugh-allan in the above post:

library(sf)
library(tidyverse)
library(osmdata)

tribble(
      ~point, ~lat, ~lon, 
      1, 32.75004, -118.951721, 
      2, 34.823302, -118.951721, 
      3, 34.823302, -117.646374, 
      4, 32.75004, -117.646374,
    ) %>% 
      st_as_sf(
        coords = c('lon', 'lat'), 
        crs = 4326
      ) %>% 
      {. ->> LA_bounds}
    
    st_bbox(LA_bounds) %>% 
      opq %>% 
      add_osm_feature(key = 'highway') %>% 
      osmdata_sf %>% 
      `[[`('osm_lines') %>% 
      {. ->> LA_streets}

If anyone knows how to get around this error with OpenStreetMaps or otherwise adjust the syntax of another package to accommodate cross streets and counties as inputs, I would greatly appreciate it.

Batch geocoding street intersections in R

0 Answers