left_join based on closest LAT_LON in R

Question

I am trying to find the ID of the closest LAT_LON in a data.frame with reference to my original data.frame. I have already figured this out by merging both data.frames on a unique identifier and the calculating the distance based on the distHaverSine function from geosphere. Now, I want to take step further and join the data.frames without the unique identifier and find ID the nearest LAT-LON. I have used the following code after merging:

v3 <-v2 %>% mutate(CTD = distHaversine(cbind(LON.x, LAT.x), cbind(LON.y, LAT.y)))

DATA:

loc <- data.frame(station = c('Baker Street','Bank'),
     lat = c(51.522236,51.5134047),
     lng = c(-0.157080, -0.08905843),
               postcode = c('NW1','EC3V'))

stop <- data.frame(station = c('Angel','Barbican','Barons Court','Bayswater'),
                lat = c(51.53253,51.520865,51.490281,51.51224),
                lng = c(-0.10579,-0.097758,-0.214340,-0.187569),
                postcode = c('EC1V','EC1A', 'W14', 'W2'))

As a final result I would like something like this:

df <- data.frame(loc = c('Baker Street','Bank','Baker Street','Bank','Baker Street','Bank','Baker 
        Street','Bank'), 
              stop = c('Angel','Barbican','Barons Court','Bayswater','Angel','Barbican','Barons Court','Bayswater'), 
              dist = c('x','x','x','x','x','x','x','x'), 
              lat = c(51.53253,51.520865,51.490281,51.51224,51.53253,51.520865,51.490281,51.51224), 
              lng = c(-0.10579,-0.097758,-0.214340,-0.187569,-0.10579,-0.097758,-0.214340,-0.187569),
              postcode = c('EC1V','EC1A', 'W14', 'W2','EC1V','EC1A', 'W14', 'W2')
              )

Any help is appreciated. Thanks.

Answer 1

All of the joining, distance calculations, and plotting can be done with available R packages.

library(tidyverse)
library(sf)
#> Linking to GEOS 3.6.2, GDAL 2.2.3, PROJ 4.9.3
library(nngeo)
library(mapview)

## Original data
loc <- data.frame(station = c('Baker Street','Bank'),
                  lat = c(51.522236,51.5134047),
                  lng = c(-0.157080, -0.08905843),
                  postcode = c('NW1','EC3V'))

stop <- data.frame(station = c('Angel','Barbican','Barons Court','Bayswater'),
                   lat = c(51.53253,51.520865,51.490281,51.51224),
                   lng = c(-0.10579,-0.097758,-0.214340,-0.187569),
                   postcode = c('EC1V','EC1A', 'W14', 'W2'))

df <- data.frame(loc = c('Baker Street','Bank','Baker Street','Bank','Baker Street','Bank','Baker 
        Street','Bank'), 
                 stop = c('Angel','Barbican','Barons Court','Bayswater','Angel','Barbican','Barons Court','Bayswater'), 
                 dist = c('x','x','x','x','x','x','x','x'), 
                 lat = c(51.53253,51.520865,51.490281,51.51224,51.53253,51.520865,51.490281,51.51224), 
                 lng = c(-0.10579,-0.097758,-0.214340,-0.187569,-0.10579,-0.097758,-0.214340,-0.187569),
                 postcode = c('EC1V','EC1A', 'W14', 'W2','EC1V','EC1A', 'W14', 'W2')
)



## Create sf objects from lat/lon points
loc_sf <- loc %>% st_as_sf(coords = c('lng', 'lat'), remove = T) %>%
  st_set_crs(4326) 

stop_sf <- stop %>% st_as_sf(coords = c('lng', 'lat'), remove = T) %>%
  st_set_crs(4326) 


# Use st_nearest_feature to cbind loc to stop by nearest points
joined_sf <- stop_sf %>% 
  cbind(
    loc_sf[st_nearest_feature(stop_sf, loc_sf),])


## mutate to add column showing distance between geometries
joined_sf %>%
  mutate(dist = st_distance(geometry, geometry.1, by_element = T))
#> Simple feature collection with 4 features and 5 fields
#> Active geometry column: geometry
#> geometry type:  POINT
#> dimension:      XY
#> bbox:           xmin: -0.21434 ymin: 51.49028 xmax: -0.097758 ymax: 51.53253
#> epsg (SRID):    4326
#> proj4string:    +proj=longlat +datum=WGS84 +no_defs
#>        station postcode    station.1 postcode.1                   geometry
#> 1        Angel     EC1V         Bank       EC3V  POINT (-0.10579 51.53253)
#> 2     Barbican     EC1A         Bank       EC3V POINT (-0.097758 51.52087)
#> 3 Barons Court      W14 Baker Street        NW1  POINT (-0.21434 51.49028)
#> 4    Bayswater       W2 Baker Street        NW1 POINT (-0.187569 51.51224)
#>                    geometry.1         dist
#> 1 POINT (-0.08905843 51.5134) 2424.102 [m]
#> 2 POINT (-0.08905843 51.5134) 1026.449 [m]
#> 3   POINT (-0.15708 51.52224) 5333.417 [m]
#> 4   POINT (-0.15708 51.52224) 2390.791 [m]



## Use nngeo and mapview to plot lines on a map
# NOT run for reprex, output image attached 
#connected <- st_connect(stop_sf, loc_sf)
# mapview(connected) + 
#   mapview(loc_sf, color = 'red') +
#   mapview(stop_sf, color = 'black')

^{Created on 2020-01-21 by the reprex package (v0.3.0)}