14  Buffers vs. Isochones

We will use the stplnar package to create buffers in meters, and the openrouteservice package to create isochrones.

14.1 Buffer

Represent a buffer of 500 m and 2000 m from IST1.

IST = st_sfc(st_point(c(-9.1397404, 38.7370168)), crs = 4326) # create a point
IST$coordinates = st_coordinates(IST)

# BUFFERist500 = st_buffer(IST, dist = 500) # non  projected - results may be weird
BUFFERist500 = geo_buffer(IST[1], dist = 500) # from stplnar, to make sure it is in meters.
BUFFERist2000 = geo_buffer(IST[1], dist = 2000)

mapview(BUFFERist500) + mapview(BUFFERist2000, alpha.regions = 0.5)

14.2 Isochrone

Isochrone from 1 point - distance

We use again the openrouteservice r package.

ISOCist = ors_isochrones(
  IST$coordinates,
  profile = "foot-walking",
  range_type = "distance", # or time
  range = c(500, 1000, 2000),
  output = "sf"
)

ISOCist = arrange(ISOCist, -value) |>  # to make the larger polygons on top of the table so the are displayed behind.
          select(-center) # unnecessary variable
mapview(ISOCist, zcol = "value", alpha.regions = 0.5)

As you can see, the distance buffer of 500m is larger than the isochrone of 500m. Actually we can measure their area of reach.

ISOCist$area = st_area(ISOCist)
BUFFERist500$area = st_area(BUFFERist500)
BUFFERist2000$area = st_area(BUFFERist2000)

ratio1 = BUFFERist500$area / ISOCist$area[ISOCist$value == 500] # 1.71
ratio2 = BUFFERist2000$area / ISOCist$area[ISOCist$value == 2000] # 1.22

The euclidean buffer of 500m is 1.7 times larger than its isochrone, and the buffer of 2000m is 1.23 times larger than its isochrone.

Isochrone from more than 1 point - time

For this purpose we will use the high schools dataset.

# import schools
SCHOOLS = st_read("../geo/SCHOOLS_basicsec.gpkg", quiet = TRUE)

SCHOOLS$coordinates = st_coordinates(SCHOOLS) # create coordinate variable

SCHOOLShigh = SCHOOLS |>
  filter(Nivel == "Secundario")  # filter the high schools
SCHOOLShigh = SCHOOLShigh |> mutate(id = 1:nrow(SCHOOLShigh)) # provide and id

# list of XY coordinates for ORS
coord_schools = data.frame(lon = SCHOOLShigh$coordinates[, 1],
                           lat = SCHOOLShigh$coordinates[, 2])
                           # id = SCHOOLShigh$id)

And proceed with the time isochrones, for a range of 20 min, with 5 min intervals.

ISOCschools = ors_isochrones(
  coord_schools,
  profile = "foot-walking",
  range_type = "time", # or distance
  range = 20*60, # 20 minutes in seconds
  interval = 5*60, # to have intervals of 5 minutes
  attributes = "area", #you can directly get area, population, and so on.
  output = "sf"
)

ISOCschools = arrange(ISOCschools, -value) |> # larger polygons on top of the table so the are displayed behind.
  select(-center)

And now merge this information with the schools’ names.

ISOCschools = ISOCschools |>
  mutate(id = group_index + 1, # because group_index starts at 0
         value = value/60) |> 
    left_join(SCHOOLShigh |>
                st_drop_geometry() |>
                select(id, INF_NOME, Alunos))
mapview(ISOCschools, zcol = "value", alpha.regions = 0.5)

And estimate areas for the 20min isochrones

summary(ISOCschools$area[ISOCschools$value==20])/1000000 # in km²
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  3.295   4.303   5.363   5.098   5.804   6.511

  1. Here I selected only the first variable because now we also have the coordinates information (unecessary for this procedure)↩︎