Create a sf from bounding box coordinates (and compare to another sf) in R

Question

I want to create a sf from bounding box coordinates and compare it to another sf to see if it intersects. When I use the code below it doesn't seem to work:

1. I can't plot the sf as I get the error Error: Invalid index: field name 'x_start' not found and

2. st_intersects gives me [,1] [,2] [,3] [,4] [1,] FALSE FALSE FALSE FALSE even though I know they should intersect.

Here is my code:

  # Read in UK shapefile
  uk.shapefile <- file.path(uk.shapefile_path) %>%
    st_read() %>%
    st_transform(crs = 27700)

  # Define the bounding box
  bbox_coords <- c(-7.84000, 49.68000, -6.52750, 51.10625)

  # Create a polygon from the bounding box
   bbox_polygon <- st_as_sf(st_sfc(st_polygon(list(matrix(c(
     bbox_coords[1], bbox_coords[2],  
     bbox_coords[3], bbox_coords[2],
     bbox_coords[3], bbox_coords[4],  
     bbox_coords[1], bbox_coords[4],
     bbox_coords[1], bbox_coords[2]   
   ), ncol = 2, byrow = TRUE))), crs = 27700))
   
   # Plot the bounding box - **returns error**
   ggplot() +
     geom_sf(data = bbox_polygon, fill = "blue", alpha = 0.3)
   
   # Find intersection - **returns FALSE**
   
   st_intersects(bbox_polygon, uk.shapefile, sparse = FALSE)

Answer 1

The easiest way from two corner points in a vector to a polygon is to name the elements and then make a bbox and make an sfc:

> bbox_coords <- c(-7.84000, 49.68000, -6.52750, 51.10625)
> names(bbox_coords) = c("xmin","ymin","xmax","ymax")

> bbp = st_as_sfc(st_bbox(bbox_coords), crs=4326)
> bbp
Geometry set for 1 feature 
Geometry type: POLYGON
Dimension:     XY
Bounding box:  xmin: -7.84 ymin: 49.68 xmax: -6.5275 ymax: 51.10625
CRS:           NA
POLYGON ((-7.84 49.68, -6.5275 49.68, -6.5275 5...

Saves all that messing of making a polygon/list/matrix from vector elements.

Answer 2

Since you are entering bounding box as defined by Longitude and Latitude (which is perfectly fine approach) consider converting it to a polygon within context of EPSG:4326 CRS (and not the British National Grid; units matter, and to paraphrase the well known meme about dollars and murders - 6 degrees is a lot, 6 meters not so much).

It will plot just fine, and in order to test intersection with your uk.shapefile object (to which we don't have access) you can easily transform it to whatever CRS the UK object is in - note how I am not hardcoding a CRS specification, but using the value returned from st_crs(uk.shapefile) - that way you are certain that your CRS definition will be aligned.

library(sf)
library(dplyr)
library(ggplot2)

# Define the bounding box
bbox_coords <- c(-7.84000, 49.68000, -6.52750, 51.10625)

# Create a polygon from the bounding box
bbox_polygon <- st_as_sf(st_sfc(st_polygon(list(matrix(c(
   bbox_coords[1], bbox_coords[2],  
   bbox_coords[3], bbox_coords[2],
   bbox_coords[3], bbox_coords[4],  
   bbox_coords[1], bbox_coords[4],
   bbox_coords[1], bbox_coords[2]   
), ncol = 2, byrow = TRUE))), crs = 4326))

# Plot the bounding box
ggplot() +
   geom_sf(data = bbox_polygon, fill = "blue", alpha = 0.3)

# Find intersection - **returns FALSE**

st_intersects(st_transform(bbox_polygon, st_crs(uk.shapefile)),
              uk.shapefile, sparse = FALSE)