Using the Universal Soil Loss Equation

Knowing where erosion occurs is vital to mitigating the problem and protecting against future soil loss. Field surveys can be time-consuming, but using GIS to compare datasets representing different factors associated with erosion can provide an encompassing snapshot of at-risk areas. In addition, land managers can use Global Mapper and free online data to apply the Universal Soil Loss equation to a site and find where erosion mitigation may be needed.

Erosion affects farming, construction projects, and cities near oceans, rivers, and terrestrial slopes. To help predict soil erosion due to water, The Universal Soil Loss Equation was developed by the US Department of Agriculture in 1965 and is hailed as one of the most significant developments in soil and water conservation in the 20th century (USDA). This formula measures the average annual rate of soil loss by erosion on a field slope based on soil type, rainfall pattern, topography, management practices, and crop system. More about the history of the Universal Soil Loss Equation can be found here on the USDA website.

The formula for applying the Universal Soil Loss Equation in a GIS is:

A = R*K*(LS)*C*P

Where:

A = Average annual soil loss (ton/ha/yr)
R = Rainfall erosivity factor. Runoff Factor for soil erosion (MJ mm ha-1 hr-1 yr-1)
K = Soil erodibility factor (ton ha hr MJ-1 ha-1 mm-1)
LS = Topographic Steepness factor (based on length and slope)
C = Land Cover Management factor
P = Erosion control practice factor

The benefit of using this formula in GIS is that the formula can easily be run over large areas in a short amount of time. Each of the values in the formula are represented by a raster layer. Rasters break the surface of the earth into grid cells that contain a single value for the area they cover using any numeric metric. Lower resolution rasters (large cell sizes) can be less precise as they have less detail, representing larger swaths of area with a single value. Keep this in mind when finding data to suit your study area.

These instructions will cover how to calculate the Revised Universal Soil Loss Equation formula in Global Mapper. It helps point you in the direction of data required for each variable, and how to run it through Global Mapper’s toolset to generate the layers and final output.

There are a few different methods based on the source of data available for your study site, and this scenario aims to cover a few of them, with the goal being that others may also easily be applied. Each variable will include a US data source and a Global data source, each with the necessary instructions.

Let’s look at how to use the universal soil loss equation in Global Mapper.

*Tip: In Global Mapper, to crop raster data to a vector feature, first select the vector feature in the workspace, then double-click on the raster layer(s) in the control center to open Raster Options. From there, find the Cropping tab; this will include the option to Crop to Selected Features. The cropped data is not deleted, the crop can be undone, but it will not be included in any analysis or export.

R – Rainfall erosivity factor

The rainfall erosivity factor represents the erosive powers of rainfall energy and runoff by measuring the intensity, quantity, and duration of rainfall. The greater the intensity and duration of the rain storm, the higher the erosion potential. This relates to a few different types of erosion, such as splash erosion, caused by rainfall, which can loosen soil particles. This data can be found premade or calculated from rainfall data.

US Data Source:

Precalculated data for the contiguous United States can be downloaded from NOAA. This website was down as of the posting of this workflow.
https://catalog.data.gov/dataset/r-factor-for-the-coterminous-united-states

Global Data Source:

Precalculated Global Data can be requested from the European Soil Data Center (ESDAC)
Global Rainfall Erosivity. This data can be downloaded as a premade raster layer with worldwide coverage. The downside is that this is low-resolution data and has been known to overestimate.

If higher resolution data is needed, download rainfall data for the area and use a formula to generate the rainfall erosivity factor. Rainfall is typically measured at specific locations and is often distributed as a point layer. In Global Mapper, you can use the Create Elevation Grid tool to interpolate these points into a solid surface raster layer using rainfall as the preset elevation value. Next, use the Raster Calculator tool to calculate the erosivity factor. Given the mean annual rainfall from 5 stations, R can be calculated using the formula (Wischmeier and Smith, 1978) R = 0.5*P*1.73, where P is the mean precipitation in the area in inches.

