SciELO - Scientific Electronic Library Online

vol.66 issue5First-order decay models to describe soil C-CO2 Loss after rotary tillageSesbania virgata stimulates the occurrence of its microsymbiont in soils but does not inhibit microsymbionts of other species author indexsubject indexarticles search
Home Pagealphabetic serial listing  

Services on Demand




Related links


Scientia Agricola

On-line version ISSN 1678-992X

Sci. agric. (Piracicaba, Braz.) vol.66 no.5 Piracicaba Sept./Oct. 2009 



Soil erosion fragility assessment using an impact model and geographic information system


Avaliação da fragilidade à erosão do solo por meio de um modelo de impacto e sistema de informações geográficas



Luiz Alberto Blanco Jorge

UNESP/FCA - Depto. Recursos Naturais - C.P. 237 - 18603-970 - Botucatu, SP - Brasil - e-mail <>




A study was taken in a 1566 ha watershed situated in the Capivara River basin, municipality of Botucatu, São Paulo State, Brazil. This environment is fragile and can be subjected to different forms of negative impacts, among them soil erosion by water. The main objective of the research was to develop a methodology for the assessment of soil erosion fragility at the various different watershed positions, using the geographic information system ILWIS version 3.3 for Windows. An impact model was created to generate the soil's erosion fragility plan, based on four indicators of fragility to water erosion: land use and cover, slope, percentage of soil fine sand and accumulated water flow. Thematic plans were generated in a geographic information system (GIS) environment. First, all the variables, except land use and cover, were described by continuous numerical plans in a raster structure.The land use and cover plan was also represented by numerical values associated with the weights attributed to each class, starting from a pairwise comparison matrix and using the analytical hierarchy process. A final field check was done to record evidence of erosive processes in the areas indicated as presenting the highest levels of fragility, i.e., sites with steep slopes, high percentage of soil fine sand, tendency to accumulate surface water flow, and sites of pastureland. The methodology used in the environmental problems diagnosis of the study area can be employed at places with similar relief, soil and climatic conditions.

Key words: water flow, hydric erosion, analytical hierarchy process


A área de estudo, com 1566 ha, abrange uma microbacia hidrográfica, posicionada em local de relevo em contato reverso-frente de cuesta-depressão periférica paulista, inclusa na Bacia do Rio Capivara, no município de Botucatu, Estado de São Paulo. Esse ambiente é frágil, podendo sofrer diferentes formas de impactos negativos, dentre as quais a erosão hídrica do solo. O objetivo principal foi desenvolver uma metodologia, com auxílio do sistema de informações geográficas ILWIS v. 3.3 para Windows, para avaliar a fragilidade à erosão do solo em diferentes posições na microbacia. Um modelo de impacto foi criado para gerar o plano de fragilidade à erosão do solo, com base em quatro planos temáticos relacionados com indicadores de fragilidade a erosão hídrica: uso e cobertura do solo, declividade, percentagem de areia fina e acúmulo de fluxo de água. Em um primeiro passo, excetuando-se o uso e cobertura do solo, as outras três variáveis foram descritas por planos numéricos contínuos em estrutura raster. O plano de uso e cobertura do solo pôde ser representado também por valores numéricos, associados a pesos atribuídos a cada classe, a partir de uma matriz de comparação por pares, considerando-se o método de hierarquia analítica. Realizou-se uma checagem de campo final, em que foram registradas evidências de processos erosivos, naquelas áreas apontadas como apresentando maiores níveis de fragilidade, ou seja, locais com declividade acentuada, solo com alta percentagem de areia fina, tendência de acúmulo de fluxo superficial de água e ocupadas por pastagem. A metodologia desenvolvida pode ser usada no diagnóstico de problemas ambientais em locais com condições similares às da área de estudo.

Palavras-chave: fluxo de água, erosão hídrica, método de hierarquia analítica




Among the descriptive and interpretive aspects of an environmental diagnosis, Gómez-Orea (2002) includes the need to understand how the land and its natural resources are utilized, considering degradations, threats, and an estimate of the fragility or vulnerability for the development of human activities.

