2.3. Annual soil erosion modeling and mapping (2005 - 2015)
The inherent soil erosion potential depends on rainfall, soil type,
LCLU, and terrain characteristics, which are represented by various
factors, such as, rainfall erosivity, soil erodibility, slope length and
slope steepness, LCLU management, and conservation practices. At the
national scale in Pakistan, in this study RUSLE (equation 1) (Renardet al. , 1991) was used to estimate soil erosion at 1 km spatial
resolution for 2005 and 2015.
\(A\ =\ R*K*L*S*C*P\) (1)
where, A = soil erosion (ton/ha/year), R = rainfall erosivity factor (MJ
mm/ha/hr/year), K = soil erodibility (ton·ha·hr/ha/MJ/mm), L =
slope-length factor (dimensionless), S = slope-steepness factor
(dimensionless), C = cover management factor (dimensionless), and P =
conservation practice factor (dimensionless).
Equation 2 (Renard & Freimund, 1994) was used to calculate rainfall
erosivity factor (R) for 2005 and 2015 using annual mean precipitation
data of the respective year. Over the diverse geographical and
topographical study areas, several researchers utilized equation 2 to
calculate rainfall erosivity factor (Howland et al. , 2018; Lamyaaet al. , 2018; Uddin et al. , 2016, 2018).
\(R\ =\ 0.0483P\hat{}1.610\) (2)
where, R = rainfall erosivity factor in Megajoule millimeters per
hectare per hour per year (MJ mm/ha/hr/year) and P = annual
precipitation in millimeters (mm).
In the FAO soil map, fourteen soil classes were found which behave
differently towards soil erosion and hence contain varying soil
erodibility (K) factor values that explain the degree of soil
erodibility of each soil type. Therefore, for this study, the K factor
values were obtained from published research articles, which were
assigned to fourteen soil classes to obtain a soil erodibility factor
map, used in the RUSLE (Table 1 and Table S1).
In this study, we used equations 3 and 4 to calculate the slope-length
(L) and slope-steepness (S) factors respectively which are adopted from
variously published studies in the South Asia region (Uddin et
al. , 2016, 2018) The combined LS-factor describes the effect of
topography on soil erosion. As temporal topography data was unavailable
for 2005 and 2015, a one-time computed L & S spatial layers were used
for annual soil erosion estimations and mapping for 2005 and 2015 over
the entire Pakistan.
\(L\ =([\lambda/22.13]\hat{}m)\) (3)
where, L = slope-length factor (dimensionless), λ =grid size (1 km x 1
km) or field slope length, and m = 0.5 for slopes > 4 %;
0.4 for 4% slope; 0.3 for slopes < 3% (Wischmeier & Smith,
1978)
\(S=\ ((0.43\ +\ 0.30\ s\ +\ 0.043\ s\hat{}2)/6.613)\) (4)
Where, S = slope-steepness factor (dimensionless), s = slope derived
from DEM in percentage (%)
In this study, MODIS’s seventeen LCLU classes defined according to the
International Geosphere-Biosphere Programme (IGBP) classification scheme
were clamped/recoded to six LCLU classes (‘Cropland’, ‘Forest land’,
‘Grassland’, ‘Other land’, ‘Settlements’ and ‘Wetlands/Snow cover’)
which were defined by Intergovernmental Panel on Climate Change (IPCC)
(Table S2). From several published articles on RUSLE, we synthesized
average cover management factor (C) values for each IPCC LCLU class,
which were assigned to 2005 and 2015 LCLU maps (Table 2 and Table S3).
In this study to determine the soil cover management (C) and
conservation practice factor (P) values, we compiled and used average
values from previously published studies in Pakistan, India, China, and
on a global scale (Table S4). Each LCLU class was assigned an average C
and P factor values given in Table 2, to get cover and use management
and conservation practice factor (P) maps for 2005 and 2015.
2.4. Quantification of soil erosion factors andannual rate of soil
erosion changes (2005 - 2015)
To generate soil erosion maps for 2005 and 2015, the spatially computed
soil erosion factors: ‘rainfall erosivity (R)’, ‘soil erodibility (K)’,
‘slope-length (L)’, ‘slope-steepness (S)’, ‘cover management (C)’ and
‘conservation practice (P)’ were inserted in equation 1.
The final annual soil erosion (ton/ha/year) maps were reclassified and
renamed into four classes; Low (<1), Medium (1–5), High
(5-20), and Very high (>20). In this study, each factor and
final annual soil erosion maps were spatially and statistically compared
between 2005 and 2015. For the annual soil erosion change assessment, we
used a cross tabular or change matrix method in which, diagonal values
show the stability of soil classes, while omission and commission values
indicate a shift in area or percentage between the four soil erosion
classes (adopted from Gilani et al., (2015) used for LCLUC assessment).
In this study, using the change matrix method, soil ‘loss’ (commission)
and gain (omission) statistically and spatially were reported between
2005 and 2015.
Latitude wise Pakistan extends from 24° N to 37.5° N. In this study
across latitude, mean soil erosion (ton/ha/year) in 2005 and 2015 were
plotted to understand eroded intensity linking to spatial distribution
of soil erosion areas.