Accessing Rangeland Analysis Platform (RAP) Gridded Products with R

library(rapr)

library(soilDB)  # get soil mapunit polygons
library(sf)      # spatial data handling
#> Linking to GEOS 3.10.2, GDAL 3.4.1, PROJ 8.2.1; sf_use_s2() is TRUE
library(mapview, warn.conflicts = FALSE) # interactive maps
library(gifski)  # create GIF from results

First we buffer 100m around a longitude/latitude coordinate (WGS84 decimal degrees) using the {sf} package.

You could change the coordinates to where you live or your favorite range spot!

You can change the buffering distance to include larger areas or remove it altogether to make sure you get just a single polygon.

x <- st_buffer(st_as_sf(
  data.frame(x = -120.26021, y = 37.99765),
  coords = c("x", "y"),
  crs = 4326
), 100)

We use {soilDB} to retrieve SSURGO mapunit polygons for our area of interest.

p <- SDA_spatialQuery(x, what = "mupolygon")

In this example the target polygon has a significant Urban land component associated with the highway and adjacent development. How much of the polygon is reflecting inter-annual variation in plant cover?

There are several housing developments and some areas have been cleared for grazing. The natural condition is closed canopy mixed oak/conifer woodland.

mapview(p, map.types = c("Esri.WorldImagery", "OpenStreetMap"))

Use {rapr} to download the Rangeland Analysis Platform “vegetation-biomass” product for 1986 to 2021 using the mapunit polygon p to define the target extent.

rap <- rapr::get_rap(p, product = "vegetation-biomass", 1986:2021, progress = FALSE)

Now we will select just the "annual forb and grass biomass" layers, iterate over them and plot symbolizing with a common range [0,5000] pounds per acre.

We write this iteration into a function called makeplot() and use {gifski} to render an animated GIF file from the R plot graphics output in each year for a total of 36 layers.

Pass the function makeplot() as the first argument to gifski::save_gif() and then specify the output file name with gif_file argument.

makeplot <- function() {
  lapply(grep("annual_forb_and_grass_biomass", names(rap)), function(i) {
    terra::plot(rap[[i]], 
                main = gsub(".*_(\\d+)_v3$", "\\1", names(rap)[i]), 
                range = c(0, 5000), 
                cex.main = TRUE)
    terra::plot(terra::vect(p), add = TRUE)
  })
}

gifski::save_gif(makeplot(), gif_file = "annual_forb_and_grass_biomass.gif", delay = 0.5)
#> [1] "/home/runner/work/rapr/rapr/vignettes/annual_forb_and_grass_biomass.gif"