Water erosion assessment, which is part of the diagnosis of environmental problems, is highly relevant since inadequate land use can accelerate the processes of erosion and deposition that occur naturally, leading to modifications related to soil conservation, water production and quality, and environmental changes in certain locations of a drainage basin. Various approaches and equations for risk assessment or predictive evaluation of soil erosion by water are available in literature. Among then the universal soil loss equation (Wischmeier & Smith, 1965, 1978) and the revised universal soil loss equation (Renard et al., 1994; Yoder & Lown, 1995). Wu & Wang (2007) proposed a model to evaluate the risk index for soil erosion by water, in which the remote sensing, GIS, the analytical hierarchy process (AHP) and modeling techniques are integrated.

The present study was undertaken in a 1,566 ha watershed situated in a cuesta relief area of São Paulo State Peripheral Depression. This environment is extremely fragile and can be subjected to different forms of negative impacts, including soil erosion. The main objective was to develop a methodology, aided by the geographic information system ILWIS version 3.3 for Windows, to assess the soil's fragility to erosion by water at the different sites.



Study area

The watershed is located in the Capivara River basin in the municipality of Botucatu, São Paulo State, Brazil (Figure 1). It comprises: i - a small area at the cuesta reverse (beginning of the S. Paulo Eastern Plateau), with altitudes of 700 to 810 m; ii - a cuesta front (a sandstone and basaltic escarpment with their originated soils ); and iii - a peripheral depression segment, with altitudes from 465 to 600 m, comprising an area of sandstones and alluvial sediments encompassing the Capivara River wetlands. The soils of the Capivara River flat wetlands are Dystric Fluvisols, Mollic Gleysols and Dystric Gleysols. Toward the front of the cuesta, but still within the peripheral depression there are frequent gently undulating uprisings (2 to 20% slope) with Albic Arenosols, Acric Ferralsols, Rhodic Ferralsols, Chromic Luvisols, Dystric Nitisols and Haplic Chernozems. The basaltic cuesta front have Lithosols on 20% to 40% slopes. Rock's outcrops are frequent on scarped relief areas, where slopes exceed 40% (Carvalho et al., 1991).

The climate has two distinct seasons, one rainy and hot (September to March) and the other dry and cold (April to August). The native vegetation partially loses its leaves in the dry cold season. It is classified as Semidecidual Seasonal Forest (Veloso, 1992), with fragments of altered forest that have undergone various levels of anthropic disturbance. Fragments of natural vegetation, a transition of forest to forested savanna (cerradão) and gallery forest, are also present. Land is used mainly as pasture. Small areas are planted with forests and cultivated as farmland.

Digital elevation model and slope

A digital elevation model (DEM) is defined as a numerical structure of data that represents the spatial distribution of the altitude of the terrain's surface (Felicísimo, 1994). The DEMs were integrated to the GIS, offering a series of applications and different options of usage. One of these options is the generation of new spatial data based on DEMs, such as slope, aspect, effects of illumination, determination of the water route, etc.

The study area database was structured (from 2005 to 2007) in the geographic information system ILWIS for Windows (International Institute for Geo-Information Science and Earth Observation, 2001) environment. The two study area planialtimetric maps were digitized (using a scanner) and georeferenced. The maps were used as the basis for the spatial data input operation. The watershed's boundaries, the contour lines, the summit elevations, and the drainage network were digitized. After digitizing the contour lines, the digital elevation model of the watershed and the thematic slope plan were created in the GIS environment.

Orthophotos and land use plan

Orthophotos were generated based on aerial photographs taken in the 2000 and 2005 years. These orthophotos served as the basis for outlining the land use and cover plan. The procedure to produce the orthophotos can be subdivided into the following steps: (i) acquisition of the aerial photography in digital form; (ii) creation of a georeference; (iii) specification of the system of coordinates and of the file of the digital elevation model utilized; (iv) recording of the fiducial markers, which consists of the specification of the focal distance of the camera and of the photo coordinates of the fiducial markers; (v) indication of the control points detected in the aerial photography and on a map (Secretaria de Economia e Planejamento do Estado de São Paulo, scale 1:10000, UTM projection).

Water flow direction and flow accumulation

