Geospatial analysis: Clipping multiple shapefiles from a single csv file

Question

I recently collected data dat, with gps cordinates. Unit of data collection was ward. On ploting I realised that some points are outside the expected wards. There are a number of wards, kindly assist with an R funtion to clip each point based on the ward shapefile

I tried to compare datasets to obtain the points within. the output came as a geometry and could merge on the original csv to obtain the result

Answer 1

Hi Omena Fulu and welcome to Stack Overflow!! You can accomplish the equivalent of an ArcGIS "clip" using the st_intersection() function in the "sf" package. I have created a reproducible example for you using example datasets available in the "spData" package (see the end). In your case, you will want to first read in the csv using either the read_csv() function or the read.table() function, and then convert the csv into a simple feature dataframe using the st_as_sf() function in the "sf" package. In my example, the "cycle_hires" dataset is a point dataset that would be analogous to your GPS point dataset in the csv. The "lnd" dataset coincidentally also consists of administrative wards in the London Metropolitan area, and should be almost exactly analogous to what you're after. Below is the example:

##Loading Necessary Packages##
library(spData)
library(sf)
library(tmap)
library(tmaptools)


#Getting Example Data##
data("lnd") #Polygon Data similar to wards
data("cycle_hire") #Point Data similar gps points

##Take a look at the raw data##
#Setting tmap to interactive viewing#
tmap_mode("view")

#drawing the map#
tm_shape(lnd)+
  tm_polygons(col="NAME", alpha=0.5)+
tm_shape(cycle_hire)+
  tm_dots(col="red", size = .01)+
tm_basemap(server=providers$Esri.WorldImagery)


##Geoprocessing the data##
Westminster<-lnd[lnd$NAME=="Westminster",] #Get just the Westminster geographical area
westminster_cycles<-st_intersection(cycle_hire, Westminster) #use the st_intersection() function to clip cycle hires to just the Westminster geographical area

##View the results##
#drawing the new map#
tm_shape(lnd)+
  tm_polygons(col="NAME", alpha=0.5)+
tm_shape(westminster_cycles)+
  tm_dots(col="red", size = .01)+
  tm_basemap(server=providers$Esri.WorldImagery)