20  Roads

Today is going to be an showing how to get some roads data.

20.1 Load libraries

To enable caching of data, set `options(tigris_use_cache = TRUE)`
in your R script or .Rprofile.
Linking to GEOS 3.12.1, GDAL 3.8.4, PROJ 9.3.1; sf_use_s2() is TRUE

Attaching package: 'dplyr'
The following objects are masked from 'package:stats':

    filter, lag
The following objects are masked from 'package:base':

    intersect, setdiff, setequal, union
knitr::opts_chunk$set(message = NA)

20.1.1 Tigris function for scraping road data.

The aptly named roads() function will pull road linestrings.

Riv_roads <- roads(state = 'CA', county = 'Riverside', year = 2023, progress_bar = FALSE) |> 
  st_transform(crs = 4326)

head(Riv_roads)
Simple feature collection with 6 features and 4 fields
Geometry type: LINESTRING
Dimension:     XY
Bounding box:  xmin: -117.5613 ymin: 33.49898 xmax: -116.3334 ymax: 33.87767
Geodetic CRS:  WGS 84
       LINEARID               FULLNAME RTTYP MTFCC
1  110458254891       Armstrong Dr Exd     M S1400
2  110458266265 Gabino Canyon Rd N Spr     M S1400
3 1102175699783           Porphyry Spr     M S1400
4  110458248901        Pub Utility Acc     M S1400
5 1103677442858        White Horse Pvt     M S1740
6 1103678156677            Bus Park Dr     M S1400
                        geometry
1 LINESTRING (-117.5613 33.83...
2 LINESTRING (-117.51 33.5231...
3 LINESTRING (-117.5396 33.87...
4 LINESTRING (-117.1296 33.50...
5 LINESTRING (-116.3343 33.75...
6 LINESTRING (-117.1577 33.49...

There are 58,448 road segments from 2023. That’s a lot of lines for a map and takes too long to load. Feel free to make a map with it if you want to hear your computer hum.

For the purpose of the road widening project, we’re mostly focusing on big roads getting widened. The documentation for roads() states that the RTTYP category provides six types.

I’m going to filter just the Interstate (I) and State (S) roads for this version.

Riv_major_roads <- Riv_roads |> 
  filter(RTTYP %in% c('I', 'S'))

There are only 81 major roads segments that meet that criteria.

  leaflet() |> 
  addTiles() |> 
  addPolylines(data = Riv_major_roads,
              color = 'black',
              weight = 3,
              dashArray = '2,6', #new feature - makes lines dashed - line lenght, space length
              label = ~FULLNAME)
Figure 20.1: Riverside roads on a leaflet map

Now let’s add Riverside County warehouses for a quick look at the proximity of those facilities to current major roadways. We did this in lesson 10.3.2

WH.url <- 'https://raw.githubusercontent.com/RadicalResearchLLC/WarehouseMap/main/WarehouseCITY/geoJSON/comboFinal.geojson'
warehouses <- st_read(WH.url) |>  
  st_transform(crs = 4326) |> 
  filter(county == 'Riverside County')
Reading layer `comboFinal' from data source 
  `https://raw.githubusercontent.com/RadicalResearchLLC/WarehouseMap/main/WarehouseCITY/geoJSON/comboFinal.geojson' 
  using driver `GeoJSON'
Simple feature collection with 9084 features and 8 fields
Geometry type: MULTIPOLYGON
Dimension:     XY
Bounding box:  xmin: -118.8037 ymin: 33.43325 xmax: -114.4085 ymax: 35.55527
Geodetic CRS:  WGS 84

Figure 20.2 shows the major roads and warehouse land-use.

palWH <- colorFactor(palette = c('orange', 'gold', 'brown'),
                     domain = warehouses$category)
palRoadType <- colorFactor(palette = c('darkblue', 'purple'),
                           domain = Riv_major_roads$RTTYP)

leaflet() |> 
  addProviderTiles(provider = providers$CartoDB.Positron) |> ## We modified this line to change to satellite view  
  addPolygons(data = warehouses,
              color = ~palWH(category),
              weight = 1,
              fillOpacity = 0.7) |> 
  addLegend(data = warehouses,
            pal = palWH,
            values = ~category) |> 
  addPolylines(data = Riv_major_roads,
            color = ~palRoadType(RTTYP),
            weight = 3,
            dashArray = '2,6',
            label = ~FULLNAME) |> 
  setView(zoom = 10, lat = 33.8, lng = -117.15)
Figure 20.2: Riverside roads and warehouses on a leaflet map

This is a good start, but there’s some missing pieces that we can show a few examples for road widening.

Let’s filter() down a data table for a couple of key roads that are approved or under review for road widening.

roads2widen <- Riv_roads |> 
  filter(FULLNAME %in% c('Cajalco Rd', 'Cajalco Expy', 'Ramona Expy', 'Gilman Springs Rd'))

This selects just a few individual county level roads. Let’s add those to a new map.

Figure 20.3 shows the result.

palWH <- colorFactor(palette = c('orange', 'gold', 'brown'),
                     domain = warehouses$category)
palRoadType <- colorFactor(palette = c('darkblue', 'purple'),
                           domain = Riv_major_roads$RTTYP)

leaflet() |> 
  addProviderTiles(provider = providers$CartoDB.Positron) |> ## We modified this line to change to satellite view  
  addPolygons(data = warehouses,
              color = ~palWH(category),
              weight = 1,
              fillOpacity = 0.7) |> 
  addLegend(data = warehouses,
            pal = palWH,
            values = ~category,
            position = 'bottomleft') |> 
  addPolylines(data = Riv_major_roads,
            color = ~palRoadType(RTTYP),
            weight = 3,
            dashArray = '2,6',
            label = ~FULLNAME) |> 
  setView(zoom = 10, lat = 33.8, lng = -117.15) |> 
  addPolylines(data = roads2widen,
               color = 'red',
               weight = 4,
              dashArray = '2,6',
              label = ~FULLNAME) |> 
  addLegend(title = 'Roads',
            colors = c('darkblue', 'purple', 'red'),
            labels = c('Interstate', 'State Route', 'Road widening'),
            position = 'bottomleft')
Figure 20.3: Riverside roads and widening roads with warehouses on a leaflet map