Estimation of Spatial and Seasonal Variability of Soil Erosion in a Cold Arid River Basin in Hindu Kush Mountainous Region Using Remote Sensing

An approach is proposed in the present study to estimate the soil erosion in data-scarce Kokcha subbasin in Afghanistan. The Revised Universal Soil Loss Equation (RUSLE) model is used to estimate soil erosion. The satellite-based data are used to obtain the RUSLE factors. The results show that the slight (71.34%) and moderate (25.46%) erosion are dominated in the basin. In contrast, the high erosion (0.01%) is insignificant in the study area. The highest amount of erosion is observed in Rangeland (52.2%) followed by rainfed agriculture (15.1%) and barren land (9.8%) while a little or no erosion is found in areas with fruit trees, forest and shrubs, and irrigated agriculture land. The highest soil erosion was observed in summer (June–August) due to snow melting from high mountains. The spatial distribution of soil erosion revealed higher risk in foothills and degraded lands. It is expected that the methodology presented in this study for estimation of spatial and seasonal variability soil erosion in a remote mountainous river basin can be replicated in other similar regions for management of soil, agriculture, and water resources.


Introduction
Soil erosion harms agriculture productivity due to the reduction of soil fertility [1][2][3][4] Besides, it affects water quality by reducing the level of dissolved oxygen [5][6][7][8]. Soil erosion is caused by different drivers including water and air [9]. A large portion of global land (1094 Mha) is altered by soil erosion due to hydrological processes, of which a significant portion (224 Mha) cannot be restored [10]. The global average of soil erosion by either precipitation or floods is 20-30 Gtyr-1 [11].
Soil erosion rates are higher in arid and semi-arid regions due to high sensitivity of such region to climate and land-use changes [12][13][14][15][16][17]. The soil erosion rates are also reported to be higher in Asian countries like China, India, and Afghanistan [18][19][20]. Therefore, the arid region of Asia is considered as highly susceptible to soil erosion. Soil erosion is a slow natural geologic phenomenon, but it costs billions of dollars for agricultural losses [21][22][23]. Hence, soil erosion study is a vital concern for the regions with agriculturalbased economies.
Afghanistan located between the south and central Asia has a predominantly arid and semi-arid climate. Soil erosion is a challenging issue for agricultural productivity in Afghanistan [20,24]. Shrestha [25] reported that it is complicated to figure out the exact amount of soil erosion in Afghanistan due to a lack of data. The findings from the Global Assessment of Human-induced Soil Degradation (GLASOD) project indicated that a large portion (75%) of the topsoil of Afghanistan is removed due to water and wind. Furthermore, excessive erosions on riverbanks lead to notable social and economic losses. Despite the high importance of the soil erosion problem in Afghanistan, a few studies have been conducted to address the issue over the different regions in Afghanistan [20,26,27].
Sahaar [27] integrated the Revised Universal Soil Loss Equation (RUSLE) with Geographic Information System (GIS) to estimate the annual soil loss of Kabul Rivers. They found that the river was significantly eroded in the rate of 19 ton·acre −1 ·year −1 . [26] evaluated the annual soil loss of the Lower Harirud watershed in Heart province using RUSLE. The research findings indicated that soil erosion was varied in the range of 0.025 to 778 Mg ha −1 ·y −1 . [20] assessed the volume of eroded soil from Gardez basin in Paktya province, Afghanistan using the Universal Soil Loss Equation (USLE) integrated with GIS. Results showed that annual soil loss was found to be in the range of 0 to more than 100 ton·ha −1 ·y −1 . The findings reported in the literature indicated that soil erosion rates are estimated to be significantly high. Thus, further studies should be conducted over the vulnerable basins to provide sustainability in agriculture and environment in the country.
Kokcha subbasin, which is located in Amu Darya River basin in Northeast of Afghanistan, is a data-scarce region with high soil erosion rate. People with agricultural income are dominated in the Kokcha subbasin. Hence, their properties are significantly dependent on the farming industry [28][29][30]. Rapid deforestation and highly frequent rainfall depth, which change the river flow pattern, provide an increasing trend in soil erosion over the Kokcha subbasin [20,[31][32][33]. Hence, it is essential to estimate the potential amount of soil erosion in the basin and provide the measures required for protecting the vulnerable lands.
A number of approaches have been developed to estimate the soil erosion such as USLE, RUSLE, Agricultural Non-Point Source (AGNPS) model, and Watershed Erosion Prediction Project (WEPP) [34]. Several factors, including precipitation intensity, soil type, physiography, land use, and anthropogenic activities, have significant impacts on soil erosion [35][36][37][38][39][40]. Hence, soil erosion modelling needs a large number of in situ measurements including hourly rainfall, land use, and digital elevation model data which are not available for the basin [25,31,33]. Satellite-based remote sensing data are widely used to assess the susceptible zones of soil erosion in data-scare regions [34,41,42]. To estimate soil erosion with remote sensing data, the RUSLE method provided better performance [43][44][45]. The parameters of the RUSLE method can be estimated using the remote sensing data collected from several sources and augmented in geographical information system (GIS). Bahrawi et al. [46] reported that remote sensing data and GIS provided a flexible and powerful tool for assessing soil erosions in data-scarce regions.
This study aims to estimate the soil erosion using the satellite-based remote sensing data in the Kokcha subbasin. In this way, the RUSLE method is employed for estimating soil erosion. The satellite-based products are used to obtain the data required for computing the RUSLE factors due to scarcity of in-situ data in Afghanistan. Cold agricultural regions share a significant portion of global food production [47]. Soil erosion poses a major threat to soil degradation and crop productivity in cold regions. Wu et al. [48] found soil erosion in cold climate is often much higher than that in other climatic regions. Snow melting at the end of winter or early spring causes high overland flow which often exceed the runoff in rainy season of other climatic zones. Reduction of soil infiltration capacity due to freezing enhance surface runoff and thus erosion [49,50]. However, studies related to soil erosion in cold climate, particularly cold mountainous climate region is very limited.
Analysis of seasonal variability and spatial distribution of soil erosion in Kokcha subbasin can contribute valuable knowledge to soil erosion in cold climatic region. The soil erosion maps generated in this study would help to understand the spatial patterns of soil erosion within the Kokcha subbasin to assist stakeholders in taking suitable mitigation measures.

Study Area
Kokcha subbasin, which is located in the north-eastern of Afghanistan, covers an area of 21,811 km 2 ( Figure 1). The basin is situated from 35.40 • to 37.40 • N and 68.90 • to 71.60 • E. Additionally, the altitude of the basin is varied from 370 m above sea level in the northwestern part to 6721 m in the southern mountains. The percentage of the basin area above 2000, 3000, and 4000 m are respectively, 75%, 50%, 30%. The area receives an annual average rainfall of 470 mm [30]. The snow covers in the basin are approximately found above 4000 m. The mountainous areas in the basin are furnished with low-density trees. The basin is characterized by cold semi-arid and warm continental climate. Rainfall events mainly occur in spring (March-May) and autumn (October-April) seasons while winter (December-February) is usually cold with heavy snowfalls. Several small glaciers are observed in the Kokcha Basin. The glaciers and snow are the primary water resources to recharge the river in the summer (June-August) [51]. Soil in most parts of the basin can be categorized as Fluvisol based on Food and Agriculture Organization (FAO) classification.

Data Sources
Most of the hydrological and climatological stations in Afghanistan are destroyed due to four decades of internal war [30,31]. Hence, adequate observation data are rarely found for a national scale study. Thus, the satellite-based data, which are freely available for large areas, are preferred. The digital elevation model (DEM) of the basin with 30 m resolution is obtained from Advanced Space Thermal Emission Radiometer (ASTER). The Global Precipitation Measurement (GPM) product (GPM_3IMERGHH, half-hourly, 0.1 • × 0.1 • ) obtained by National Aeronautics and Space Administration (NASA) is used to provide the rainfall erosivity factor (R) data. The Normalized Difference Vegetation Index (NDVI) provided by Moderate Resolution Imaging Spectroradiometer (MODIS) is employed to extract the cover and management factor (C). Slope length and steepness (LS) factors are derived from DEM of Afghanistan, while soil erodibility factor (K) is obtained from the Soil Types Map from ISRIC World Soil Information. All the data are obtained at a spatial resolution of 0.1 • × 0.1 • .

Methodology
The RUSLE method has been widely used for predicting the rate of inter-rill and rill erosions [52]. The RUSLE integrated with GIS is used in the present study to predict soil erosion in the Kokcha subbasin. The RUSLE factors are calculated and mapped using remote sensing data and GIS techniques. The RUSLE method is defined as follows: where A is the potential annual soil erosion (t ha −1 year −1 ), R is rainfall/runoff erosivity factor, K is soil erodibility factor, LS is the slope length/steep factor, C is vegetation cover factor, and P is support practice factor.

Rainfall Erosivity Factor
The rainfall erosivity factor (R) reflects the effect of rainfall intensity on soil erosion [38,53]. A long period rainfall time series is required to compute the accurate R factor; however, that is not possible for many regions like Kokcha subbasin. To solve this limitation, the halfhourly GPM product is used to provide the rainfall data required in the present study.
The formula proposed by Wischmeier and Smith [53] is used to compute the R factor as the following equation: where R denotes the mean annual rainfall erosivity (100 foot-ton in acre −1 h −1 year −1 ), m indicates the number of storms in the N-year period, I i is the mean intensity of the i-th storm (in h −1 ), I30 i is maximum 30-min intensity of the i-th storm (in h −1 ).

Soil Erodibility Factor
The soil erodibility factor (K) quantifies the soil erosion caused by runoff [54]. Several parameters, namely soil texture, organic content, permeability and soil type, have significant impact on K factor [55]. Herein, the Erosion Productivity Impact Calculator (EPIC) is used to compute the K factor (100 acre) −1 ft −1 Mg −1 in −1 as follows [56]: where Cl, Sa and Si are the fraction (%) of clay, sand and silt respectively; and C is the soil organic carbon content (%).

Slope Length and Steepness Factor
The impacts of the slope length and slope gradient on soil erosion is quantified using the slope LS factor [54]. The accuracy of the LS value depends on the spatial resolution of the DEM, availability of computing resources, and the scale of the study area. DEMs with spatial resolution coarser than 100 m do not accurately capture the flow network of a catchment [44]. Hence, the LS factor is computed based on a DEM with a resolution of 30 m using the following formula: where Q a is flow accumulation grid; S is grid slope in percentage; M is grid size (x × y), and y is a dimensionless exponent in the range of 0.2 to 0.5.

Land Cover Management Factor
The Land Cover Management Factor (C) is used to express the effect of cropping and management practices on soil erosion [52]. The factor is most sensitive to human activities [57]. The Normalized Difference Vegetation Index (NDVI) is used to calculate the C factor in intraseasonal and interannual timescales [58]. The C factor can be computed as the following expression [59,60]: where α and β are the dimensionless parameters. van der Knijff et al. [60] proposed the best values of α and β were 2 and 1, respectively.

Support Practice Factor
Conservation or Support Practice (P) factor is quantified the practices that affect the runoff rate and thus, the soil erosion rate [52]. It denotes the ratio of soil erosion provided after a practice (i.e., contour farming, cross-slope cultivation, and strip cropping) to the regular farming slope [52]. The P-factor is considered 1 due to a lack of proper data of conservation and support activities in the basin.

Validation of Soil Erosion Model
The suspended sediment loads estimated by Hydrology Department of Ministry of Energy and Water of Afghanistan at Khawajaghar station, which is located at the basin outlet, are used as an observed dataset for the validation of RUSLE estimated soil erosion. The monthly erosion pattern of observed and modelled soil erosion is compared to show the ability of RULSE model to estimate the seasonal variability of observed soil erosion.

Spatial Pattern of the RUSLE Factors
To estimate the soil erosion in the Kokcha subbasin, the maps of RUSLE factors (e.g., R, K, LS, and C) are prepared using the ArcMap software ( Figure 2). The map of R factor is provided based on the mean annual rainfall data (Figure 2a). It is characterized based on the five classes in the range of 174 to 317 (MJ mm ha −1 ·h −1 ·y −1 ). The highest erosivity is observed in the southern of the basin while the northwestern region exposes the lowest values. The central parts of the basin represent the rainfall erosivity from 244 to 266 MJ mm ha −1 ·h −1 ·y −1 . A large area of the basin from north to south attains the erosivity in the range of 266 to 285 MJ mm ha −1 ·h −1 ·y −1 .

Spatial Pattern of Soil Erosion
Having provided the RUSLE factors maps, the annual soil erosion pattern in Kokcha can be obtained using Equation (1) through the raster calculator of spatial analyst tool in ArcMap software (Figure 3). Jenk's optimum classification algorithm (JOCA) [61] is used for classification of estimated soil erosion values to prepare the map. Based on JOCA, the soil erosion (t ha −1 ·y −1 ) in basin was classified into five classes including slight [0-10), moderate , high , severe [150-500), and very severe (>500). From Figure 3, it is evident that slight erosion is significantly observed over the study area, while moderate erosion is mainly found in the center area of the basin. The other erosion classes, including high, severe, and very severe, are rarely observed over the Kokcha subbasin.  Table 1 demonstrates the percentage of soil erosion in the classes defined in the present study. As shown in the table, it is clear that the slight erosion (0-10 t·ha −1 ·y −1 ) is found in 71.34%. The moderate erosion (10-50 t·ha −1 ·y −1 ), which has the second rank in vastity, covers 25.46% of the basin. The other erosion classes include 3.2% of the Kokcha subbasin in the mountains/foothills and degraded lands. The high slope steepness, very poor vegetation and no conservation practices are the most prominent causes of soil erosion over the Kokcha subbasin. There is no universal rule for classification of soil erosion values to prepare erosion map. Therefore, different classification schemes were adopted in different studies. Classification of soil erosion class is subjective which can significantly affect the areas covered by different soil erosion zones. JOCA can be used to avoid subjectivity in classification of soil erosion. JOCA classifies values in such a way that the variability of the values within a class is least but among the class is maximum. This ensures the best arrangement of erosion values in different classes. JOCA has been widely used for classification of environmental data for which any generalized classification scheme is not available [62].

Validation of the RUSLE Model
The available observed soil erosion data ( Table 2) is used to validate the results. Table 2 shows that 26 suspended sediment samples for different months are available for the period of 2014-2017.  Figure 4 represents the monthly pattern of soil erosion over the Kokcha subbasin. Figure shows the highest erosion rate in summer (June, July and August) due to snow melting from high mountains. This finding has good agreement with those obtained from the Khawajaghar station where the maximum sediment loads are observed in June and July months. The total annual soil loss estimated for the study area (22,116,974 t) is slightly higher than the mean annual sediment load observed at Khawajaghar station (17,847,711 t). The RUSLE model provides either the underestimated soil erosions due to ignoring the gully erosion, channel erosion, bank erosion, and mass wasting events such as landslides or the overestimated soil losses due to missing the deposition and sediment routing [44]. Hence, the difference between the measured and estimated erosion is acceptable.

Assessing the Soil Erosion over the Different Land Use and Land Cover
To assess the soil erosion over different land use classes of the Kokcha subbasin, the land use/land cover map of the study area is intersected with the soil erosion map. Table 3 demonstrates the percentage of soil erosion for different land use types. It is clear from the table that Rangeland (52.2%) followed by rainfed agriculture (15.1%) and barren land (9.8%) provides the highest percentage of erosion. The erosion in snow-covered area was found 3.4%. The erosion of snow-covered area might be due to an increase in temperature over Afghanistan [31]. Additionally, slight erosion is found in regions with fruit trees, forest and shrubs, and irrigated agriculture land.  Table 4 illustrates the soil erosion classes (slight, moderate, high, severe, very severe) of different land use existed in the Kokcha subbasin. As presented in the table, it is obvious that the highest erosion rates with different classes are found in the rangeland which followed by rainfed agriculture land and barren land. Due to the importance of rainfed agriculture land, the government should take urgent actions for soil conservation, sustainable land use practices for controlling soil erosion risks.

Conclusions
This study proposes an approach to estimate the mean annual soil erosion over the data-scare Kokcha subbasin in Afghanistan. In this way, the RUSLE model, coupled with GIS, is used to estimate soil erosion using the satellite-based data. The following conclusions are obtained from the present study. First, the main fraction (71.34%) of the basin is characterized based on the slight erosion class, while the small area (0.01%) of the study area is induced by very severe soil erosion risk. Second, the central parts of the Kokcha subbasin provide the highest erosion rate. Third, the soil erosion risk is found to be higher in regions with steep slopes, very poor vegetation, and inadequate conservation practices. Fourth, the severe and very severe soil erosion risks are mainly found in mountains/foothills and degraded lands. Fifth, the results of the study have adequate consistency with those obtained from measured data in Khawajaghar station. Sixth, the rangeland has the highest percentage of soil erosion followed by rainfed agriculture land and barren land. Comparison with observed data showed that RUSLE model is able to simulate the seasonal variability of observed soil erosion in the study area. However, a slight overestimation by RUSLE is noticed which may be due to ignoring the deposition and sediment routing. High erosion from agricultural land can significantly reduce soil quality and agricultural sustainability in the region. Different cost-effective adaptation measures including natural confinements of agricultural land, conservation tillage and plant residue management should be adopted for mitigation of soil erosion from agricultural land. Satellite-based precipitation and NDVI data used in this study have considerable uncertainty in estimation. In future, precipitation and NDVI data from different satellite products can be used to assess uncertainty in simulated soil erosion.

Data Availability Statement:
The data that support the findings of this study are available from the corresponding author upon justifiable request.