ESDAC R layer loaded and cropped in Global Mapper. Transparency settings were used here for display purposes. The sample study area for this post is the Missouri River along the South Dakota border.

*Tip: To view cell values in a raster, select the Feature Info tool and click on a cell. Values that Global Mapper is interpreting as elevation will appear in the Elevation section.

K – Soil erodibility (soil characteristics)

Soil erodibility describes the susceptibility of a soil type to erosion, measured by the soil’s physical properties such as texture, structure, and permeability. Soil types respond to an erosive event differently based on their particle size. The larger the particle size, the more soil will be lost. For example, sandy soil (large particles) will erode far more quickly than clay soil (small particles).

Values of K range from the lowest erodibility, 0.02, to the highest, 0.69. All other factors being equal, the higher the K value, the greater the susceptibility of the soil to rill and sheet erosion by rainfall.

US Data Source:

Downloaded data from the USDA NRCS Web Soil Survey includes precalculated K values. Because of the high volume of attributes, they are stored separately from the shapefile and must be joined in Global Mapper. In order to join a text file to a shapefile, there must be a common attribute that Global Mapper can use to associate each attribute to its specific location. Here we are going to use a third file that has attributes to match both the shapefile and the K attribute to bridge them together.

  1. Acquire data:
    First, download data for your desired study area from the Web Soil Survey: http://websoilsurvey.sc.egov.usda.gov/App/HomePage.htm
    You don’t need to specify any particular type of data. Use the Download Data tab to get the data itself instead of generating a map (via the cart). More instructions on how to do that are here.
  2. Format comp.txt:
    From the Tabular folder, open comp.txt in Excel. Use the Data >Text to Column to delineate the features by pole (|). For ease of use, delete everything other than the last two columns. These are MUKEY and COKEY, in that order. (More information on what all of the attributes are can be found here.) Label them in a new row. Save it as a new file.
  3. Format chorizon.txt:
    Next, open chorizon.txt in Excel. After delineating, columns DH and DI in Excel are the K factors. DH is Kw (adjusted for the effect of rock fragments), and DI is Kf (not adjusted).
    Delete all columns other than your chosen K factor and the second to last column, which is COKEY. Label them in a new row. Save this as a new file.
  4. Join the files in Global Mapper:
    In Global Mapper, open soilmu_a_aoi.shp. Join your edited comp.txt to soilmu_a_aoi using the settings below (Based on MUKEY). This adds COKEY as an attribute. Then, join Join chorizon.txt (based on COKEY).

  5. Create the raster:
    Open Vector options for the shapefile. In the Elevation tab, set your K value as the elevation attribute. Use the Create Elevation Grid tool to create a raster named K. This tool will create a raster of whichever attribute is set as elevation in a dataset.
  6. You are now ready for the next variable, LS.
K as a finished raster.

Global Data Source:
Globally, the K variable can be downloaded or streamed directly from the United Nations Agriculture Program (UNAP). By streaming the data directly into Global Mapper, you can limit the data to your study area without downloading the entire world and cropping it down.

Here are the settings to streaming K into the Online Data tool from the UNAP

To add this as a custom source in the Online Data tool, click Add New Source, choose WMS, then fill out the dialog as shown below. More information on this process can be found here.
https://data.apps.fao.org/map/gsrv/gsrv1/geonetwork/wms
Layer name dsmw_14116

Now that the data is loaded in Global Mapper, you can calculate K. For more information on how to use this dataset, see this video with step-by-step instructions from GeoDev.

(LS) – Slope Length and Steepness

L is slope length, and S is slope steepness. Also called the Topographic Factors, these factors together alter runoff speeds. Steeper and longer slopes can cause an increase in water speeds, which can lead to an increase in soil loss. To meet the limitations of GIS, LS is calculated by the USPED (Unit Stream Power Erosion and Deposition) method, which uses a raster calculation between flow accumulation and slope. First, generate the layers needed for L and S, then use the Raster Calculator tool to create the final LS layer.

