# Extending' the RUSLE with the Monte Carlo error propagation technique to predict long-term average off-site sediment accumulation.

ABSTRACT: To evaluate if the adaptation of the basically two-dimensional Revised Universal Soil Loss Equation (RUSLE) to a three-dimensional reality is appropriate for predicting off-site sediment accumulation, it was extended with the Monte Carlo error propagation technique. This technique generates the true probability distribution of model output and gives the possibility to explain whether the difference between the model output and the field observations is largely due to the uncertainty of the model input or is mainly due to the uncertainty and limitations of the model itself It was found that the RUSLE was able to accurately predict off-site sediment accumulation in the water reservoir of a study area. The value of the measured sediment input was within the 68% confidence interval around the predicted value, with a difference of only 1.4%. Therefore, the error propagation explained this difference as mainly due to the uncertainty of the model input parameters. Consequently, it can be concluded that the topographic factor of the RUSLE model also can be considered as a measure of the sediment transport capacity of the overland flow, although it was originally developed for situations where detachment limits the sediment load.Key words: Error propagation, Monte Carlo, RUSLE, water erosion.

Abundant rains in the winters of 1993-94 1994-95 resulted in substantial on-site and off-site erosion problems in Belgium. Therefore, the Regional Land Management Board of the hilly region in the south of West-Flanders requested a scientific study of these problems. This Board was installed by Ministerial Decree to preserve and protect the landscape of this area, which is considered to be of an exceptional value. The motivation for this request came from the need for information to support "land management agreements" between farmers and this Board.

Although process-based erosion models, such as the Water Erosion Prediction Project (WEPP) (Flanagan and Nearing 1995) are being developed to replace the empirical models (Laflen et al. 1991), the Revised Universal Soil Loss Equation (RUSLE) (Renard et al. 1996) was selected as a basis to develop an erosion expert system. The latter should be able to assess both on-site soil losses and off-site sediment accumulations. This system was implemented in ANSI C (American National Standards Institute), intended to be used inside a Geographic Information System (GIS) environment. The RUSLE was chosen because it requires only a limited amount of data to perform a field-scale erosion analysis for large areas, compared to other process-based models.

To predict the off-site sediment accumulation in rivers or water reservoirs, the equations in the erosion model describing the hydrological processes must be a measure of the transport capacity of the overland flow. Only when this is realized can the amount of sediment leaving a field be estimated. The RUSLE is a factor-based erosion model designed to predict longterm average soil losses carried by runoff from specific field slopes in specified cropping and management systems (Renard et al. 1996). The factor that summarizes the hydrological components of the water erosion process is the topographic factor LS. Foster and Wischmeier (1974) stated that the LS equation, derived from unit field plots with a length of 22.13 m (72.6 ft) and 1.80 m (5.9 ft) width, applies to situations where detachment limits the sediment load, and is not a transport-capacity equation. However, recent rill erosion experiments conducted by Nearing et al. (1997) indicate that transport capacity in eroding tills is already reached with in a sample length of 2.5 m (8.2 ft) for slopes ranging between 3 to 28%. These results support the statement of Moore and Wilson (1992) that the equations of the topographic factor in the RUSLE model also are a measure of the sediment transport capacity of overland flow. This implies that RUSLE can estimate the sediment actually leaving a field and does not account for deposition as colluvium (intrabasinal storage).

The end product of a model always is the result of operations and computations performed on uncertain data. In GIS studies, every layer of information has its associated uncertainty caused by different sources of variance (Heuvelink et al. 1989). Some of these sources can be unavoidable (e.g., the intrinsic variability of the climate or the uncertainty associated with interpolation methods). If the results of a model are not in agreement with field observations, it is important to know whether this is due to the model itself or to the uncertainty of model input. This can be evaluated with the Monte Carlo error propagation technique (Wesseling and Heuvelink 1993). With this technique, the model output is generated at least a few hundred times, but instead of using the parameter values, their stochastic distributions are used. This allows a modeler to determine the stochastic properties of model output.

