1

I want to do spatial join to an OSM line data and GPS point data both in geographic coordinate system, the result which I want is a road layer(has a unique ID) with the attributes of of road and the point file. I don't want to use buffer: sp::over (I get NA values). I tried using gDistance() it only gives me the distance but does not join the attribute tables. I want to replicate NNJoin of QGIS in R since I can't use batch processing for NNJoin.

I don't think I am satisfied by the answer here https://stackoverflow.com/questions/47675571/r-spatial-join-between-spatialpoints-gps-coordinates-and-spatiallinesdatafra.

longitude <- c(10.86361, 10.96062, 10.93032, 10.93103, 10.93212)        
latitude <- c(44.53355, 44.63234, 44.63470, 44.63634, 44.64559)
longlat <- data.frame(longitude, latitude)  
my_spdf<- SpatialPoints(longlat, proj4string=crs("+init=epsg:4326"))  
shapefile_roads<- spLines(cbind(longitude+.1, latitude), cbind(longitude-.1, rev(latitude)), cbind(longitude-.1, latitude+1), crs=crs("+init=epsg:4326"))
shapefile_roads <- spTransform(shapefile_roads, CRS("+proj=utm +zone=43 ellps=WGS84"))
my_spdf <- spTransform(my_spdf, CRS("+proj=utm +zone=43 ellps=WGS84"))
result <- over(shapefile_roads, my_spdf)
9
  • 1
    Time to try the sf package. It has a function called st_join which does exactly hat you want when you specify st_nearest_feature as the join.
    – TimSalabim
    Commented Jan 11, 2019 at 15:01
  • I get this error: Error in st_join(shapefile_roads, my_spdf, join = st_nearest_feature) : inherits(x, "sf") && inherits(y, "sf") is not TRUE
    – ADIR
    Commented Jan 11, 2019 at 16:12
  • Even when I convert them into projected Coordinate system.
    – ADIR
    Commented Jan 11, 2019 at 16:21
  • 1
    Why exactly are you "unsatisfied" with the linked answer? It is correct and efficient. You are getting the error using the sf st_nearest_feature function because you are not coercing your data into a sf class object. You really need to learn the basics of R spatial classes before you start rejecting solutions. Commented Jan 11, 2019 at 16:29
  • Actually i m new to R. I did not understand the error.
    – ADIR
    Commented Jan 11, 2019 at 16:31

1 Answer 1

3

Here's a solution using library(sf):

library(sf)
library(rgdal)
library(sp)
library(raster)

longitude <- c(10.86361, 10.96062, 10.93032, 10.93103, 10.93212)        
latitude <- c(44.53355, 44.63234, 44.63470, 44.63634, 44.64559)
longlat <- data.frame(longitude, latitude)  
my_spdf<- SpatialPoints(longlat, proj4string=CRS("+init=epsg:4326"))  
shapefile_roads<- spLines(cbind(longitude+.1, latitude), cbind(longitude-.1, rev(latitude)), cbind(longitude-.1, latitude+1), crs=CRS("+init=epsg:4326"))
shapefile_roads <- spTransform(shapefile_roads, CRS("+proj=utm +zone=43 ellps=WGS84"))
my_spdf <- spTransform(my_spdf, CRS("+proj=utm +zone=43 ellps=WGS84"))

my_sfdf = st_sf(id_pt = 1:length(my_spdf), geom = st_as_sfc(my_spdf))
shapefile_roads_sf = st_sf(id_ln = 1:length(shapefile_roads), 
                           geom = st_as_sfc(shapefile_roads))

# get nearest point for each line
shapefile_roads_sf_w_pts = st_join(shapefile_roads_sf, 
                                   my_sfdf, 
                                   join = st_nearest_feature)

shapefile_roads_sf_w_pts
Simple feature collection with 3 features and 2 fields
geometry type:  LINESTRING
dimension:      XY
bbox:           xmin: -4358172 ymin: 7326658 xmax: -4226019 ymax: 7428067
epsg (SRID):    32643
proj4string:    +proj=utm +zone=43 +ellps=WGS84 +units=m +no_defs
  id_ln id_pt                           geom
1     1     2 LINESTRING (-4350438 732675...
2     2     1 LINESTRING (-4348820 735297...
3     3     5 LINESTRING (-4243179 742480...

Please, when posting a reproducible example, make sure that it does indeed run and also include the packages that need to be loaded in order to execute the code. We cannot know where functions come from if the execution fails! Also, if you say you want IDs, please make sure that your example includes them!

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.