The data needed for this step is a Digital Elevation Model (DEM).

US Data source:

For higher-resolution data, download a DEM from USGS’s Earth Explorer. Use the shortcut found in the Online Data tool under Terrain to automatically open the explorer in our default browser, loaded to match the zoom extent and location in Global Mapper.

Global Data Source:

If you have an Intermap subscription, Next Map from the online data tool is a great resource for high-resolution terrain data.
The Online Data tool also has free resources under the Terrain section, including SRM Worldwide Elevation Data and ASTER GDEM.

1. L factor:

The L factor (slope length factor) is the soil loss ratio from a slope length relative to the standard erosion plot length of 22.1m.

L= [(FA * cell size)/22.13]^m
Where:
FA = Flow Accumulation
Cell size = resolution of the Digital Elevation Model (DEM)
m ranges from 0.2-0.6.

1.2 Calculate Flow Accumulation (FA) in Global Mapper

With your DEM open in Global Mapper, open the Watershed tool. This tool maps the movement of water across the terrain, downhill from cell to cell. We will take advantage of a side function in this tool, Create Flow Points with Symbol, to create a point layer where each point contains Flow Accumulation as an attribute. The point layer can then be gridded into a raster.

In the watershed tool, check the Create Flow Points with Symbol box in the lower right corner.

Many of the settings in this tool aren’t relevant to our purpose, but be sure to set the Resolution on the right side to match the resolution of your DEM (or at least the lowest resolution layer in your dataset) if not automatically populated.

This will create two layers. Delete the unnecessary line layer. The second layer, the point layer, will only appear when you zoom in very closely, but it is always there. Change its Vector options > Elevation tab to the Flow Accumulation attribute, then use the Elevation grid tool to create a raster of those values.

Flow accumulation layer: streams will have higher values than the surrounding areas.

Now that we have the necessary layer for the L factor, we can calculate it simultaneously with the S factor.

2. S Factor

The S factor (slope steepness factor) is the ratio of soil loss relative to a 9% slope, which
is the standard slope that experiment plots use.

S = [(sinB * 0.01745)/0.09]^n

Where:
B = the slope angle in percent
n ranges from 1.0-1.3

2.1 How to Calculate Slope in Degrees in Global Mapper

Global Mapper will automatically display the slope of a raster when the Shader is changed to Slope Shader. If you change it from the dropdown, the change will apply to the entire workspace, and that would affect the other factor layers. Instead:

  1. Right-click on the DEM and open Raster Options.
  2. In the Display tab, change the Shader dropdown from Default to Slope Shader. This layer will use the slope shader no matter what the global shader is changed to.

3. Calculate L*S

In Global Mapper, when using a formula in the Raster Calculator, instead of layer names, you will use B variables such as B1, B2, etc. Once the formula has been completed and you click OK, you will receive a prompt to associate each B variable with a Raster Layer. If the desired layer does not appear, it means the layer is either turned off in the workspace or excluded from the original tool. Close Raster Calculator and start over, making sure that all desired layers are checked.

To calculate LS in Raster Calculator

  • Open Raster Calculator – Apply Formula from the Analysis dropdown menu
  • Click Add Custom Formula.
  • Copy and paste this into the Custom Formula box:

(((B1*70)/22.1)^0.4)*(((SIN(B2*0.01745))/0.09)^1.4)

  • Click Add Band. The formula will have been added as a band.
  • Click OK
  • Band 1 (B1) is flow accumulation
  • Band 2 (B2) is the DEM slope layer
The LS factor is higher in streams and areas with steeper slopes.

C = Cover factor

The Cover factor is a measure of plant and/or crop cover in an area. The existence of vegetation in an area reduces soil loss, so more ground cover typically means there’s less soil loss.
There are two ways of measuring ground cover, normalized difference vegetation index (NDVI) using Landsat imagery, and Land Use/Land Cover (LULC). This workflow will cover LULC.