The planes of water flow direction and flow accumulation were generated in the ILWIS environment, based on the digital elevation model (DEM). The latter plan is related to the concept of contribution area. Any watershed may contain a large variety of geological, topographical, soil, vegetation and land use features that will affect the relation between precipitation and flow. There are many points, which show a similar hydrologic behavior, with similar water balance and runoff characteristics. The hydrologic similarity index was developed as part of the rainfall-runoff model named the Topmodel (Beven, 2001). One of the basic premises of Topmodel has to do with the dynamics of saturated regions of the watershed that can be approximated considering successive representations of the saturated zone of an area "a" drained to a point of the slope (contribution area).

Soil attribute

Soil sampling points were defined starting from a regular grid allocated so that the interval between points in the X and Y directions was 500 m, using the navigation capacity of the GPS to reach each point. Sixty-three points were sampled. At each point, a sample was taken at the 0 to 20 cm soil depth. The soil samples were analyzed to obtain information about the granulometry, with emphasis on the percentage of fine sand.

Attributing the soil sampling points to their respective plan coordinates obtained by GPS, a connection was drawn between the location of these points and the variable of the soil, whose results were tabulated into a table of attributes. The spatial correlation module of ILWIS provided experimental semivariogram about the spatial behavior of the variable of the soil property. A semivariogram is a variance graph of paired sample measures. It provides a means to quantify the commonly grouped samples observed tendency to have more approximate values than more widely separated samples. The semivariogram described by Isaaks & Srivastava (1989), which provides the spatial dependence, is expressed by:

where: γ(h) = semivariance for the distance h; xi and xi + h = sampling sites separated by the distance h; Z(xi) and Z(xi + h) = measured values of the variable.

In addition to the experimental semivariogram, the spatial correlation module of ILWIS also allows one to obtain the values of Moran's I (1948) and Geary's C (1954) statistics, which are classical methods to evaluate the existence or absence of spatial autocorrelation. These methods allow one to verify if the values of a variable at a given site are independent of the values of the variable at neighboring sites. A positive spatial autocorrelation refers to a pattern in which points with similar values tend to cluster, while a negative spatial autocorrelation indicates a pattern in which points with similar values are scattered throughout the map. When there is not a statistically significant spatial autocorrelation, the spatial distribution pattern is considered random. The interpretation of the two statistics can be summarized as: 0 < C < 1 and I > 0 strong positive autocorrelation; C > 1 and I < 0 strong negative autocorrelation; C = 1 and I = 0 random distributions of the values.

The power, exponential and Gaussian models were tested to adjust the estimated semivariogram in the GIS environment. The appropriate model was selected to describe the semivariance of the variable as a function of the coefficient of determination (R2), of the sum of squares of residues (SSR) and of the mean percentage deviations (MPD). Entering the estimators of the selected model into the ILWIS kriging module, using the ordinary kriging method, the values obtained at the sampled points were interpolated, deriving the corresponding continuous plan in a raster structure.

Soil erosion fragility

The soil erosion fragility layer was drawn from the four thematic plans generated in the geographic information system environment: land use and cover, slope, percentage of fine sand, and water flow accumulation. The landscape elements were chosen as indicators of soil erosion fragility based on the field survey and information analysis, combined with professional expert judgment. In the first step, all the variables except land use and cover were described by continuous numerical plans in a raster structure.

The land use and cover plan was also represented as numerical values associated to weights attributed to each class, based on a pairwise comparison matrix, considering the analytical hierarchy process-AHP. The significances of land use and cover classes for fragility to soil erosion by water were compared mutually and dissimilarities between their levels of importance were quantified, by the numbers 1 to 9, to indicate escalated differences, while reciprocals of those numbers were specified to their comparisons with an inverse order. Specifically, number 1 means that levels of importance for both classes are equal while number 9 represents that the level of importance for the first one is extremely higher than for the second one, according to the methodology presented by Saaty (1977). Using the AHP, the consistency of the allocation of weights to the classes of criteria involved was evaluated. The latent root and the consistency index (CI) (CI = latent root - n) / (n - 1); n = number of lines / columns in the comparison matrix) were calculated. One then obtains the consistency ratio (CR) by dividing the consistency index by the random index (RI), which has a fixed value (Bantayan & Bishop, 1998).