The factor of the RUSLE model that poses the most problems in the error propagation process is the LS factor. Desmet and Govers (1996) developed a method to calculate the LS factor on topographically-complex landscapes within a GIS, based on the unit contributing area. For this study, an alternative method to calculate the LS of a field was developed. This alternative method can be used in the Monte Carlo analysis and is more closely related to the linear structure of the RUSLE model and the linear micro-topography that can be observed in the field. The linear micro-topography is created by cultivation techniques: crops are cultivated in rows and the use of heavy machinery creates linear furrows, which initiates rill erosion. This induces a parallel flow pattern and prevents the concentration of the Horton runoff into bigger rills and gullies. These features cannot be captured by the unit contributing area because the resolution of most

Digital Elevation Models (DEM) is too low to describe this linear micro-topography.

This paper has two main objectives: 1) to evaluate if the RUSLE can be used to predict long-term off-sire sediment accumulation, which is equivalent to checking if the RUSLE LS equations are a measure of the transport capacity, and 2) to perform a Monte Carlo error propagation to determine the uncertainty of the calculated on-site soil losses and off-site sediment accumulation.

Materials and methods

The procedure described below was applied to the Kemmelbeek Watershed (Figure 1), located in the south of West-Flanders, within Belgium's loess belt. It covers an area of 1,075 ha (2,655 ac) and feeds a drinking-water reservoir for the city of Ieper. The highest elevation is 151 m (495 ft), which drops to 23 m (75 ft) at the reservoir inlet. The average slope steepness is 4.6%, with a maximum slope, of 71%, although 99% of the slopes are less than 30%. This watershed was chosen as a pilot test area to determine the on-site soil losses and off-site sediment accumulation in the reservoir, using the adapted RUSLE expert system. The model predictions could be validated, since data on the sediment input in the reservoir are available. The reservoir has a sediment-trapping efficiency of nearly 100%. This is based on the fact that 1) most incoming sediment is deposited directly behind the reservoir inlet; 2) the daily amount of water used for drinking water production is a volume of 4,000 [m.sup.3] (141,200 [ft.sup.3]), which corresponds with a waterlayer of only 1.15 cm (0.45 in) for a total reservoir area of 34.88 ha (86 ac); and 3) excessive water is only pumped into the downstream drainage system when the reservoir exceeds a certain critical level. (The pumps are located 700 m (2,296 if) from the reservoir inlet).