US Data Source:

LULC data can be streamed directly from the Online Data tool under NLCD US National Land Cover Dataset. In the US, it’s produced by the USGS. This dataset quantifies how land is being used into preset categories such as urban, shrubland, agriculture, forest, etc. The full legend will appear with the data.

The National Land Cover dataset with a legend for each cover type.

Global Data Source:

If land Cover data isn’t available for your study area, you can use NDVI. Landsat 7 imagery can be accessed freely from the USGS Earth Explorer’s website for free. From that data, Bands 3 and 4 are used to calculate NDVI. For more information on how to calculate NDVI from Landsat 7 imagery, see the link attached (Calculating NDVI in ArcGIS using Landsat 7 satellite image – YouTube).

Calculate Cover Factor

The NLCD data streamed from the online data tool has labeled the land based on land use.

Each cover type has a value and a description (for example open water =11). We will use the Raster Reclassify tool to assign C values to each land cover type.

  1. Connect to NLCD data from the Online Data tool.
  2. Export this layer as a global mapper package file from the layer menu, and reimport.
    1. Change the resampling method in Raster Options to No Resampling before export.
    2. GMP files are the best option when exporting online data as they only include one zoom level.
  3. From the Analysis Menu, open Raster Reclassification.
  4. Choose Load Palette from Layer
  5. In the Raster Reclassification Rules, click on a row to add a new classification for each land cover type based on the C values in the table below.
  6. Here, the colors were reassigned with classification values different from the standard, so the land classification type was determined based on the color.
  7. Click OK to reclassify the data.

P = Management Factor

This value represents practices farmers implement to reduce soil loss, such as contour cultivation to slow downhill water movement where furrows are perpendicular, not parallel, to the slope of a hill. This helps slow runoff during rainfall, which allows more water to soak into the ground. In areas where the practice varies or is unknown, the default value for P is “1” where straight-row farming practices are assumed.

Final Calculation
To generate the final calculation layer, multiply the layers together in Raster Calculator. The Universal Soil Loss equation is:

A = R*K*(LS)*C*P

  1. In Raster Calculator, copy and paste this formula:
    B1*B2*B3*B4
  2. Click Add Band. Click OK.
  3. Associate each band with the R, K, LS, and C layers in any order.

The final layer will highlight which areas in the terrain are most at risk for erosion based on the data.

Here, sloped urban areas are at higher risk than the surrounding forested areas.

Knowing which areas are at higher risk for erosion can help guide mitigation plans. Berms, Are you interested in adjusting the terrain to change the water flow? Try modeling the mitigation in the terrain with the Terrain Paint tool. Would you like to measure how much soil has been lost over time? Subtract two elevation datasets of the same area captured at different times to measure the difference in elevation. Use the Combine/ Compare Terrain Layers to compare elevation rasters like DEMs, and the Compare Point Clouds tool to measure lidar differences.



WORK MADE EASY WITH GLOBAL MAPPER PRO

Want to try Global Mapper? Sign up for a 14-day free trial. You can also request a demo from one of our experts to see this workflow or other Global Mapper processing abilities.


Learn More


References:

Kim, Y. (2014, December 5). Soil erosion assessment using GIS and revised universal soil loss ... Soil Erosion Assessment using GIS and Revised Universal Soil Loss Equation (RUSLE). Retrieved March 30, 2023, from caee.utexas.edu/prof/maidment/giswr2014/ProjectReport/Kim.pdf.

Oppong, J. (2022, June 12). How to use ArcGIS pro to estimate soil erosion from a catchment basin. GIS Lounge. Retrieved March 28, 2023, from gislounge.com/arcgis-pro-soil-erosion-catchment-basin/.

Wischmeier, W.H. and Smith, D.D. (1978) Predicting Rainfall Erosion Losses. A Guide to Conservation Planning. The USDA Agricultural Handbook No. 537, Maryland.

Companies using Blue Marble’s geospatial technology