An impact model, considering the first impact layer (land use weights) and the second impact layer (w21 * (slope * fine sand) + w22 * flow accumulation), was proposed to generate the erosion fragility map:

where: fragility = soil erosion fragility layer; impact1 = first impact layer; impact2 = second impact layer; w1, w2, w21, w22 = assigned weights.

The highest values obtained after the layer was normalized indicate points of greater fragility, i.e., points where the land cover has a smaller protection capacity and the impact is greater (areas of higher slope grade and percentage of fine sand, greater concentration of water flow during rainfall). The layer was classified into five levels of fragility.

The Cohen's kappa coefficient, used as a measure of model accuracy and validation, was derived from the cross-tabulation of the mapped levels of soil erosion fragility against those observed in the 60 ground control points. The kappa statistic is expressed by:

where: Po = observed agreement; Pe = expected agreement;

Kappa is a measure of the difference, standardized to lie on a -1 to 1 scale, where 1 is perfect agreement, 0 is exactly what would be expected by chance, and negative values indicate agreement less than chance. Interpretation of the statistic according to Landis & Koch (1977) is shown in Table 1.




The sites with the highest slope values of the watershed are located at the cuesta front (Figure 2). Seventeen land use and cover categories (Table 2) were observed in the study area. The values (weights) attributed to the different classes of land use and cover can be considered satisfactory, since the calculated consistency ratio (CR) is lower than 0.10 (Table 3, underneath). The highest weights correspond to the land use classes that offer the least protection against soil erosion by water (Table 3 and Figure 3).





In the ILWIS environment, the DEM allowed for the derivation of the layer of accumulated water flow (contribution area). The Neperian logarithm was applied to the numerical plan in an attempt to reduce the wide range of values related to the contribution areas on each pixel of the raster map. For purposes of presentation, the values were normalized (amplitude of 0 to 1), and the division into 5 classes of accumulated flow enabled us to produce the thematic map (Figure 4).



Entering the percentage of fine sand values (Figure 5) into the ILWIS spatial correlation module resulted in the experimental semivariogram from the spatial behavior of the variables. The theoretical semivariogram was also adjusted (Table 4). The semivariogram of the variable percentage of fine sand did not show stabilization of the semivariance, considering the maximum distances between pairs of sampled points. The power model adjusted adequately to experimental values, which can be ascertained by analyzing the coefficient of determination - R2, the sum of squares of residues - SSR and the mean percentage deviations - MPD (Table 4).



In addition to the Table 5 data the values of the statistics I of Moran (1948) and Geary's C (1954) were obtained, which allow one to evaluate the existence of autocorrelation or random distribution of the variable. An analysis of the values of I and C (Table 5) reveals that the variable percent of fine sand (0 - 20 cm depth) shows a positive autocorrelation for pairs of values that are up to 2,000 meters distant (pattern in which points of similar values tend to cluster). There is a transition band from positive to negative autocorrelation between 2,000 and 2,500 m, and the negative autocorrelation for pairs of values between 2,500 and 4,500 m (pattern in which points of similar values are scattered over the map).

After adjustment of the estimated semivariogram by entering the coefficients of the model into the GIS kriging model, the continuous plan was derived for the variable percentage of fine sand in the 0 to 20 cm depth (Figure 6 and 7).





The four thematic plans related to water erosion fragility landscape elements were used to produce the erosion fragility map, according to the proposed modeling approach (Figure 8). The weights were assigned in the model after the field surveys. The impact layer was classified into five levels of fragility (Figure 9). The kappa coefficient value of 0.61 indicated a substantial agreement between the model and the ground control points.





Fragility levels 3, 4 and 5 allow one to detect the areas in the landscape where the phenomenon of erosion is the most intense. In the final check after generating the soil erosion fragility map, the main erosive processes occurring in the watershed under study were recorded using photographs taken at each site. The photographs (Figure 9) is a sample of the comparison between the highest levels of fragility shown in the thematic plan and the reality found in the field.