In a GIS environment, a model can be written as (Wesseling and Heuvelink (1993):

M = F([Z.sub.1], [Z.sub.2],..., [Z.sub.n], [a.sub.1], [a.sub.2],..., [a.sub.m]) [1]

where the resulting map M is obtained by applying expression F on input maps [Z.sub.i] and model coefficients [a.sub.i]. Because most model input is subject to uncertainty, not all parameters of F are exactly known. Therefore, the parameters of F must be represented by probability distributions rather than by deterministic quantities. A Monte Carlo error propagation technique can then be used to determine the distribution of the model output. This technique repeatedly (e.g., 500 times) runs the model with input values that are sampled from their distributions. If the number of runs is sufficiently large (depending on the model complexity), the distribution obtained from the runs will approximate the true distribution of the model output.

Before an error propagation can be performed, the following must be determined (Wesseling and Heuvelink 1993):

1. the properties of function F;

2. which parameters are stochastic and which are deterministic;

3. the probability distribution of the stochastic parameters;

4. the correlation between the different parameters at the same location; and

5. the correlation between the spatial parameters at different locations.

In this study, the function F in Equation 1 is the RUSLE model. To test if the elaborated methodology of the RUSLE expert system is capable to predict off-site sediment accumulation, it is hypothesized that the RUSLE model itself induces no error. The uncertainty of the output is then only induced by the uncertainty of the model input parameters: LS, R, K, C, and P. Therefore, the stochastic properties of these input parameters must be determined.

The topographic factor (LS)

The RUSLE topographic factor describes the combined effect of slope length (L) and slope steepness (S). Because a terrain element downslope gets more runoff water than a terrain element near the water divide, Foster and Wischmeier (1974) subdivided a slope into a number of uniform segments. The LS of a flowline in the landscape can be calculated by:

LS = [[sigma].sub.i] [([[L.sup.[M.sub.i]+1].sub.i] - [[L.sup.[M.sub.i]+1].sub.i-1] * [S.sub.i]/[L.sub.tot] * [22.13.sup.[M.sub.i]] [2]

where [L.sub.tot] is the total length of a flowline [m], [L.sub.i] is the length from the top of the slope to the foot of segment i [m], [L.sub.i-1] is the length from the top of the slope to the top of segment i [m], [S.sub.i] is the slope steepness factor for segment i [-], and [M.sub.i] is the slope-length-exponent [-]. The slope-length-exponent can be written as (McCool et al. 1989; Renard et al. 1996):

M = [beta]/1+[beta] [3]

[beta] = [sin([alpha]) / 0.0896/3.0 [(sin([alpha])).sup.0.8] +0.56] * r [4]

where [beta] is the rill/interrill ratio [-] and [alpha] is the slope steepness [radians]. When field conditions favor rill erosion r = 2 [e.g., on ridged potato (Solanum tuberosum) fields]; when field conditions favor sheet erosion r = 0.5 (e.g., a field that remains bare for a long time); and for in-between conditions r = 1. The main crops in the study area are beets (Beta vulgaris), maize (Zea mays), potatoes (Solanum tuberosum), and wheat (Triticum aestivum). The field conditions for these crops favor rill erosion. Therefore, an r-value of 2 was used, except for pasture fields and forested areas where r was set to 1.

The slope steepness factor (for a slope length longer than 4.56 m, or 15 ft) is given by this equation (McCool et al. 1987; Renard et al. 1996):

S = {10.8 sin([alpha]) + 0.03 if (slope [less than] 9%)

S = {16.8 sin([alpha]) - 050 if (slope [greater than or equal to] 9%) [5]

Due to the small parcels in the study area (the average field area is only 1.40 ha or 3.46 ac) and the dense ditch system, every field can be considered as a separate hydrological unit concerning Horton flow. Consequently, the flowlines in the landscape start at the upper field boundaries and end at the lower field edges, where the runoff water flows into the drainage system. The calculation process for the LS factor is shown in Figure 2.

For example, suppose a square field or hydrological unit has a length and width of be 80 m (262 ft). In a 10-m (32.8-ft) resolution DEM, this results in a grid with 8 rows and 8 columns. Figure 2a gives the drainage directions for each elementary cell of the grid. These drainage directions are used to construct linear flowlines in the field (Figure 2b). The flowline matrix (FLN), given in Figure 2c, indicates the number of flowlines running through a cell. The representative area for each segment of a flowline is the cell area divided by the FLM value for that segment.

The soil erosion [t . [yr.sup.-1]] along a flowline, [A.sub.f], can then be written as:

[A.sub.f] = [[[sigma].sup.N].sub.i=1] [([[L.sup.[M.sub.i]+1].sub.i] - [[L.sup.[M.sub.i]+1].sub.i-1])*[S.sub.i]/[L.sub.tot] * [22.13.sup.[M.sub.i]] *[a.sub.i]*[R.sub.i]*[K.sub.i]*[C.sub.i]*[P.sub.i]] [6]

where N is the number of segments in a flowline, [a.sub.i] is the representative area of a segment in a flowline [ha], and [R.sub.i], [K.sub.i], [C.sub.i], and [P.sub.i] are the other RUSLE factors, respectively, the rain erosivity [MJ * mm * [ha.sup.-1] * [h.sup.-1] * [yr.sup.-1]], the soil erodibility [t*ha*h*[ha.sup.-1]*[MJ.sup.1]*[mm.sup.-1]], the cover management factor [-] and the support practice factor [-]. If dx and dy are the cell dimensions [m], the representative area [ha] of a segment in a flowline can be calculated by:

[a.sub.i] = dx * dy/10,000 * [FLM.sub.i] [7]

The total soil erosion [t * [field.sup.1] * [yr.sup.1]] of a field or hydrological unit, [A.sub.field], results from:

[A.sub.field] = [[[sigma].sup.F].sub.f=1] [A.sub.f] [8]

where F is the number of flowlines in a field or hydrological unit.

There are two sources of variance in the LS algorithms that influence the uncertainty on the predicted sediment loss:

1. The uncertainty on the elevations in the DEM. For every iteration in the error propagation process, an error surface was created and added to the original DEM. If the error between the DEM elevation and the real elevation for a certain point is high, the error also will be high for positions in the neighbourhood of that point. Therefore, the error surfaces must be autocorrelated over a short distance. This requires the construction of a random generator, which can create autocorrelated error surfaces.

2. The rilling pattern that can be observed in a field is never the same for every storm event. To simulate this randomness of the Horton flow, a stochastic flow routing model was used to determine the flow directions.

Construction of autocorrelated surfaces. The fractional Brownian motion (fBm) and fractional Gaussian noise (JGn) serves as the basis for many models for natural fractal shapes such as landscapes (Polidori 1991; Peitgen et al. 1992). In its one-dimensional form, the fractional Brownian motion model is defined as a continuous function, fix), of the independent spatial variable, x, having the following properties (Molz and Liu 1997; Peitgen et al. 1992):

1. The increments of/are stationary. This means that for all values of x and a fixed increment, h

E[f(x + h) - f(x)] = [C.sub.1](h) [9]

E[[(f(x + h) - f(x)).sup.2]] = [C.sub.2](h) = [gamma](h) [10]

where E[X] expected value of the random variable X, [gamma](h) is the variogram, and [C.sub.1] and [C.sub.2] are functions of h. Because the increments of f(x) are defined to be stationary, one can define:

n(x,h) =f(x + h) -f(x) [11]

with the statistical properties of n depending only on h.

2. The variable n has a Gaussian distribution with:

E[n(x,h)] = 0 and E[[(n(x,1)).sup.2]] = [[sigma].sup.2] [12]

Thus, in equations [9] and [10], [C.sub.1](h) = O for all h, and [C.sub.2](1) [[sigma].sup.2].

3. The increments, n(x,h), are statistically invariant with respect to an affine transformation. This means that the random variable n(x,rh) and [r.sup.H] n(x,h), with r and H constants (0 [less than or equal to] r [less than] [infinity] and 0 [less than]H [less than] 1), have the same Gaussian distribution. Thus:

E[n(x,rh)] = E[r.sup.H] n(x,h)] = 0 [13]

E[[(n(x,rh)).sup.2]] = E[[([r.sup.H] n(x,h)).sup.2]] = [r.sup.2H] E[[(n(x,1)).sup.2]] [14]

Using Equations [12], [13], and [14], one can write:

E[[(n(x,h)).sup.2]] = [h.sup.2H] E[[(n(x,1)).sup.2]] = [h.sup.2H] [[sigma].sup.2] = [gamma](h) [15]

In order to define fBm = f(x), it was necessary to define the properties of its increments, n(x,h). This last function is known as fGn. The concepts of fBm and fGn are generalizations of the classical concepts of Brownian motion and Gaussian noise, denoted as cBm and cGn. These functions are obtained by setting H = 0.5. The parameter, H, is the so-called Hurst coefficient and is an indicator of the surface complexity. H is related to the fractal dimension of the surface and can be estimated by fitting a power function, y = a*[x.sup.b], on the variogram. According to Equation 15, a [approximately equals] [[sigma].sup.2], b [approximately equals] 2H, when y is the semivariance and x is the lag distance. The fBm can be divided into three distinct categories: H [less than] 0.5, H = 0.5 and H [greater than] 0.5. The case H = 0.5 is the ordinary Brownian motion, which has independent increments. For H[greater than] 0.5, there is a positive correlation between the increments. For H [less than] 0.5, there is a negative correlation between the increments and the curves or surfaces seem to oscillate more erratically. As an illustration, Figure 3 gives 3 fBm series with a Hurst coefficient of 0.2, 0.5, and 0.8 with their respective variograms. Random surfaces with a given Hurst coefficient can be created by the "random midpoint displacement" algorithm (Peitgen et al. 1992).

Flow-routing algorithm. Several flow-routing models exist to determine the flow directions in a grided DEM: the deterministic 8-neighbors method (D8) (O'Callaghan and Mark 1984; Jenson and Domingue 1988); the stochastic 8-neighbors method (Rho8) (Fairfield and Leymarie 1991); multiple-direction methods (Quinn et al. 1991; Freeman 1991); and DEMON (Costa-Cabral and Burges 1994). For this application, the Rho8 method, which is basically a stochastic extension of method D8, was chosen based on its stochastic character to determine flow directions.

In method D8, each pixel discharges into one of its 8 neighbors. The drainage direction is determined by the direction of the largest weighted elevation drop (LD), calculated by:

LD = ([H.sub.C] - [H.sub.n]) * [rho] [16]

where [H.sub.c] and [H.sub.n] is the elevation of the center cell and the neighbor cell, respectively, and [rho] is the weight factor for the direction and equals: 1/dx in the x direction, 1/dy in the y direction, and 1/sqrt([dx.sup.2] + [dy.sup.2]) in the diagonal directions.

If dx equals dy, [rho] is equivalent to 1 for the cardinal directions, and 1/sqrt(2) for the diagonal directions. This method gives maximum errors in the flow direction of 22.5 degrees in planar areas which are not aligned with the grid orientation (Fairfield and Leymarie 1991). To solve this problem, method Rho8 gives the weight factor [rho] a stochastic character. For the cardinal directions, [rho] equals 1; for the diagonal directions, [rho] ranges between [0.5,1.0] with a mean value of 1/sqrt(2), according to its cumulative distribution function (cdf):

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] [17]

Let r be a uniformly distributed random variable between [0, 1], then the inverse of the cdf, [rho] = 1/(2 - r), generates random weight factors for the diagonal directions, with a mean value of 1/sqrt(2). Using this flow-routing model, every iteration in the Monte Carlo error propagation process generates another flow pattern.

To calculate the LS factor, a DEM of the study area was interpolated from the contour lines, digitized from the Belgian topographic maps with a 1:10,000 scale. The grid resolution was chosen as 10 m (32.8 if). The field boundaries were digitized from the 1:10,000 orthophotos. The random midpoint displacement method, used in the Monte Carlo simulations, requires two parameters: the Hurst coefficient, indicating the complexity of the DEM error surfaces, and a standard deviation, indicating the dimension of the errors. To determine the Hurst coefficient, the elevation was measured along a transect down a hillslope. These measurements were compared with the elevations in the DEM.

Figure 4a gives the trace of the errors along the transect, and Figure 4b shows the respective variogram of this error trace. Fitting a power function to the variogram results in a Hurst coefficient of 0.83, indicating a positive correlation between the errors at successive points. The interval of the digitized contour lines, used to interpolate the DEM, was 2.5 m (8.2 ft). Consequently, the maximum error of the elevations in the DEM is 2.5 m, which can be considered to be normally distributed. The standard deviation, used in the random midpoint displacement method, is then approximately one-sixth of the contour interval, or 0.42 m (1.38 if), according to (with f(x) the normal probability density function):

[[[integral of].sup.+3*[sigma]].sub.-3*[sigma]] f(x) * dx = 0.997 [18]

Figure 5 gives the results of the Monte Carlo error propagation: the mean LS value of a field and the uncertainty of these estimations, expressed by the standard deviation.

Rain erosivity (K)

Precipitation over a 27-yr span was used to calculate the rain erosivity index. Over this period, rain intensity was recorded every 10 min. Equations used to calculate erosivity can be found in Renard et al. (1996). For the Kemmelbeek Watershed, the mean yearly rain erosivity is 724 [MJ*mm*[ha.sup.-1]*[h.sup.-1]*[yr.sup.-1]]. The natural variability of the precipitation characteristics (intensity and amount) is very large. Consequently, the standard deviation of the R-factor also is very large: 224 [MJ*mm*[ha.sup.-1]*[h.sup.-1]*[yr.sup.-1]]. Because the RUSLE model predicts the mean soil loss over long time periods (a few decades), it is not necessary to take this uncertainty into account in the model calculations. Therefore, the R-factor was considered to be deterministic, with a value of 724 [MJ*mm*[ha.sup.-1]*[h.sup.-1]*[yr.sup.-1]].

Soil erodibility (K)

The soil erodibility factor can be calculated by the following equation (Renard et al. 1996):

K=[2.1*[(S*(100-C)).sup.1,14]*[10.sup.-4]*(12-OM)/100]

0.1317 [19]

where K is the soil erodibility [t*ha*h*[ha.sup.-1]*[MJ.sup.-1]*[mm.sup.-1]], S is the textural fraction between 2 and 100 [micro]m [%], and OM is the organic material content [%].

The Belgian soil map contains only qualitative information, so it is not appropriate to convert this into quantitative information. Therefore, 153 locations were sampled within the watershed to determine the textural fractions and the organic material content. Using equation [19], the erodibility of the samples was calculated and block-kriged (Webster and Oliver 1990; Van Meirvenne 1991) using the variogram given in Figure 6. Blocks 10 m (32.8 ft) square were used, resulting in an erodibility grid (Figure 7a) of the same resolution as the DEM. The uncertainty of this interpolation is expressed by the kriging variance (Figure 7b). The K factor can be considered to be normally distributed with the kriged value as mean and the kriging variance as a measure of the spread of the estimation error.

The cover management factor (C)

The dimensionless C factor, which has a range between 0 and 1, expresses the degree of protection of the soil surface by the crops or vegetation. The information for calculating the C factor was obtained from a detailed inquiry with farmers (Ghekiere 1997). However, there was not enough information (e.g., soil biomass, crop residues, soil consolidation, etc.) to use the RUSLE methodology. Instead, the Universal Soil Loss Equation (USLE) was used to calculate the C factor for every crop rotation.

For every crop, the growing season was subdivided in 6 stages. However, not every crop cultivated in western Europe can be found in the USLE crop database. The C value for these crops must be estimated using the data of similar crops. Therefore, for every growing stage, a possible minimum and maximum C value was chosen from the USLE C factor table [see Table 5 in Agriculture Handbook No. 537 (Wischmeier and Smith 1978)]. These values were weighted with the erosivity value of that period, which resulted in two maps (Figure 8) that represent the minimum C factor and the maximum C factor for a field. The C-factor of a field can be considered to be uniformly distributed between this minimum and maximum C value.

The support practice factor (P)

The dimensionless support practice factor, P, takes into account the effect of special management practices, such as strip cropping and terraces. Because no farmer applied such soil conservation practices, the P factor has a deterministic value of 1 for the entire study area.

Soil erosion (A)

Once the stochastic distribution of every RUSLE parameter is determined, the RUSLE model can be represented by Equations 6 and 8. Because the stochastic RUSLE parameters (LS, K, and C) are mutually independent, there are no correlation terms that should be taken into account in the error propagation. To assess the soil loss of every field or hydrological unit and the uncertainty on these results, a Monte Carlo simulation, with 500 runs per field, was performed. For every run, a new DEM was created for the field, and the K and C values for every pixel in the field were randomly sampled according to their probability distribution. If the number of iterations is high enough (at least a few hundred), the soil loss will be distributed normally. This resulted in two maps (Figure 9) that show the soil loss per field and the uncertainty of these estimates, expressed by the standard deviations.

Validation

Using the soil erosion map and the standard deviation map on the soil loss, one can calculate (using a Monte Carlo simulation) the average yearly sediment input in the water reservoir and the uncertainty on that value. The calculated yearly sediment input was 4,376 [t*[yr.sup.-1]] with a standard deviation of 75 [t*[yr.sup.-1]]. Given the uncertainty of the model input parameters, the predicted average yearly sediment input, and its 68% confidence interval, was 4,376 [plus or minus] 75 [t*[yr.sup.-1]]. In 1982, 204,141 t sediment was dredged out of the reservoir (Gabriels 1985). The previous dredging dated from 1936. This means that in the 46 years between the two dredging operations 204,141 t sediment was deposited in the reservoir, which represents a yearly input of 4,438 [t*[yr.sup.1]]. This value lies within the 68% confidence interval; the model output and field amount differed by only 1.4%. It can be concluded that the RUSLE model, in combination with the methods presented in this paper, p articularly the adaptation of the linear RUSLE model to a three-dimensional reality, are capable of predicting both on-site soil losses (sediment leaving each field) and off-site sediment accumulation within acceptable accuracy. Note that the RUSLE model cannot be used to estimate intrabasinal sediment storage as colluvium.

Discussion and conclusions

Extending the RUSLE with an error propagation was very important in evaluating the difference between what can be predicted by a model and what exists in the field. If field measurements coincide reasonably with model output, then the difference between the model output and the field truth can be explained mainly by the uncertainty of the model input. Otherwise, the difference must be mainly due to the model itself. For example, in regions where gully erosion is prominent, the RUSLE will underestimate soil losses and sediment accumulation. Also, an incorrect assignment of the r-value of a field in Equation 4 can be responsible for considerable deviations between model output and field truth. However, it was not possible to estimate the amount of error induced by the model itself. This requires supplementary statistical information about the model regression equations (confidence and prediction interval equations) from the model designers.

Scheinost et al. (1997) and Sinowski et al. (1997) investigated the error contribution of a model (a pedotransfer function to predict soil water retention) and the error contribution of the spatial interpolation. Probably due to error self-compensation, the overall error was substantially smaller than the sum of both single components. The small deviation between measured and predicted sediment accumulation in this study may indicate a similar error compensation rather than the contribution of errors from the input data only.

Computing sources of error can be done by pretending that all parameters have no error, except for the parameters that are traced. Comparing these results with the original model output shows the (relative) error contribution of that particular parameter. Percentile error maps of each parameter can be calculated by dividing the variance by the total variance and multiplying by 100 (Spiegel and Meddis 1982; Wesseling and Heuvelink 1993). Figure 10 gives the percentile error maps for the LS, K, and C factors. In the case that the uncertainty induced by the data is too large, one can select the input parameters and locations that must be sampled more precisely or at a higher spatial resolution. The major error contribution is from the LS factor, indicating the need of high-quality DEM data in erosion studies.

This study indicates the possible power of the RUSLE model when applied in agricultural watersheds, and when used within the boundary limits of the model. Until now, there was no physical model capable of accurately predicting off-site sediment accumulation with the same amount of input data. Because the calculated and measured mean yearly sediment input in the water reservoir differed only by 1.4%, and gully erosion is not significant in the study area, the RUSLE LS equation proved to be, in this study, a reliable estimate of the sediment transport capacity of the overland sheet flow and small till flow.

Modeling environmental processes is very complex. No matter how complex the model, it always is a generalization of real processes, with output always the result of computations and operations on uncertain input data. Therefore, model calculations should be performed using probability distributions of the input parameters rather than their deterministic values.

Jan Biesemans is research scientist, Marc Van Meirvenne is professor of spatial information processing, and Donald Gabriels is professor of soil and water conservation for the Department of Soil Management and Soil Care at University Gent, Gent, Belgium. Their research is financed by the Flemish Institute for the Promotion of Scientific and Technological Research in the Industry (IWT) and the Province of West Flanders.

REFERENCES CITED

Costa-Cabral, M.C., and S.J. Burges. 1994. Digital elevation model networks (DEMON): A model of flow over hillslopes for computation of contributing and dispersal areas. Water Resources Research 30(6):1681-1692.

Desmet, P.J.J., and G. Govers. 1996. A GIS procedure for automatically calculating the USLE LS factor on topographically complex landscape units. Journal of Soil and Water Conservation 51(5):427-433.

Fairfield, J., and P. Leymarie. 1991. Drainage networks from digital elevation models. Water Resources Research 27(5):709-717.

Flanagan, D.C., and M.A. Nearing. 1995. USDA-Water Erosion Prediction Project (WEPP). West Lafayette, NSERL Report No. 10.

Foster, G.R., and W.H. Wischmeier. 1974. Evaluating irregular slopes for soil loss prediction. Transactions of the ASAE. 17:305-309.

Freeman, TG 1991. Calculating catchment area with divergent flow based on a regular grid. Comput. Geosci. 17(3):413-422.

Gabriels, D. 1985. Erosion and pollution (in Dutch). Water 4:139-144.

Ghekiere, G. 1997. Assessment of the actual erosion risk in the Kemmelbeek watershed (in Dutch). Master's thesis, University Gent: Belgium.

Heuvelink, G.B.M., P.A. Burrough, and A. Stein. 1989. Propagation of errors in spatial modelling with GIS. International Journal of GIS. 3:303-322.

Jenson, S.K., and J.O. Domingue. 1988. Extracting topographic structure from digital elevation data for geographic information system analysis. Photogramm. Eng. Remote Sens. 54(11):1593-1600.

Laflen, J.M., L.J. Lane, and G.R. Foster. 1991. WEPP, a new generation of erosion prediction technology. Journal of Soil and Water Conservation 46(1):34-38.

McCool, D.K., L.C. Brown, G.R. Foster, C.K. Mutchler, and L.D. Meyer. 1987. Revised slope steepness factor for the Universal Soil Loss Equation. Transactions of the ASAE. 30(5):1387-1396.

McCool, D.K., G.R. Foster, C.K. Mutchler, and L.D. Meyer. 1989. Revised slope length factor for the Universal Soil Loss Equation. Transactions of the ASAE. 32(5):1571-1576.

Molz, F.J., and H.H. Liu. 1997. Fractional Brownian motion and fractional Gaussian noise in subsurface hydrology: a review, presentation of fundamental properties, and extensions. Water Resources Research 33(10):2273-2286.

Moore, I.D., and J.P. Wilson. 1992. Length-slope factors for the Revised Universal Soil Loss Equation: simplified method of estimation. Journal of Soil and Water Conservation 47(5):423-428.

Nearing, M.A., L.D. Norton, D.A. Bulgakov, G.A. Larinov, L.T. West, and K.M. Dontsova. 1997. Hydraulics and erosion in eroding tills. Water Resources Research 33(4):865-876.

O'Callaghan, J., and D.M. Mark. 1984. The extraction of drainage networks from digital elevation data. Comput. Vision Graphics Image Process. 28:323-344.

Peitgen, H., H. Jurgens, and D. Saupe. 1992. Chaos and Fractals, new frontiers of science. Springer-Verlag: New York, NY.

Polidori, L. 1991. Description of terrain as a fractal surface, and application to digital elevation model quality assessment, Photogramm. Eng. Remote Sens. 57(10):1329-1332.

Quinn, P., K. Beven, P. Chevalier, and O. Planchon. 1991. The prediction of hillslope flow paths for distributed hydrological modelling using digital terrain models. Hydrol. Processes. 5:59-79.

Renard, K.G., G.R. Foster, G.A. Weesies, D.K. McCool, and D.C. Yoder. 1996. Predicting soil erosion by water: a guide to conservation planning with the Revised Universal Soil Loss Equation (RUSLE). USDA-ARS, Agriculture Handbook 703.

Scheinost, A.C., W. Sinowski, and K. Auerswald. 1997. Regionalirarion of soil water retention curves in a highly variable soilscape, I. Developing a new pedotransfer function. Geoderma 78:129-143.

Sinowski, W., A.C. Scheinost, and K. Auerswald. 1997. Regionalization of soil water retention curves in a highly variable soilscape, II. Comparison of regionalization procedures using a pedotransfer function. Geoderma 78:145-159.

Spiegel, M.R., and R. Meddis. 1982. Schaum's outline of theory and problems of Probability and Statistics, SI (metric) edition. McGraw-Hill: New York, NY.

Van Meirvenne, M. 1991. Characterization of soil spatial variation using geosratisties. Ph.D. thesis, University Gent: Belgium.

Webster, R., and M.A. Oliver. 1990. Statistical methods in soil and land resource survey, Spatial Information Systems. Oxford University Press: Oxford, MA.

Wesseling, C.G., and G.B.M. Heuvelink. 1993. ADAM, an error propagation tool for geographical information systems. University of Utrecht, Department of Physical Geography, second beta version, 52 p.

Wischmeier, W.H., and D.D. Smith. 1978. Predicting rainfall erosion losses: a guide to conservation planning. Agricultural Research Service: Washington, DC. Agriculture Handbook 537. U.S. Department of Agriculture.

Printer friendly Cite/link Email Feedback | |

Author: | Biesemans, J.; Meirvenne, M. Van; Gabriels, D. |
---|---|

Publication: | Journal of Soil and Water Conservation |

Geographic Code: | 1USA |

Date: | Jan 1, 2000 |

Words: | 5971 |

Previous Article: | REMM: The Riparian Ecosystem Management Model. |

Next Article: | Are green lots worth more than brown lots? An economic incentive for erosion control on residential developments. |

Topics: |