Landslide Exercise

Pre-requisites:

Setup
Basics
Data Hygiene

Objectives:

Prep data for a spatial analysis
Perform a spatial analysis utilizing both vector and raster data
Visualize results from the spatial analysis

Exercise: Landslide Susceptibility Model

Landslide Susceptibility Model

Overview

This case study has been chosen to demonstrate the concepts being highlighted by the bootcamp. We have chosen a landslide susceptibility analysis due to the diversity of inputs and a general curiousity into the methodology. The study area has been chosen in response to the Oso Mudslide[^1] which took place in Washington State on March 22nd, 2014. The methodology for this analysis is based off of work done by the California Geological Survey in Conjunction with the USGS[^2].

1. Load the Geology Layer

Open QGIS and first set the project projection to ESPG:2927. This will ensure we have the correct projection after importing layers. It's recommended to download and unzip the wa_geology layer due to its larger file size.

Load wa_geology.shp through the iRods plugin:
/iplant/home/shared/aegis/Spatial-bootcamp/spatial-analysis/landslide-exercise/wa_geology.shp
Or download, unzip, and Add Vector Layer :
wa_geology.zip

QGIS: Add geology vector layer into map view

2. Inspect Washington Geology Attributes

If you open the Washington geology attributes you'll notice there's a column 'rock_type' with values as string. Our landslide susceptibility model requires that we have a numeric value representing 'rock_type'. Notice the range of 'rock_type': soft, medium, hard. Wow that's an easy one! It's easy as 1, 2, 3 (literally). You're welcome for the heavy lifting. See 'geo_unit' for actual geologic unit types.

We are going to add a new field to the wa_geology shapefile to represent the strength of the geologic unit. This new field data type will be number, not a string. This will allow us to perform arithmetic opertions later on.

Think about it: the value '1' represented as a string does not equal the numeric value of 1. In this case, an integer.

3. Add Rock strength index

Open field calculator in Attribute Table and enter Edit Mode:

Reminder: Making changes to your data in Edit Mode are permanent. It's good practice to store source data, or the original files in a separate directory so in the event you make a catastrophic mistake you don't have to ask your colleague, boss, or whomever for another copy of the data.
1. Right-click wa_geology (Layers List) > Open Attribute Table > Toggle editing mode (top left)
2. While in edit mode,click Open field calculator (top right)
  
  Whil in Edit Mode the Attribute table icons will be enabled: Save, Delete, New Column.
  
  Notice the shapefile has turned red on the map canvas. You can also edit the shapefile nodes (e.g. coordinates = point = node), changing the geometry of the file, since all we are really dealing with are numbers. We are simply visualizing these numbers (lat/long) and their geometry.
Add and calculate new field:
Be sure to match your Field Calculator window to the one below for creating a new field strength with a range of values of 0-3; if rock_type is not soft, medium, or hard then assign 0 to the value. Be sure to match it exactly as it's displayed

Note: double-quotes are evaluated as a column whereas single-quotes are evaluated as column values. Save your changes once finished.

Configure Field Calculator inputs as follows:
Create a new field: (checked)
Output field name: strength
Output field type: Whole number (integer)
Output field width: 1 (we only need single-digit integers)
Expression:

Copy and paste will work for entering the code into the Expression field.

case when "rock_type" = 'soft' then 3 when "rock_type" = 'medium' then 2 when "rock_type" = 'hard' then 1 else 0 end

This field calculation will take a couple of seconds. We're evaluating 57,657 rows.
Confirm field calculation of 'strength':

Be sure to Save Edits then Disable Editing before continuing. Your geology layer will no longer be red (Edit Mode).

4. Visualize the geology by strength

Let's visualize the new field we have just calculated.

Open wa_geology Style properties: Right-click wa_geology (Layer list) > Properties > Style
Change style type from Single Symbol to
Categorized .
Categorize by strength and click Classify to add classes.
Match your colors with the example below or get fancy with your own styling.

There's a category with no Value. You can leave this one there or highlight and delete it.

Be sure to APPLY to make changes then OK to close the Properties window

Drum roll...

Reminder: Red = soft geology (land) = higher risk of landslides

5. Rasterize the rock strength layer

To be able to perform the landslide susceptibility model we must intersect rock strength with slope (derrived from a DEM). So the best (highly suggested) method is converting our geology shapefile into a raster. A simple process thanks to GDAL_rasterize.

Because GDAL is more powerful in the commandline, there are certain parameters that cannot be set within QGIS.

Open gdal_rasterize: Menu Bar > Raster > Conversion > Rasterize (Vector to raster)
Configure inputs as follows:
Input file: wa_geology
Attribute filed: strength (strength value assigned to pixels)
Output file for rasterized vectors (raster): wa_geo_coded
"The output file doesn't exist. You must set up the output size or resolution to create it:" OK
Raster size in pixels: 3280 x 3280 (1km resolution)
Enable editing the gdal_rasterize code by clicking the Edit icon near the block of code. This allows to add custom parameters.

Add the line of code: -at -a_nodata 0

The -at flag enables the ALL_TOUCHED rasterization option so that all pixels touched by lines or polygons will be updated, not just those on the line render path, or whose center point is within the polygon. Defaults to disabled for normal rendering rules.[^3]

-a_nodata 0 The above snippet specifies that 0 should be recognized as a missing or null value, which prevents pixels with 0 value from being visualized or used in analysis in the output raster.[^3]

Drum roll...

Notice the light blue (water) above. We declared 0 as no_data so water is no longer included in our raster dataset.