Photograph 1 (Figure 9) depicts a linear erosive process (fragility levels 3 and 4) found within a forest fragment that regenerated from a pastureland. Photograph 2 shows gully erosion in an area (level 3 of fragility in the surroundings and levels 4 and 5 at the sites of the erosive process) situated at the head of a drainage channel, in an area of pasture with steep slope and a soil with high percentage of fine sand. Photograph 3 depicts a situation (level 3 of fragility) occurring at a site that encompasses the three aforementioned landscape elements: the land use is pasture, the slope is steep, and there is a high percentage of fine sand. The methodology employed to evaluate the erosion fragility can be used for the detection of processes that have already occurred and/or are evolving, but it also seeks to identify areas with a potential for soil loss.



The final check, performed in the field to record evidence of erosive processes, indicated a good agreement with the levels of fragility to water erosion as calculated for the thematic plan. Therefore, the methodology used in the environmental problems' diagnosis of the study area can be employed at places with similar relief, soil and climatic conditions.



BANTAYAN, N.C.; BISHOP, I.D. Linking objective and subjective modelling for landuse decision-making. Landscape and Urban Planning, v.43, p.35-48, 1998.         [ Links ]

BEVEN, K.J. Rainfall: runoff modeling; the primer. Chichester: John Wiley, 2001. 360p.         [ Links ]

CARVALHO, W.A.; PANOSO, L.A.; MORAES, M.H. Levantamento semidetalhado dos solos da Fazenda Experimental Edgardia - Município de Botucatu. Botucatu: UNESP/FCA 1991. 467p. (Boletim Científico, v.2).         [ Links ]

FELICÍSIMO, A. Modelos digitales del terreno: introducción y aplicaciones en las ciencias ambientales. Oviedo: Pentalfa, 1994. 118p.         [ Links ]

GEARY, R.C. The contiguity ratio and statistical mapping. The Incorporated Statistician, v.5, p.115-145, 1954.         [ Links ]

GÓMEZ-OREA, D. Ordenación territorial. Madrid: Mundi-Prensa, 2002. 704 p.         [ Links ]

INTERNATIONAL INSTITUTE FOR GEO-INFORMATION SCIENCE AND EARTH OBSERVATION - ITC. ILWIS 3.0 for Windows: user´s guide. Enschede: ITC, 2001. 530 p.         [ Links ]

ISAAKS, E.H.; SRIVASTAVA, R.M. An introduction to applied geoestatistics. New York: Oxford University Press, 1989. 560p.         [ Links ]

LANDIS, J.R; KOCH, G.G. The measurement of observer agreement for categorical data. Biometrics, v.33, p.159-174, 1977.         [ Links ]

MORAN, P.A.P. The interpretation of statistical maps. Journal of the Royal Statistical Society. Series B, v.37, p.243-251, 1948.         [ Links ]

RENARD, K.G.; FOSTER, G.R.; YODER, D.C.; McCOOL, D.K. RUSLE revisited: status, question, answers, and the future. Journal of Soil and Water Conservation, v.49, p.213-220. 1994.         [ Links ]

SAATY, T.L. A scaling method for priorities in hierarchical structures. Journal of Mathematical Psychology, v.15, p.234-281. 1977.         [ Links ]

VELOSO, H.P. Manual técnico da vegetação brasileira. Rio de Janeiro: IBGE/Departamento de Recursos Naturais e Estudos Ambientais, 1992.         [ Links ]

WISCHMEIER, W.H.; SMITH, D.D. Predicting rainfall erosion losses from cropland east of the Rocky Mountains. Washington, DC: USDA/ARS, 1965. 47p. (Agricultural Handbook, n.282).         [ Links ]

WISCHMEIER, W.H.; SMITH, D.D. Predicting rainfall erosion losses: a guide to conservation planning with the universal soil loss equation (USLE). Washington, DC: USDA, 1978. 58p. (Agricultural Handbook, 537).         [ Links ]

WU, Q.; WANG, M. A framework for risk assessment on soil erosion by water using an integrated and systematic approach. Journal of Hydrology, v.337, p.11-21. 2007.         [ Links ]

YODER, D.; LOWN, J. The future of RUSLE: inside the new revised universal soil loss equation. Journal of Soil and Water Conservation, v.50, p.484-489, 1995.         [ Links ]



Received January 18, 2008
Accepted February 16, 2009

Creative Commons License All the contents of this journal, except where otherwise noted, is licensed under a Creative Commons Attribution License