11.2 Case study: Landslide susceptibility

This case study is based on a dataset of landslide locations in Southern Ecuador, illustrated in Figure 11.1 and described in detail in Muenchow, Brenning, and Richter (2012).A subset of the dataset used in that paper is provided in the RSAGA package, which can be loaded as follows:

data("landslides", package = "RSAGA")

This should load three objects: a data.frame named landslides, a list named dem, and an sf object named study_area.landslides contains a factor column lslpts where TRUE corresponds to an observed landslide ‘initiation point’, with the coordinates stored in columns x and y.⁶²

There are 175 landslide points and 1360 non-landslide, as shown by summary(landslides).The 1360 non-landslide points were sampled randomly from the study area, with the restriction that they must fall outside a small buffer around the landslide polygons.

To make the number of landslide and non-landslide points balanced, let us sample 175 from the 1360 non-landslide points.⁶³

# select non-landslide points
non_pts = filter(landslides, lslpts == FALSE)
# select landslide points
lsl_pts = filter(landslides, lslpts == TRUE)
# randomly select 175 non-landslide points
set.seed(11042018)
non_pts_sub = sample_n(non_pts, size = nrow(lsl_pts))
# create smaller landslide dataset (lsl)
lsl = bind_rows(non_pts_sub, lsl_pts)

dem is a digital elevation model consisting of two elements:dem$header, a list which represents a raster ‘header’ (see Section 2.3), and dem$data, a matrix with the altitude of each pixel.dem can be converted into a raster object with:

dem = raster(
  dem$data, 
  crs = dem$header$proj4string,
  xmn = dem$header$xllcorner, 
  xmx = dem$header$xllcorner + dem$header$ncols * dem$header$cellsize,
  ymn = dem$header$yllcorner,
  ymx = dem$header$yllcorner + dem$header$nrows * dem$header$cellsize
  )

Figure 11.1: Landslide initiation points (red) and points unaffected by landsliding (blue) in Southern Ecuador.

To model landslide susceptibility, we need some predictors.Terrain attributes are frequently associated with landsliding (Muenchow, Brenning, and Richter 2012), and these can be computed from the digital elevation model (dem) using R-GIS bridges (see Chapter 9).We leave it as an exercise to the reader to compute the following terrain attribute rasters and extract the corresponding values to our landslide/non-landslide data frame (see exercises; we also provide the resulting data frame via the spDataLarge package, see further below):

• slope: slope angle (°).
• cplan: plan curvature (rad m−1) expressing the convergence or divergence of a slope and thus water flow.
• cprof: profile curvature (rad m-1) as a measure of flow acceleration, also known as downslope change in slope angle.
• elev: elevation (m a.s.l.) as the representation of different altitudinal zones of vegetation and precipitation in the study area.
• log10_carea: the decadic logarithm of the catchment area (log10 m2) representing the amount of water flowing towards a location.
  Data containing the landslide points, with the corresponding terrain attributes, is provided in the spDataLarge package, along with the terrain attribute raster stack from which the values were extracted.Hence, if you have not computed the predictors yourself, attach the corresponding data before running the code of the remaining chapter:

# attach landslide points with terrain attributes
data("lsl", package = "spDataLarge")
# attach terrain attribute raster stack
data("ta", package = "spDataLarge")

The first three rows of lsl, rounded to two significant digits, can be found in Table 11.1.