Rivers contribute 95% of sediment entering the ocean. A river’s sediment load is related to many watershed-based processes generally consisting of erosion and deposition. Previous global sediment studies use coarse resolution (0.5 degree) spatial data and empirical models based on area, relief and temperature to derive estimates. Using global datasets of soils, rainfall, elevation and vegetation at a 30 arc-second (approximately 1km) spatial resolution, we were able to use the RUSLE model to predict a global erosion potential map. By summing the erosion potential by basin, we derived estimates of annual sediment load for each river. The estimates for 121 rivers were compared to observed sediment loads and were found to predict within 2 orders of magnitude in all but 3 cases. Techniques were explored to adjust the estimates for the effects of dams, reservoirs and other depositional factors but these did not increase the predictive power of the model.
The Universal Soil Loss Equation was developed by the USDA Agricultural Research Service in 1965 as an empirical method to compute soil loss in agricultural fields. The Revised Universal Soil Loss Equation (RUSLE) is based on the same empirical theory as USLE but has been improved in many areas such as methods for calculating rainfall and slope parameters(Brooks, 2003). The RUSLE model is detachment capacity limited (Mitasova, ); it only takes into account erosional factors and does not consider the effects of sediment deposition. Annual erosion is calculate by combining spatially explicit raster maps according to this formula:
E = R K LS C P
E = Annual Erosion
R = Rainfall Erosivity Factor
K = Soil Erodibility Factor
LS = Slope Factor
C= Vegetative Cover
P = Management Factor
All data referred to in subsequent sections have been preprocessed in order to normalize the spatial referencing. Because the model relies on area calculations, it was essential to use an equal-area projection. Seven regions (roughly corresponding to the continents) were defined and each uses a custom Lambert Azimuthal Equal Area projection to minimize distortion. Thus global model is run in seven steps and requires each input dataset to be reprojected into each custom projection and clipped to the region's extents. In order to minimize resampling errors, each region uses a different cell size; the approximate size of a 30 arc second cell at the center lattitude of the projection.
The SRTM30 Dataset (****) was used to calculate the slope factor (LS). Since the DEM was clipped to the sea mask and slope calculations require a 3x3 window, it was necessary to add a 5 pixel buffer with a value of zero around the existing DEM to prevent edge effects. Flow accumulation was calculated using the D-infinity algorithm. Slope percentages were calculated using a 3x3 moving window derivative approach (****). LS was calculated from slope and flow accumulation according to the modified method for complex terrain (Mitasova et. al. 1996, 1999)
Vegetative Management Factor
The MODIS land cover dataset (http://duckwater.bu.edu/lc/mod12q1.html), classified according to the IGBP scheme, was used to calculate the vegetative management factor (CP). Each class was assigned a CP factor; certain values were known ( ) and others were interpolated by their relative influence on erosion. The IGBP raster was then reclassified according to these values to derive a CP map.
|IGBP Vegetation Class||Vegetative Management Factor (VM)|
|Permanent Wetlands||0.000 |
|Snow and Ice||0.000 |
|Evergreen Needleleaf Forest||0.005 |
|Evergreen Broadleaf Forest||0.007 |
|Mixed Forests||0.009 |
|Deciduous Needleleaf Forest||0.010 |
|Deciduous Broadleaf Forest ||0.012 |
|Closed Shrublands ||0.050 |
|Woody Savannas ||0.050 |
|Open Shrublands ||0.075 |
|Grasslands ||0.075 |
|Savannas ||0.080 |
|Urban and Built-Up ||0.100 |
|Cropland/Natural Vegetation Mosaic ||0.140 |
|Barren or Sparsely Vegetated ||0.200 |
|Cropland ||0.320 |
K was calculated based on the USDA global soil suborder dataset, available through the NRCS (CITE). Each soil order was assigned a factor based on their described soil properties (CITE http://soils.usda.gov/technical/classification/taxonomy/ ) using the soil nomograph (USDA, 1980). The soil suborders were collapsed into orders and reclassified by their K factor.
| Shifting Sand||0.99|
Using interpolated monthly rainfall data (***) spanning 40 years, the average rainfall erosivity index (R) was calculated. The modified fourniers index (MFI) has shown to be a good predictor of R. The annual MFI for each of the 40 years was calculated by
sum( montly precipitation / total annual precipiatation) over each month
Raster algebra was used on the annual MFI grids to derive an average annual MFI. Finally the grid was interpolated to the optimal regional resolution using regularized spline tension (Mitasova, 1993).
Using the Database of World Rivers and their Sediment Yields (cite), we calculated the total sediment load (Yield x Area) for each river system. These rivers were observed over varying time spans and the observed values represent long term average sediment loading. In cases where multiple observation stations existed on a single river, the station furthest downstream was used to represent the sediment load at the coast. We joined the observed sediment data with the gtn basins (cite) yielding a vector polygon dataset of 123 rivers with known sediment values.Using starspan, a raster-on-vector zonal statistics tool (Rueda, 2005), we calculated the total basin-wide erosion potential by taking sum of the RUSLE values within each polygon.
R used to correlate, calculate deviation (modvobs)
Attempts to improve sediment model by incorporating auxillary data
Using the esri global dams dataset, we performed a spatial intersection with our validation basins and the total dam volume per basin was calculated. Neither the total dam volume nor the volume-to-area ratio were significantly correlated with our models deviations. The Global Lakes and Waterbodies database (GLWD) was also used to calculate the total volume of lakes and resevoirs per basin. This dataset also showed no significant correlation with our model deviation. Despite the fact that dams and resevoirs have a known effect on sediment trapping, neither global dataset had a significant effect on the model.
An additional problem is that sediment transport in rivers is a highly variable phenomenon with possible deviations of 2 orders of magnitude between years. Therefore, an annual average sediment load approach is not likely to predict the sediment loading for a given year.
Rueda, C.A., Greenberg, J.A., and Ustin, S.L. StarSpan: A Tool for Fast Selective Pixel Extraction from Remotely Sensed Data. (2005). Center for Spatial Technologies and Remote Sensing (CSTARS), University of California at Davis, Davis, CA.
Mitasova, H., Mitas, L., 1999. Modeling soil detachment with RUSLE 3d using GIS. http://skagit.meas.ncsu.edu/~helena/gmslab/erosion/usle.html. Helena Mitasova, Lubos Mitas, University of Illinois at Urbana-Champaign.Mitasova, H., J. Hofierka, M. Zlocha, L.R. Iverson, 1996, Modeling topographic potential for erosion and deposition using GIS. Int. Journal of Geographical Information Science, 10(5), 629-641. (reply to a comment to this paper appears in 1997 in Int. Journal of Geographical Information Science, Vol. 11, No. 6)
Mitasova, H. and Mitas, L. (1993) Interpolation by Regularized Spline with Tension: I. Theory and Implementation, Mathematical Geology, 25, 641-655.