You really are forced to define a spatial extent to bound the data. You can draw a bounding box on a plot screen using raster::drawExtent()
. One issue with this data is that it comes as MULTIPART geometry, which is one attribute row representing multiple spatial features. This is why a simple indexing approach or a geometry intersection operator will not work using an sp class object.
I would recommend exploding the polygon geometry into SINGLEPART.
Add libraries, data and subset
library(raster)
library(sp)
library(rgdal)
library(spatialEco)
library(rgeos)
incheon <- rgeos::readOGR(getwd(), "gadm36_KOR_1")
incheon <-incheon[incheon$NAME_1 == "Incheon",]
Here you can explode the MULTIPART geometry. If you look at dim(incheon)
you will see 1 feature with 10 attributes and 111 features with 10 attributes after applying explode
.
incheon <- spatialEco::explode(incheon, sp = TRUE)
Here is where we create an extent SpatialPolgons
object.
e <- as(raster::extent(125.6698, 126.8242, 36.76623, 37.90437),
"SpatialPolygons")
And, finally intersect the data to create the spatial subset.
incheon <- rgeos::gIntersection(incheon, e, byid=TRUE)
plot(incheon)