It seemed that I couldn’t find a good visual representation of what some of the operations in the sf package did, so I went exploring what was available. Many of these tasks might be done in typical GIS but I was eager to see how they might translate to doing geographic analysis in R.
I made up some fake data that approximates what would be encountered in doing work with census data. We have a city boundary and a number of geographical units. The 2 layers should be aligned, so that the lines and vertices of the polygons overlap.
First, let’s see what I made:
ggplot(geo_units) + geom_sf(size = 2) + geom_sf(data = city_bnd, colour = "red", fill = NA, size = 1)
From the vingette we find:
“The commands st_intersects, st_disjoint, st_touches, st_crosses, st_within, st_contains, st_overlaps, st_equals, st_covers, st_covered_by, st_equals_exact and st_is_within_distance return a sparse matrix with matching (TRUE) indexes, or a full logical matrix”
Thanks to Kyle Walker (here) we can use this matrix to subset those units that fall within certain locations relative to the city boundary.
x1 <- st_intersects(geo_units, city_bnd)
x2 <- map_lgl(x1, function(x) {
if (length(x) == 1) {
return(TRUE)
} else {
return(FALSE)
}
})
ex1 <- geo_units[x2,]
ggplot() +
geom_sf(data = geo_units, lwd = .25, fill = NA) +
geom_sf(data = ex1, fill = "slategray2") +
geom_sf(data = city_bnd, fill = NA, color = "red", size = 1) +
ggtitle("sf geometrical operations", subtitle = "st_intersects")
We can see that intersecting includes when a ploygon falls on the line of the city bounday too.
st_disjoint appears to be the opposite of st_intersects, good to know.
Touches appears to mean on the outside of the city boundary only.
Perhaps this shows that none of the unit polygons completely crosses the city boundary? Not sure about this one.
These polygons are within the city boundary, though the upper left one should meet this criteria too. Maybe is has somthing to do with the city boundary hitting mid-line on the lower side.
Perhaps none of the unit polygons contain the city?
Jagged lines mean that the city boundary is overlapping a unit polygon. The upper left probably should not be selected.
In this case ‘covered_by’ equals ‘within’.
These are some of the operations that might play into a spatial analysis. Being able to do these in R can allow for more GIS functionality without having to go out to another application. Nice!
(I am by no means an expert on this, so if you see something that’s off please do reach out - data@sogletree.com. Thanks!)