**5.2 K factor**

For K factor layer, a digital soil map (IBGE, 2001 and 2007) was considered. For determining the value of erodibility we used a soil profile database provided for whole Brazilian territory (Cooper et al., 2005). The erodibility was calculated indirectly through the method proposed by Boyoucos (1935), called clay-ratio method (Equation 3).

$$\text{Erodibility} = \left[ (\text{sand} + \text{silt}) \;/ \; (\text{clay}) \right] / 100 \tag{3}$$

Where: erodibility expressed in t.h.MJ-1.mm-1 and the textural attributes expressed in g.kg-1.

Driest areas occur in northeastern region, where annual rainfall amount is approximately 600 mm. Northern region encompasses areas whose annual rainfall amount and erosivity are expressive (annual rainfall amount is usually over 2,500 mm), for example Amazon region (Silva, 2004). In most of Brazilian territory, annual rainfall depth ranges from 1,000 to

Fig. 4. Annual rainfall amount map of Brazilian territory (in mm.y-1). Source: Silva (2004).

For each one NPE factor a single, separated layer was elaborated (Figure 5) and for each

For R factor layer, we used the digital map of rainfall erosivity (Silva, 2004). In this study, the author considered eight majors Brazilian regions covering the whole of the territory of Brazil, and for each region, one adapted equation was applied using pluviometric records obtained from 1,600 weather stations with continuous database of at least twenty

For K factor layer, a digital soil map (IBGE, 2001 and 2007) was considered. For determining the value of erodibility we used a soil profile database provided for whole Brazilian territory (Cooper et al., 2005). The erodibility was calculated indirectly through the method proposed

Where: erodibility expressed in t.h.MJ-1.mm-1 and the textural attributes expressed in g.kg-1.

Erodibility = [(sand + silt) / (clay)] / 100 (3)

**5. The natural potential for soil erosion map** 

factor, a specific database was elaborated, as described below.

by Boyoucos (1935), called clay-ratio method (Equation 3).

2,000 mm (Figure 4).

**5.1 R factor** 

consecutive years.

**5.2 K factor** 

Fig. 5. Soil loss prediction through USLE approach and using GIS, modified from Mongkolsawat et al. (1994).

This equation is based only in textural characteristics of the soil. Due this simplicity and feasibility of obtaining, it is still largely used (Waswa et al, 2002, Lopes-Assad et al., 2009). We also used a complementary database regarding soil erodibility that occurs specifically along São Paulo State (Silva & Alvares, 2005). Soil erodibility values were classified into five classes, as shown in Table 1 (Giboshi, 1999).


Source: Giboshi (1999).

Table 1. Values and interpretation classes for soil erodibility.

Official Brazilian soil map (IBGE, 2001) was crossed with the soil erodibility map, and thus we obtained average values of the K factor as large groups of Brazilian soils. For this, we used the tool "Zonal statistics as table" in GIS ArcMAP 10 (Theobald, 2007).

Natural Potential for Erosion for Brazilian Territory 11

The annual rainfall erosivity ranges from 3,116 to 20,035 MJ mm ha-1 h-1 year-1 (Silva, 2004). The region with the lowest values is represented by the northeastern region and an occurrence in southeastern region. Highest values are found in the northern region, mostly in the Amazon region (Figure 6). Predominant class was "> 12,000 MJ mm ha-1 h-1 y-1" with

Spatial distribution of the rainfall amount and erosivity are irregular (Figures 4 and 6). In some regions the annual erosivity normally is incipient and others are extremely high (almost ten times more erosive than the areas with lowest erosivity). Maps elaborated by Rao et al. (1996) and showed in Figure 7 confirm this information. Such maps show that the trimester with major or minor contribution over seasonal distribution of amount of rain

Fig. 6. Annual erosivity map (MJ mm ha-1 h-1 y-1). Source: Silva (2004) – reclassified.

Table 2. Percentage of occurrence of each class of the R factor along Brazilian territory.

Erosivity (MJ mm ha-1 h-1 y-1) % < 4,500 3.0 4,500-7,000 25.0 7,000-9,500 13.0 9,500-12,000 22.0 > 12,000 37.0

**6. Results and discussion** 

37.0% of occurrence (Table 2).

along Brazilian territory is also changeable.

**6.1 R factor** 

Experimental semivariograms were determined until approximately 50% of the geometric camp, since after this value the semivariogram did not seem correct (Guerra, 1988), i.e., its accuracy was reduced due to a smaller number of possible pairs to calculate the semivariance at this distance. A geometric camp of 16 degrees (geographic coordinates) with partition groups (lags) of 1 degree was considered, as these lags are the estimators of the experimental semivariograms (Deutsch & Journel, 1998). Theoretical models considered, such as spherical, exponential, Gaussian and linear, were described by Guerra (1988) and Andriotti (2003).

Only this theoretical semivariogram group was considered because it usually covers the general dispersion situation of soil science spatial events (Burrough & McDonnell, 1998; Soares, 2006). Correlation coefficient of selected models were obtained through Cross Validation routine of the geostatistical software GS+, version 9. The spatial dependence index (SDI) was used according to Zimback (2001), which measures a sample's structural variance effect on total variance (sill). SDI comprises the following interpretation break: weak SDI ≤ 25%, moderate SDI between 25% and 75% and strong SDI ≥ 75%. This index is a complement of the traditional method recommended by Cambardella et al. (1994) in which the nugget weight effect (randomness) on total variance is evaluated. Through structural parameters obtained from experimental semivariograms, maps of some properties were created using GIS ArcMap v.10 (ESRI, 2010). A punctual ordinary kriging estimator was used for geostatistic interpolation.

#### **5.3 Topographic factor**

For LS layer (computed jointly), the Digital Elevation Model for Brazilian territory was obtained from SRTM project (Shuttle Radar Topography Mission) (Farr & Kobrik, 2000), that it is in the fourth version (Jarvis et al., 2008). The LS map was generated through the algorithm available in Wischmeier and Smith (1978):

$$\text{LS} = (\lambda/22.1)^{\text{m}\_s} \left( 0.065 + 0.045 \,\theta + 0.0065 \,\theta^2 \right) \tag{4}$$

Where λ = slope length (m); θ = slope gradient (%); and m = 0.5 if the percent of slope is 5 or more, 0.4 on slopes of 3.5 to 4.5 percent, 0.3 on slopes of 1 to 3 percent, and 0.2 on uniform gradients of less than 1 percent.

The values of λ and θ were derived from DEM (ESRI, 2010). For determination of λ value we used the method proposed by Moore & Burch (1986).

$$
\lambda = \text{(Flow Acceleration}^\* \text{ Cell Size)} \tag{5}
$$

Where: Flow Accumulation is a grid theme of flow accumulation expressed as number of grid cells (readily derived from watershed delineation processing steps) and Cell Size is the length of a cell side (m). Flow Accumulation was derived from the DEM, after conducting Fill and Flow Direction processes in ArcGIS 10 (Theobald, 2007; ESRI, 2010).

#### **5.4 NPE map**

Using Equation 2 and approach shown at Figure 5, the NPE layer was created. The final map was reclassified into interpretative classes. We analyzed which feature(s) influence(s) more expressively the spatial variability of the values of NPE along the study area. Hence, we interpreted the map considering the possibilities of aggravation of erosive process by present land use patterns which can be changed.
