Subset Spatial Points to Extract those inside of a Polygon (Countries Borders)

1.7k Views Asked by At

I have a data.frame with lat long coordinates:

df<-data.frame(
     lat=c(40, 30, 40.864),
     lon=c(0, 20, 1.274)
    )

And the border of a country (Spain),

library(raster)
border <- getData("GADM",country="Spain",level=0)

I would like to select for df only the points inside the border.

How can I do this?

Note: in my reproducible example, df, the first entry point is inside, the second is clearly outside and the third is outside but close to the coast.

1

There are 1 best solutions below

0
On BEST ANSWER

While I believe the post that I referred to answers this question, I think some clarification is needed; hence, posting an answer.

To find the intersection of two geographic feature, we need to have the same projection. library helps us to do so. Make sure to put longitude and latitude in the right place:

sp::SpatialPoints(c(my_point$long,my_point$lat),proj4string=CRS(proj4string(my_raster)))

Using library we can check if there is an intersection between two spatial dataset/feature:

rgeos::gContains(my_raster,my_projected_point)

So, here is how it works for the OP's example:

library(sp)        #for projection
library(raster)    #for getting the border data
library(rgeos)     #for finding intersection
library(tidyverse) #for illustration only

#data
border <- getData("GADM",country="Spain",level=0)
df <- data.frame(
 lat=c(40, 30, 40.864),
 lon=c(0, 20, 1.274)
                )

#this is the part that actually check if a point is inside the border
#adapted from https://stackoverflow.com/questions/21971447
sapply(1:3,function(i)
  list(id=i,
       intersect= gContains(border,SpatialPoints(df[i,2:1],proj4string=CRS(proj4string(border))))))

#           [,1] [,2]  [,3] 
# id        1    2     3    
# intersect TRUE FALSE FALSE

#a map for better understanding
ggplot()+
  geom_polygon(data=border, aes(x=long, y=lat, group=group), 
                       fill=NA, color="grey50", size=0.25) +
  geom_point(data=df,aes(x=lon,y=lat), color="red", size=1)

enter image description here