6. Style the new Rock Strength Raster

It’s good to visually check your results, if any large issues may occur. Sometimes you must look at your data values to locate errors. Let’s style our new geology raster.

Open wa_geo_coded style properties.
Render type: Singleband pseudocolor
Load actual values:
- Load min/max values: Min/max
- Extent: Actual (slower)
- Click LOAD to load actual values
This loads actual values into the classification. QGIS Approximates the Minimum and Maximum values by default.
Change Mode to Equal Interval with 3 classes.

Remember, we only have three rock strength classes: soft, medium, hard. The color maps are not stored into the raster, so the raster will load with the default color map (black to white) each time its opened in a new project.
Click Classify to add categories.
Click APPLY to confirm changed and OK to exit properties.

7. Prepare slope for analysis:

Now it's time to perform an analysis. We're going to reclassify our slope into something useful.

Import wa_slope.tif

This slope raster was created in the Spatial Functions. If you still have it, great. Otherwise the slope raster is available here:

iRods access: /iplant/home/shared/aegis/Spatial-bootcamp/spatial-analysis/landslide-exercise/wa_slope.tif

Or download file: wa_slope.tif

In order to reclass our slope layer into something useful, we will use the grass r.reclass tool

The r.reclass tool reads text files which define the classification rules (rule files). We must prepare a text file with the following contents using a text editor not Microsoft Word

Create slope-rule.txt:

Copy and paste the code below into a new file and save it as slope-rule.txt

0 thru 2.99 = 1 3 thru 4.99 = 2 5 thru 9.99 = 3 10 thru 14.99 = 4 15 thru 19.99 = 5 20 thru 29.99 = 6 30 thru 39.99 = 7 40 thru 90 = 8

Reclassifying inputs definition - how to interpret slope-rule.txt:

Slope

Category	Slope
1	< 3°
2	3-5°
3	5-10°
4	10-15°
5	15-20°
6	20-30°
7	30-40°
8	> 40°

Once we have created the slope-rule.txt file, we can select the r.reclass tool from the Processing Toolbox.

Open the Processing Toolbox: Menu Bar > Processing > Toolbar

And enable Advanced Interface by clicking the dropdowin list at the bottom of the Processing Toolbox:

From the Advanced Processing Toolbox, select GRASS commands > Raster > r.reclass
```
If you receive this error, download the reclassified slope here: wa_reclass_slope.tif


Or to enable GRASS Provider: Menu Bar > Processing > Options > Providers > (enable GRASS)


Otherwise continue on.


This means you don't have GRASS installed or bad configuration.
```
Continuing on...

We will use the following settings for the reclassify tool:
- Input Raster: WA_slope
- File containing reclass rules: slope-rules.txt
- Output raster layer: slope-reclass.tif (NOTE: select Save to file... or QGIS will not create a permanent output)
  
  Your output should look similar to the input. However, your legend should read different categories.
  
  The output layer it titled "Output raster layer" but is still saved to the hard drive. Rename this output raster to wa_slope_reclass.

Susceptibility

We will have to do some raster calculations to create the output data set. The classifications are taken from the California study[^2] and are as follows:

Since we need to cross reference these categories, we need to create a raster which represents each unique combination of the slope and rock strength classes. A method to create unique combinations is to promote the rock strength categories a decimal place. This would make the rock classes: 10,20, and 30. This is done in Raster Calculator.

We will be using both wa_geo_coded.tif, and wa_slope_reclass.tif

Open the Raster Calculator: Menu Bar > Raster > Raster Calculator

Our raster calculator expression will be:

("wa_geo_coded@1" * 10) + "wa_slope_reclass@1"

Output layer: slope_geo_sect.tif

View the new slope_geo_sect.tif
We can apply the landslide susceptibility classification to slope_geo_sect.tif using r.reclass and the following rule file:

Create landslide-rule.txt:

Copy and pase the code into a new file and save it as landslide-rule.txt

11 12 13 21 31 = 0 14 = 3 22 23 = 5 15 = 6 16 32 33 = 7 17 18 24 = 8 25 thru 28 34 = 9 35 thru 38 = 10

Open r.reclass: Processing Toolbox > GRASS > Raster > r.reclass

Configure as follows:
Input raster layer: slope_geo_sect
File containing reclass rules: landslide-rule.txt
Output raster layer: wa_landslide_suscept.tif

Notice the projection of slope_geo_sect is EPSG:2286. This is close to EPSG:2927. Read here about the difference. This is not a published analysis so we will move on.
Our final product of landslide susceptibility within the entire state of Washington.

Style wa_landslide_suscept:

Download style rules here: wa_landslide_suscept.txt

Open wa_landslide_suscept (layer not text file) style properties and import the wa_landslide_suscept.txt style rules:
We have just located areas susceptible to landslides within the state of Washginton!

References:

http://www.usgs.gov/blogs/features/usgs_top_story/landslide-in-washington-state

http://www.conservation.ca.gov/cgs/information/publications/Documents/MS58.pdf

http://www.spatialreference.org/ref/epsg/2927

https://lta.cr.usgs.goc/GTOPO30

http://www.dnr.wa.gov/ResearchScience/Topics/GeosciencesData/Pages/gis_data.aspx

http://wdfw.wa.gov/conservation/gap/land_cover_data.html

http://www.esrl.noaa.gov/thredds/catalog/Datasets/cpc_us_precip/catalog.xml#Datasets/cpc_us_precip/RT

http://gadm.org/country