Water Erosion Risk Assessment for Conservation Planning in the East Hararghe Zone, Ethiopia

Abstract

Water erosion is accelerating soil loss rates in the East Hararghe Zone due to inappropriate human activities and their complex and intertwined interactions with natural factors, particularly in sensitive agroecosystems that lack soil and water conservation (SWC) measures. Although these dynamic processes cause prolonged impacts, a comprehensive assessment of the risk of soil erosion has not yet been undertaken at the zonal level. To bridge this gap, we employed the revised universal soil loss equation (RUSLE) prediction model, along with remote sensing and geographic information systems (GIS), to estimate annual soil erosion rates, analyze the temporal-spatial patterns of erosion risk, and evaluate the potential of standard conservation practices to reduce soil loss in croplands. Total soil erosion (in millions of tonnes/year; Mt yr⁻¹) was estimated to be 9 in 1990, 14 in 2000, 12 in 2010, and 11 in 2020, with average rates of 33, 50, 44, and 39 t ha⁻¹ yr⁻¹, respectively. This suggests an overall 18% increase in soil erosion from 1990 to 2020. Over 75% of the area showed a tolerable soil loss rate (<10 t ha⁻¹ yr⁻¹) and low susceptibility to erosion risk. A mountainous landscape in the northwest presents extremely high erosion (>120 t ha⁻¹ yr⁻¹), which accounts for more than 80% of soil loss, making SWC planning a priority. Analysis of land-use land-cover change (LULCC) confirmed a higher increase in soil loss for LULCC that involved conversion to croplands, with average rates of 36.4 t ha⁻¹ yr⁻¹ (1990–2000), 70 t ha⁻¹ yr⁻¹ (2000–2010), and 36 t ha⁻¹ yr⁻¹ (2010–2020). The results have further revealed that implementing supportive practices such as terracing, stripping, and contouring could reduce average soil erosion by approximately 87%, 65%, and 29%, respectively, compared to the baseline model’s prediction. Therefore, a rigorous cost–benefit analysis is essential to design and implement optimal location-specific practices that maximize investment returns in SWC efforts and ecological restoration. However, we acknowledge the limitations of this study, associated with an empirical model that does not account for all forms of erosion, as well as reliance mainly on secondary data, which may affect the accuracy of the predicted outcomes.

1. Introduction

Water erosion is a pressing global issue with profound geoenvironmental, social, and economic implications [1,2,3,4,5,6,7]. This process involves the detachment, transport, and deposition of soil and rock debris over a single event or extended period, which can result in the loss of fertile topsoil layers and the exposure of underlying materials [4,5,6,7,8,9]. Waterborne erosion is influenced by natural factors, such as surface topography, geological and pedological conditions, dynamic climate patterns (including precipitation and temperature), hydrological and ecological processes, and alterations in groundwater fluctuations [7,10,11]. Unfortunately, inappropriate human activities have sped up the rate of soil loss, leading to soil degradation on different scales that sometimes surpass the rates of natural soil formation [12,13,14]. This, in turn, reduces soil fertility and potential productivity [13,15], ultimately increasing the cost of agricultural production [16]. Accelerated erosion also poses several on-and off-site hazards that vary considerably in space and time. These threats include but are not limited to flooding, habitat loss, limited soil depth for root growth, little water and nutrient storage, loss of soil structure and nutrient movement, sediment accumulation in dams and reservoirs, pollution, freshwater depletion, and their significant impact on biodiversity and natural ecosystems [17,18,19,20,21,22,23,24,25,26,27,28,29].

In 2001, the global annual average potential soil erosion caused by water was estimated to be 35 Pg yr⁻¹ [20]. However, a study conducted in 2015 [1] found an increase in this figure to ${43}_{-7}^{+9.2}$ Pg yr⁻¹, indicating a significant intensification of soil erosion risk over time [2]. Nonetheless, the extent of soil erosion varies worldwide, influenced by agroecological, biodiversity, and microenvironmental conditions [18]. Tropical developing countries suffer from the most significant soil erosion rates, primarily due to heavy precipitation and the profound impact of rapid population growth on land resources [22]. Frequent and severe hydrometeorological extremes (e.g., droughts and floods), high soil erodibility, deforestation, intensive agriculture, and poor conservation practices are additional factors that explain why tropical countries are prone to severe risk of soil erosion [11,23]. Moreover, recent estimates suggest that the rate of soil erosion will likely increase due to changes in land use and land cover (LULC) and future climate scenarios [1]. Regions of the global south, particularly Africa, may encounter high soil erosion risks. To address these multifaceted problems, adopting and implementing comprehensive pathways, such as agroforestry, conservation agriculture, climate-smart agriculture, and agroecology, are critical to protecting, preserving, and sustainably promoting natural resources [24].

Sub-Saharan Africa (SSA), where agriculture is the mainstay of the economy and supports more than 50% of the population, is one of the regions in the world most affected by land degradation [20,25]. Overgrazing and agricultural intensification/expansion are the main drivers of widespread land degradation in SSA. These practices involve the conversion of grasslands and deforestation, leading to increased erosion [17,26,27], which is responsible for more than two-thirds of cropland degradation in the region and attributed to poor soil quality [25,27]. In addition, essential nutrients are removed from agricultural soils during harvest, risking sustainable production. Intensive measures such as fertilizers, manure, compost, and other soil quality improvements must be executed [30]. These challenges are compounded by variable/unpredictable moisture availability, fluctuating precipitation patterns, and recurrent droughts that trigger crop failures and worsen food insecurity [30]. Interestingly, this condition is the same or worse in East African countries, including Ethiopia, where a large part of the population depends on agriculture for their livelihoods [4,6,11,12,13,18,22,23,31,32,33,34,35,36].

In recent years, researchers in Ethiopia have increasingly focused on studying soil erosion, prompted by the need for more precise information on erosion risk’s extent and temporal variations across regions. Scholars have concluded that water erosion that substantially alters landscape characteristics contributes to land degradation, threatening ecological diversity [4,23,34,36,37,38,39,40]. Ethiopia’s net soil erosion by water has been estimated to be around 410 × 10⁶ t yr⁻¹, representing 22% of the gross soil loss rate of 1.9 × 10⁹ t yr⁻¹ [23]. Moreover, it has been found that 25% of the total area’s soil erosion exceeds the national soil loss tolerance value of 10 t ha⁻¹ yr⁻¹ [39]. Croplands contribute 50% of the total, equivalent to 36.5 t ha⁻¹ yr⁻¹, covering 23% of the geographical area [23]. This condition threatens the sustainability of the agricultural sector, which is the pillar of the economy but faces underproductivity due to the lack of advanced technologies [40].

Given that a significant proportion of the country’s population relies on subsistence agriculture and is susceptible to frequent droughts and fluctuations in rainfall, the adverse consequences of soil erosion and land degradation have a cascading effect on livelihoods, food security, and economic growth [4,41,42], both directly and indirectly. In response to this problem, intensive efforts have been made since the 1970s and 1980s, through soil and water conservation (SWC) and restoration of abandoned fields, to reduce the severity of land degradation and improve agricultural productivity [43,44]. However, despite the proactive measures taken over the past decades, erosion-induced land degradation has persisted, suggesting more targeted conservation planning and sustainable land management (SLM) strategies to ensure millions of people’s long-term food security and livelihoods. This issue underscores the importance of reliable decision-supporting evidence to advance the implementation of SLM practices. Therefore, evaluating past and present soil erosion rates and their variations/patterns can improve our understanding of the extent to which the land is susceptible to soil erosion risk, quantitatively and geographically [42,45].

Understanding soil erosion processes in the context of reciprocal influences of natural and human-caused drivers on the occurrence, origin, and severity of erosion risk, on the one hand, and estimating soil material generation, on the other, requires a more systematic approach [10,18]. Generally, soil erosion assessment can be performed using field measurements or spatial modeling [7,45,46]. On-site measurements involve the physical observation of natural areas and yield precise results, but they are time-consuming and ineffective for larger areas [35]. In contrast, models that predict spatially distributed soil erosion rates have gained wide usage and can be classified into three major groups based on their data requirements, complexity, algorithm setups, and simulation procedures: (i) physically based models, (ii) stochastic models, and (iii) empirical models [21,37,47,48].

Of the empirical models, the revised universal soil loss equation (RUSLE) [49], a derivative of the universal soil loss equation (USLE) model [50], in particular is the most operational prediction model in which long-term impact caused by natural and human factors on sheet-and-rill erosion is considered. The model has helped with quantitative soil loss estimates and erosion risk mapping on temporal (that is, daily, monthly, seasonal, and annual basis, depending on data availability [3]) and geographic scales from a plot or minor catchments (watersheds) to more large hydrological basins under various climatic conditions, soils, and management strategies [20,23,31,34,38,51]; Hurni has also adopted it for Ethiopia’s context [52]. The RUSLE is a “parametric model” [3] that has increasingly been applied in different regions around the globe, including agriculture-dominated landscapes in Ethiopia [4,18], Poland [7], and Italy [9], Lake Kivu basin in DR Congo-Rwanda to assess soil erosion resulting from deforestation [12], the Mediterranean conifer forest to account for the response of erosion to forest disturbances [53], and protected areas (Greece) [54]. In China, the RUSLE model has been applied to appraise the impact of afforestation measures on reducing the risk of soil erosion in the Three Gorges Reservoir Region [48]. These examples establish the practicability of the model.

Furthermore, with cutting-edge developments in geographic information systems (GIS) tools and state-of-the-art remote sensing technology, spatial modeling approaches have become increasingly important in environmental hazard analysis [4,38,55]. However, soil erosion estimates based on models have uncertainties or remain challenging to validate, typically in cases where field-measured references or relevant quality data are inadequate or inaccessible [1,10,39,46]. Despite the challenges, evaluating soil erosion rates at local and global levels is feasible using reliable geospatial datasets on biophysical and hydrological parameters with an operational predictive model [4,20]. To this end, the accessibility of high-resolution Earth observation datasets, open-access tools, and cloud computing platforms (e.g., Google Earth Engine, GEE) has made large-scale soil erosion assessment feasible [9,53,54]. This development has created a favorable condition for more accurate erosion risk identification, mapping, and prioritization for SWC planning.

The East Hararghe Zone (EHZ) in Eastern Ethiopia is home to some of the most vulnerable landscapes that experience severe soil erosion and land degradation, leading to serious environmental and socioeconomic problems [18,56]. Due to an ever-growing population, the demand for food, fuel, fiber, and animal feed has increased, which has had an interrelated impact, resulting in extensive degradation at different scales [4]. Similar to various other ecoregions in the country, the combined effects of the changing climate and LULC, coupled with intense rainfall and unsustainable activities such as deforestation and cultivation on steep slopes, have caused widespread runoff and soil erosion [38,57,58]. These have led to a loss of nutrients and soil fertility, reduced productivity, siltation, damage to irrigation canals, and shorter reservoir lifespan [33,44,59].

In past decades, several studies have estimated soil loss rates at the catchment scale [60,61,62,63,64,65]; More recent work demonstrated the potential of combining prediction models with spatial analysis to identify areas with high soil erosion risk for the planning of SWCs on the watershed/subbasin scale [4,18]. However, prior study has focused primarily on small geographical areas and neglected spatiotemporal patterns of erosion risk in the whole landscape at the zonal level. Here, we aim to address knowledge gaps by examining the extent of soil erosion in the context of the long-term impacts of natural and anthropogenic factors. Our specific objectives are to: (i) estimate soil loss rates for the years 1990, 2000, 2010, and 2020; (ii) analyze the prolonged evolution of soil erosion risk in the 1990–2020 period; (iii) investigate the effects of LULC, rainfall, and topography; and (iv) evaluate the implications of conservation practices in reducing soil loss in cropland areas.

2. Study Area

The EHZ (latitude: 07°39′–09°43′ N; longitude: 07°30′–09°47′ E) is located in the Oromia Regional State and has a total surface area of approximately 25,040 km², of which croplands cover nearly 50%. The zone area is divided into 20 woredas—mid-level administrative entities equivalent to districts, as shown in Figure 1a–c. EHZ, being an enclave of the Harari Region (a masked portion within the boundary of the studied landscape, Figure 1c), is located spatially adjacent to the Dire Dawa in the north, Fafan in the east, southeast, and northeast, Erer in the south, the East Bale in the southwest, the West Hararghe in the west, and the Siti in the northwest. Significant altitudinal variations were observed in the study area, with the northern part characterized by a complex topography of deep gorges, rugged mountains, and plateaus. The southern part has relatively flat terrain. The altitudinal range extends from 520 m at Meyu Muleke in the southernmost extent to approximately 3386 m at the peak of Gara Muleta in the northern part, shared by four districts: Girawa, Bedeno, Kersa, and Kurfa Chele. The average altitude is 1505 m above sea level.

The rainy season, locally known as Ganna, occurs from June to August and has an average maximum precipitation of 1200 mm yr⁻¹ (Figure 2). The mean lowest temperature (T_min) and maximum temperature (T_max) are 10 °C and 33 °C, respectively [66]. Sorghum and maize are the main crops grown in most districts and are grown as sole crops or intercropped with chat (Catha edulis), fava beans, common beans, sweet potatoes, field beans, and potatoes [67,68]. Livestock rearing is also an important economic activity, and small-holder farming communities rely on chat and coffee (Coffea arabica L.) as sources of cash income.

3. Materials and Methods

3.1. Geospatial Data Sources

This study used multisource remote sensing datasets to generate input parameters for the RUSLE model. Multitemporal spatial and spectral resolution Landsat satellite images (Level 2 products; surface reflectance, SR; Table S1) were used for LULC classification and change analysis. Landsat-5 TM (Thematic Mapper) images were used for 1990, 2000, and 2010, while Landsat-8 OLI/TIRS (Operational Land Imager/Thermal Infrared Sensor) images were used for 2020. To minimize cloud contamination and phonological variations, we selected and acquired five Landsat scenes that cover the study area (path/row: 166/53; 166/054; 166/055; 167/53; 167/054) for the dry period of each year (December, January, and February), with less than 5% cloud cover. All Landsat images were downloaded from the United States Geological Survey (USGS) online service at the EROS Center [69]. The United Nations (UN)-Food and Agriculture Organization (FAO)-Educational Scientific and Cultural Organization (UNESCO) Digital Soil Map of the World (DSMW) and the consistent attributes of sand, silt, clay, and organic carbon content of the soil classes were accessed in the form of an Environmental System Research Institute (ESRI) shapefile [70]. The Land Process Distributed Active Archive Center (LPDAAC) provided the Shuttle Radar Topographic Mission (SRTM) void-filled digital elevation model (DEM) with a spatial resolution of 30 m × 30 m [71]. The DEM data generated thematic maps for topographic features, including slope length and steepness factor. To create a rainfall erosivity map, we accessed monthly rainfall data from 1981 to 2020, which had a spatial resolution of 0.05 × 0.05°, from the Climate Hazards Group Infrared Precipitation with Stations (CHIRPS) [72].

The sample points representing LULC types in the study area were collected using the stratified random sampling method from high-spatial-resolution images in Google Earth Pro (Google Inc., Amphitheatre Parkway, Mountain View, CA, USA) archives. Of the total samples, 70% of the data points were used to train and run the classification algorithm. In comparison, the remaining 30% were reserved for post-classification validation to calculate the accuracy of the classified images for each year. Moreover, we referred to prior studies (e.g., [4,18,63,65]) and consult with local experts with a deep knowledge of the study area for further information.

3.2. Methods

3.2.1. Land-Use and Land-Cover Classification

The accessed Landsat satellite images were calibrated, corrected, and georeferenced without requiring radiometric atmospheric corrections or rectifications. Essential preprocessing steps, such as multispectral band composite and mosaicking, were performed, and images covering the area of interest (AOI) were extracted. Subsequently, successive image classification frameworks were developed to produce LULC maps. This involved preliminary training to select the spectral signature with polygon layers from composite images and actual classification [73]. The first step was assigning reference points to LULC classes using a stratified random sampling method [74]. The spectral signature representing each LULC class was collected using the spatial location of the training points. Then, Landsat TM/OLI images were classified into seven LULC types (bare land, built-up, cropland, forestland, grassland, shrubland, and waterbody) for 1990, 2000, 2010, and 2020 using a supervised classification method based on a support vector machine (SVM) algorithm [74,75] (Table S2). Lastly, transition matrices were developed by overlaying classified LULC images, showing the conversion between different LULC classes within three ten-year windows: 1990–2000, 2000–2010, and 2010–2020.

Accuracy assessment is the most intuitive task for optimizing LULC map applications derived from remotely sensed data. It is commonly achieved by comparing the classified LULC image with field-established references or other sources, which provides confidence regarding the robustness of classification algorithms. In this study, due to the absence of field data, the classification accuracy of the LULC map was validated using point data marked in predefined areas from Google Earth Pro imagery. Fine-resolution satellite images also helped with LULC classification and accuracy verification in regions where field truth data were scarce but demonstrated promise for reducing misinterpretations of remote sensing images [18]. This study assessed accuracy by comparing validation points to their respective locations in classified images.

We used standard statistical metrics such as user accuracy (UA), which determines the accuracy of existing data for class i; producer accuracy (PA), which evaluates how well the model predicts new information for class j; overall accuracy (OA), which indicates the performance of the classifier (model) overall; and Kappa (K^) coefficients [76,77,78]. Moreover, the percentages of commission errors (CE), which are unintended exclusions from a given category and can be calculated as 100 minus the user’s accuracy (UA), and omission errors (OE), which refer to included references to an incorrect class and can be calculated as 100 minus the producer’s accuracy (PA), were acquired from the confusion matrix [74]. Equations (1)–(3) were used to compute UA, PA, and OA [79] for every year’s LULC map:

$$U{A}_i=\frac{n_{ii} }{n_{i+}},$$

(1)

$$P{A}_j=\frac{n_{jj}}{n_{+j}},$$

(2)

$$OA=\frac{{{\displaystyle \sum }}_{i=1}^n{n}_{ii}}{n},$$

(3)

In a confusion matrix, n_ij refers to the value at the intersection of row i and column j. K^ was calculated using Equation (4) as follows [80,81]:

$${\mathrm{K}}^{\wedge }=\frac{\mathsf{{\rm N}}{{\displaystyle \sum }}_{\mathrm{i}=1}^{\mathrm{r}}{\mathsf{\chi }}_{\mathrm{ii}}-{{\displaystyle \sum }}_{\mathrm{i}=1}^{\mathrm{r}}\left({\mathsf{\chi }}_{\mathrm{i}+}\times { \mathsf{\chi }}_{+\mathrm{i}}\right)}{{\mathrm{N}}^2-{{\displaystyle \sum }}_{\mathrm{i}=1}^{\mathrm{r}}\left({\mathsf{\chi }}_{\mathrm{i}=1}\times {\mathsf{\chi }}_{+\mathrm{i}}\right)},$$

(4)

where N represents the total number of observations, r indicates the number of rows or columns in the matrix, ${\mathsf{\chi }}_{\mathrm{ii}}$ represents the number at the intersection of row i and column i, ${\mathsf{\chi }}_{+\mathrm{i}}$ represents the sum of values in row i, and ${\mathsf{\chi }}_{\mathrm{i}+}$ shows the sum of values in column i.

3.2.2. Derivation of the RUSLE Model Factors

The RUSLE model was selected for estimating soil loss rates in the present study area. The primary factors that drove this decision were the model’s simple process and minimal geospatial and numerical datasets requirement. In addition, the model simplifies a complicated system to a reasonably basic one while defining the fundamental governing principles inducing soil erosion by water. When integrated with GIS-spatial analysis tools, its compatibility allows for practical soil loss estimation at different scales [9,20], and increasing predictive capabilities is an ideal choice for this research. For the RUSLE application, we derived thematic factors using the spatial analyst toolset in ArcGIS Desktop software version 10.8 (Environmental Systems Research Institute (Esri), Inc., Redlands, CA, USA). Each thematic layer developed in a raster format to ensemble the model was harmonized to a 30 m × 30 m spatial resolution and reprojected to a standard reference of the World Geodetic System 1984 spheroid with the Universal Transverse Mercator (UTM) and Adindan Zone 37 N. Lastly, the input factors of the RUSLE model (A), which contribute to soil erosion, such as rainfall erosivity, soil erodibility, topography, cover management, and conservation practices, were multiplied by Equation (5) [49].

$$A=R\times K\times LS\times C\times P,$$

(5)

where A is the average annual soil loss rate; R is the rainfall-runoff erosivity factor; K is a soil erodibility factor; LS is a slope length-steepness factor (dimensionless); C is a cover management factor (dimensionless); P is a support practice factor (dimensionless).

We used the natural breaks method, or ‘Jenks classification’ [82], to categorize the severity of the erosion risk. This method groups data into classes based on their ‘natural’ groupings, resulting in an optimized arrangement of values into ‘natural’ categories [83]. Additionally, class breaks are determined to similar group values while maximizing differences between classes to ensure accuracy and minimize potential errors in the classification process [82,83]. Consequently, the study area was characterized into eight soil erosion risk classes, following the literature [4,18]: very-low (0–5 t ha⁻¹ yr⁻¹), low (5–10 t ha⁻¹ yr⁻¹), low-medium (10–20 t ha⁻¹ yr⁻¹), medium (20–40 t ha⁻¹ yr⁻¹), high-medium (40–60 t ha⁻¹ yr⁻¹), high (60–80 t ha⁻¹ yr⁻¹), very-high (80–120 t ha⁻¹ yr⁻¹), and extremely high (>120 t ha⁻¹ yr⁻¹). Furthermore, the administrative districts prone to soil erosion risk at different scales were prioritized to help plan for developing suitable precaution or prevention measures against land degradation, improving productivity and ecosystem sustainability.

Rainfall Erosivity (R) Factor

Factor R, a significant driving agent of soil erosion [11,16], accounts for about 80% of soil loss [50]. It represents the erosion intensity of rainfall-runoff to detach and transport the soil materials due to the combined forces of rainfall volume, kinetic energy, duration, and potential [10,45,84]. In this study, CHIRPS datasets [72] were first aggregated into an annual time series based on the composite average monthly rainfall data within four-time frames, namely 1981–1990, 1991–2000, 2001–2010, and 2011–2020. We then separately calculated the factor R for the study years 1990, 2000, 2010, and 2020 based on the annual precipitation values. We used the R factor computation method (Equation (6)) proposed by Lo et al. [85] as a robust system widely applied in recent studies in East Africa, confirming good prediction results [4,11,22].

$$R=38.46+3.48\times P,$$

(6)

R is the rainfall erosivity (MJ mm ha⁻¹ h⁻¹ y⁻¹), and P is the average annual rainfall.

Soil Erodibility (K) Factor

The soil erodibility factor, K, determines the resistance of the aggregate soil particles to detachment and transport by rainwater [11]. The K factor represents a numerical value of how susceptible a particular soil is to erosion, depending on inherent soil profile characteristics and properties involving texture, organic matter, structure, permeability, and human-caused impacts [18,34]. It accounts for the resultant soil loss rate determined for erosion index unit at a standard experimental plot (“72.6 ft long on a 9% steepness”) condition [49]. The K factor was calculated using the erosivity-productivity impact calculator (EPIC) tool [86] based on the DSMW [70] with the values (accessed from the Excell document “Generalized_SU_Info.xls”) that reflect the relative ratios of sand (SAN), silt (SIL), clay (CLA), and organic carbon (C) contents for each topsoil class (Equation (7a–d)) as:

$$K_{USLE }=F_{csand}\times F_{cl-si}\times F_{orgC}\times F_{hisand},$$

(7)

where:

$$F_{csand}=\left\{0.2+0.3exp\left[-0.0256SAN\left(1-\frac{SIL}{100}\right)\right]\right\},$$

(7a)

$$F_{cl-si}={\left[\frac{SIL}{CLA+SIL}\right]}^{0.3},$$

(7b)

$$F_{orgc}=\left[1.0-\frac{0.0256C}{C+exp\left(3.72-2.95C\right)}\right]$$

(7c)

$$F_{hisand}=\left[1.0-\frac{0.70 SN1}{SN1+exp\left(-5.51+22.9SN1\right)}\right],$$

(7d)

SN1 is sand content subtracted from 1 and divided by 100.

Slope Length and Slope Steepness (LS) Factor

The topographic factor, LS, signifies a soil loss rate for a particular site condition with a slope length (L) of 22.13 m and a steepness (S) of 9%, clear vegetation, and in a seedbed condition [18]. The LS factor for the present study area was determined using NASA’s SRTM-DEM [66] with a 30 m cell size. The Arc-Hydro extension and raster calculator were used as geoprocessing tools in the ArcGIS Desktop software version 10.8. As per the method outlined in various prior studies, the LS factor consists of two sub-factors: L and S [87,88,89], computed as follows (Equations (8)–(12)):

$$\mathrm{L}={\left(\frac{\lambda }{22.1}\right)}^m,$$

(8)

where λ is the field slope length in meters, and m is the variable L exponent related to the value of the slope gradient.

Literature suggests an m value of 0.5 for a slope gradient of 4.5% or more, 0.4 for slopes 3–4.5%, 0.3 for slopes 1–3%, and 0.2 for slopes ≤1% [47]. We subdivide the S sub-factor into segments for the up-slope drainage parts [88,90].

$$L_{i,j}=\frac{ {\left(A_{i,j-in}+D^2\right)}^{m+1}-A_{i,j-in}^{m+1}}{D^{\mathrm{m}+1}·x_{\mathrm{i},\mathrm{j}}^m·{22.13}^m},$$

(9)

$$\mathrm{m}=\left(\frac{\mathsf{\beta }}{1+\mathsf{\beta }}\right),$$

(10)

$$\beta =\frac{Sin \theta /0.0896}{3{\left(Sin \theta \right)}^{0.8}+0.56},$$

(11)

$$S_{i·j}=\left\{\begin{array}{c}10.8sin{\mathsf{\theta }}_{i.j}+0.03, tan {\mathsf{\theta }}_{i.j}<9\% \\ 16.8sin{\mathsf{\theta }}_{i.j}-0.50, tan {\mathsf{\theta }}_{i.j} \ge 9\%\end{array}\right.$$

(12)

where $A_{i,j-in}$ is the contributing area at the grid cell inlet ($i,j$) in m²; D is the grid cell size; $x_{\mathrm{i},\mathrm{j}}$ is sin $a_{i,j}$ + ${\cos }_{i,j}$; $a_{i,j}$ is the aspect direction of the grid cell ($i,j$); m is a variable slope length exponent related to the ratio of β of rill erosion produced by flow and inter-rill erosion caused by raindrop impacts, and θ is the slope angle (gradient).

Cover Management (C) Factor

In the context of the RUSLE application, the cover management factor, C, indicates the sum of soil loss from land with specific vegetation cover to soil loss under clean, continuous tilled fallow or managements to reduce soil loss via water erosion [18]. The factor C value ranges from 0 for non-erodible conditions in areas with significant vegetation cover and insignificant anthropogenetic impacts to 1, which resembles a more considerable magnitude of soil loss through extensive tillage, which leaves a very smooth surface that creates much rainfall-runoff and makes the soil susceptible to erosion [49]. This shows that vegetated areas have a lower C value; they reduce the eruptive force of rainfall and stormwater runoff [16]. The C coefficients associated with the LULC types were obtained from scholarly literature [28,60,87,91,92,93].

Support Practice (P) Factor

The P factor represents the conservation measures implemented in a particular AOI to protect soil from degradation [49]. These interventions allow for the preservation of soil as a biological habitat, returning organic matter to the soil to play an essential role in establishing a healthy and rich ecosystem, while maintaining its value for the future [94,95]. Other contributions of soil conservation include (i) controlling runoff, (ii) protecting highly exposed bare surfaces and landscapes, such as steep slopes, and (iii) keeping downstream areas from sedimentation and pollution [96,97].

The P factor is characterized as the cumulative soil loss following erosion control practices to soil loss after conservation practices [49,98], which are critical mechanisms in reducing the soil loss caused by runoff velocity and concentration, drainage patterns, and hydraulic drives on the soil surface [4,97]. The P factor value ranges from zero to one, with zero showing fields with effective conservation practices and one implying areas with little or no effective conservation practices [57].

There are extensive SWC practices and management efforts to control soil erosion in the study area; due to the non-availability of a digital map of conservation practices, we correlated the P factor to its consistent coefficients for cropland and the non-cropland areas based on the LULC maps (from 1990, 2000, 2010, and 2020) and the slope classification based on a 30 m × 30 m grid cell DEM. Accordingly, the P factor values were determined as follows: for the cropland areas with slope 0–5%, P = 0.1; 5–10%, P = 0.12; 10–20%, P = 0.14; 20–30%, P = 0.2; 30–50%, P = 0.25; >50%, P = 0.3; non-cropland classes P = 1, adhering to the suggestions of Wischmeier and Smith [50]. Additionally, we conducted separate modeling to estimate the soil erosion rates in croplands for the year 2020. We considered the P factor and its associated characteristics, derived from the referred sources [11,12,15,22], for three management supports practices: contouring, stripping, and terracing, assuming that these conventional practices could be implemented in the cropland areas.

4. Results and Discussion

The application of remotely sensed data and GIS, in combination with an approach based on RUSLE modeling, has provided estimates of soil erosion by water on a pixel-by-pixel basis (i.e., at a spatial resolution of 30 m × 30 m) in the East Hararghe landscape of Ethiopia. The following sections present spatiotemporal trends and dynamics of soil loss rates and erosion risks for various slope gradients, rainfall patterns, and LULC changes. Additionally, we have identified hotspot districts susceptible to higher soil erosion risks as priority areas for future intervention in SWC measures.

4.1. Overview of the RUSLE Factors

Five map layers of input parameters (R, K, LS, C, and P) were utilized based on the RUSLE to model the spatially distributed soil erosion rate during four reference periods, i.e., 1990, 2000, 2010, and 2020. These timeframes enabled us to analyze the potential impact of policy shifts and land management over the past three decades.

The R factor varied from 1318 to 4050 MJ mm ha⁻¹ h⁻¹ y⁻¹ in 1990 (Figure 3a), 1073 to 3744 MJ mm ha⁻¹ h⁻¹ y⁻¹ in 2000 (Figure 3b), 1227 to 4029 MJ mm ha⁻¹ h⁻¹ y⁻¹ in 2010 (Figure 3c), and 1306 to 4308 MJ mm ha⁻¹ h⁻¹ y⁻¹ in 2020 (Figure 3d). The average values for the same periods were 2429, 2160, 2378, and 2416 MJ mm ha⁻¹ h⁻¹ y⁻¹, respectively. There are distinct variations in rainfall patterns between districts, with the northern and northwest landscapes receiving relatively higher rainfall than the southern and southeastern parts, resulting in higher R values in the uplands of the north and northwest and the lowest R values in the lowlands of the south and southeastern sections (Figure 3a–d).

Furthermore, the landscape that has been studied is characterized by eleven soil classes, illustrated in Figure 4a, including Haplic Xerosols (Xh), Calcaric Regosols (Rc), Haplic Yermosols (Yh), Eutric Cambisols (Be), Dystric Cambisols (Bd), Eutric Nitosols (Ne), Eutric Regosols (Re), Humic Cambisols (Bh), Calcic Xerosols (Xk), Chromic Vertisols (Vc), and Cambic Arenosols (Qc). These soil classes contributed to the zonal area in varying proportions, with Xh representing the most significant proportion, followed by Rc and Bd (Table S3). The erodibility of the soil classes varied considerably, with K values ranging from 0.3 t h MJ⁻¹ mm⁻¹ for Ne (which had a lower organic C content of <11% and were located in the northeast and southeast) to 1 t h MJ⁻¹ mm⁻¹ for Be in the northwest. The K value for each soil class can be found in Table S3, and a spatial variation in erodibility between different soil types is shown in Figure 4b. The average K-value for all soil classes is approximately 0.96 t h MJ⁻¹ mm⁻¹ with a standard deviation of 0.3 t h MJ⁻¹ mm⁻¹.

The LS factor is a crucial parameter in estimating soil erosion rates, considering the impact of terrain and surface topography on erosion. Steeper and longer slope gradients typically result in higher surface runoff and flow speed, leading to significant soil loss due to the greater influence of topography. However, gentle slopes usually have lower erosion rates. The study area exhibited a range of slope gradients, with a considerable share comprising gentle slopes with a rise of less than 7%, covering about 44% of the total area (Figure 4c). On the contrary, steeper terrains with slopes of more than 27% were predominant in regions with a high LS factor, covering almost 19% of the area. Thus, the northeast, northwest, and southeast parts had the highest LS values and steeper gradients (Figure 4b).

Factor C is an essential component that reflects the impact of soil cover, such as the vegetation canopy and crop residues, surface roughness, and soil moisture conditions, on soil loss resulting from land management practices [11,33]. We adopted the value of the C factor for LULC classes from the coefficients published in prior studies from Ethiopia and other East African countries, which ranged from 0.08 to 0.1, with a standard deviation of 0.1 (Figure 5a–d; Table S4). Meanwhile, the P factor reflects the implication of different conservation practices in reducing soil loss via water erosion, as shown in the thematic maps (Figure 6a–d), compensating for the lack of available digital data. Additionally, the P-factor values for supporting practices such as contouring, strip-cropping, and terracing are provided in Table S5. Previous studies (e.g., [11,12,15,22]) have identified these practices as practical intervention measures to reduce cropland soil loss. Thus, a comprehensive evaluation of the potential of range support practices on soil loss reduction could help to apply the most efficient measures.

4.2. Soil Erosion Dynamics over 30 Years (1990–2020)

The total amount of soil lost via water erosion was 9 × 10⁶ t in 1990, 14 × 10⁶ t in 2000, 12 × 10⁶ t in 2010, and 11 × 10⁶ t in 2020, as we reported in Tables S6–S9. The resulting spatial patterns, depicting locations with varying levels of erosion risk, are shown in Figure 7. The figure illustrates that areas with lower LS, R, and C factors in the southern, southeast, and northeast regions were less affected during in each study period. At the same time, steeper slopes and higher rainfall increased erosion risk. This highlights the importance of topography and rainfall in influencing soil erosion. Our statistical results and closer visual analysis of the erosion risk maps have revealed significant temporal-spatial variations in soil erosion caused by running water. These findings are consistent with a similar study conducted in Ethiopia [38,58,99,100], including the catchment/watershed or sub-basin scale [4,18], which has found heterogeneity in soil loss rates.

The Ethiopian highlands—where human settlement and agriculture have been into practice for at least three millennia, resulting in modifications to an equilibrium ecosystem, covering around 44% of the land area, and supporting the livelihoods of 87% of the population—are the most vulnerable regions to severe soil erosion risk [42,56,101]. The average annual soil loss rate in highland areas is estimated to be 300 t ha⁻¹ yr⁻¹. Related to this problem, Woldemariam and Harka [4], Haregeweyn et al. [44], and Akale et al. [59] reported that multiple human-induced and natural factors contribute to a higher erosion rate in most of the highland regions. These factors include rugged terrain, low vegetation cover, poor SWC practices, and intense precipitation, particularly unsustainable human activities. This has accelerated soil erosion and caused on-site and off-site impacts. Moreover, this problem is exacerbated by increasingly frequent and severe weather events such as droughts and floods, making it challenging to meet the food needs of the rapidly growing population and improve their standard of living [55,102].

Modeled erosion rates (annual average; t ha⁻¹ yr⁻¹) in 1990, 2000, 2010, and 2020 were 33, 50, 44, and 39, respectively (Tables S6–S9). The average soil loss estimates for the four reference periods studied are much higher than the average rate of 29 t ha⁻¹ yr⁻¹ estimated for the West Hararghe Zone, which ranges from 1 t ha⁻¹ yr⁻¹ to over 400 t ha⁻¹ yr⁻¹ [103]. Similarly, our estimates were two to three times higher than the national average rate of 16.5 t ha⁻¹ yr⁻¹, which ranges from 0 to 200 t ha⁻¹ yr⁻¹, reported in Ethiopia [23]. Furthermore, the modeled average rates are considerably higher than the “tolerable” soil loss rate (10 t ha⁻¹ yr⁻¹) predicted for the country [39] and the standard soil loss acceptances, which range from 5 to 11 t ha⁻¹ yr⁻¹ [50]. The results also exceed the average soil loss values reported in other East African countries [11,22]; however, they are relatively consistent with the average soil loss rate in Ekiti State, Nigeria [104].

The estimated soil loss in the present study area ranged from 0 to 1502 t ha⁻¹ yr⁻¹, 0 to 1005 t ha⁻¹ yr⁻¹, 0 to 1290 t ha⁻¹ yr⁻¹, and 0 to 1151 t ha⁻¹ yr⁻¹, respectively. The results revealed that in 1990, 2000, 2010, and 2020, 81%, 80%, 80%, and 82% of the study area had a tolerable soil loss rate, respectively. The annual soil loss estimate range was low compared to previous researchers’ findings for Ethiopia’s highland areas, which reported annual soil loss rates ranging from 1248–23,400 million tons [56]. In 1990, 19% of the zone’s surface area was occupied by sites that experienced an average soil erosion rate of more than 10 t ha⁻¹ yr⁻¹. However, this figure increased to 20% in 2000 and 2010 before decreasing to 18% in 2020, representing a 5% decrease from 1990 to 2020. Conversely, our estimates are lower than those reported in other local-scale studies, such as Belayneh et al. [60] in the Gumara Watershed, Eniyew et al. [33] in the Telkwonz Watershed, and Wagari and Tamiru [13] in the Fincha Catchment. Our analysis shows that various factors, such as temporal-spatial scale and uncertainties in prediction models, may contribute to the variation in soil erosion rates modeled in this study and previous works.

Each cell’s annual soil erosion rates were categorized into eight levels according to the severity of the erosion risk (Figure 7).

Table 1. Mean soil loss rates per erosion risk class and contribution to total annual soil erosion, mean slope, and mean rainfall in 1990, 2000, 2010, and 2020.

Year	Erosion Risk Class	Area (km2)	% of the Area Covered	Mean Soil Loss (t ha−1 yr−1)	% of Total Soil Loss	Mean Rainfall (mm yr−1)
1990	Very low	19,769	78.9	0.1	0.3	691
	Low	581	2.3	7	0.5	656
	Low Medium	1067	4.3	15	1.9	635
	Medium	1179	4.7	29	4.1	648
	High medium	621	2.5	49	3.7	668
	High	381	1.5	69	3.2	685
	Very high	408	1.6	97	4.8	703
	Extremely high	1034	4.2	650	81.5	727
	Total	25,040	100	916	100	687
2000	Very low	19,253	76.9	0.2	-	615
	Low	361	1.4	4	0.1	566
	Low Medium	496	2	7	0.3	564
	Medium	1049	4.2	15	1.2	543
	High medium	1094	4.4	29	2.5	559
	High	899	3.6	57	4.1	588
	Very high	408	1.6	98	3.2	618
	Extremely high	1480	5.9	757	88.6	665
	Total	25,040	100	966	100	610
2010	Very low	19,106	76.3	0.1	-	679
	Low	313	1.3	4	0.1	640
	Low Medium	528	2.1	7	0.4	630
	Medium	1104	4.4	15	1.5	608
	High medium	1204	4.8	29	3.1	627
	High	997	4	57	5.2	655
	Very high	438	1.7	98	3.9	677
	Extremely high	1350	5.4	698	85.8	702
	Total	25,040	100	907	100	672
2020	Very low	19,381	77.4	0.3	0.1	689
	Low	553	2.2	3	0.2	642
	Low Medium	606	2.4	7	0.5	647
	Medium	1019	4	15	1.5	622
	High medium	1072	4.3	29	3.2	635
	High	871	3.5	57	5.1	660
	Very high	371	1.5	98	3.7	686
	Extremely high	1167	4.7	718	85.7	739
	Total	25,040	100	926	100	683

provides statistical details on the severity of soil erosion risk, showing that most areas had very low levels of risk, with coverage of 79%, 77%, 76%, and 77% in 1990, 2000, 2010, and 2020, respectively. These findings are consistent with studies from eastern Ethiopia’s landscapes, such as Woldemariam et al. [18] for the Gobele Watershed and Woldemariam and Harka [4] for the Erer Sub-Basin, who also reported a significant area of very low risk. In addition, Chala and Dadi [103] found similar results for the West Hararghe Zone, indicating that 72% of the area had low to moderate soil loss rates, while 28% had high to very high rates.

As shown in Table 1, the erosion loss rates varied between different slope classes and years, indicating that all these factors affect soil erosion rates. Analysis of temporal erosion risk reveals a decrease in the proportion of areas with very low soil loss between 1990 and 2020, with a reduction in the low erosion risk class during the same period (

Table 2. The rates of spatial changes per erosion risk class in the East Hararghe Zone, Ethiopia.

Erosion Risk Class	1990–2000		2000–2010		2010–2020		1990–2020
	km2	%	km2	%	km2	%	km2	%
Very low	−516	−3	−147	−0.8	275	1	−388	−2
Low	−220	−38	−48	−13	240	77	−28	−5
Low Medium	−571	−54	32	7	78	15	−461	−43
Medium	−130	−11	55	5	−85	−8	−160	−14
High medium	473	76	110	10	−132	−11	451	73
High	518	136	98	11	−126	−13	490	129
Very high	-	-	30	7	−67	−15	−37	−9
Extremely high	446	43	−130	−9	−183	−14	133	13

). On the contrary, the most pronounced reduction was observed in soil loss of medium intensity, classified as low-medium risk (−43%), followed by the medium class. Similarly, the very high risk of erosion showed a marked reduction. In contrast, the high erosion risk showed the most significant increase, followed by high-medium and extremely high classes in areas with high rainfall levels and rugged terrain. The share of erosion risk classified as low, low-medium, medium, and high-medium declined significantly, along with their contribution to annual soil loss.

In contrast, the high erosion risk accounted for approximately 3.2% of soil losses in 1990 but increased to 5.1% in 2020. This represents a 59.4% increase in its contribution to total soil losses over the studied period, as shown in Table 2. On the other hand, the percentage of areas affected by high-to-medium erosion increased from 8.7% to 11.8% of the total area, along with an aggregated mean annual soil loss rate. However, the contribution of erosion risk to total soil loss decreased in the same period.

Furthermore, an extremely high erosion risk, with an average soil loss rate of more than 120 t ha⁻¹ yr⁻¹, accounted for 4.2% of the landscape area in 1990, 6% in 2000, 5.4% in 2010, and 4.7% in 2020. However, the corresponding annual contribution to soil loss was 81.5%, 88.6%, 85.8%, and 85.7%, respectively, in highland ecosystem districts in the northwest. Overall, the results show that erosion risk classes with higher mean soil loss rates contribute more to total annual soil erosion, highlighting the importance of controlling soil erosion in areas with extremely high risk, which are more susceptible to significant soil loss.

Erosion control measures have been implemented in the studied landscape since 2000 (Figure 8). These measures have likely improved the state of vegetation, resulting in enhanced overall forest and shrub covers and reducing soil erosion. The findings indicate an average soil loss rate increase of 17 t ha⁻¹ yr⁻¹ during the first period (1990–2000), accompanied by an estimated total annual soil loss increase of 4.9 × 10⁶ t at the zonal level. However, the total and average soil loss rates decreased during the second (2000–2010) and third (2010–2020) periods. Our analysis shows that the total annual and average erosion decreased by 1.8 × 10⁶ t and 6 t ha⁻¹ yr⁻¹ in the second period and 1.4 × 10⁶ t and 5 t ha⁻¹ yr⁻¹ during the third period, respectively. As a result, the magnitude of soil erosion decreased slightly from 2000–2020. Additionally, reduced rainfall erosivity values (Figure 3a–d) likely contributed to decreased soil loss, particularly during 2000–2010.

Our findings are consistent with those of the watershed-level results from the northeast and southeast parts of the present study area [18]. They estimated a 33% reduction in average soil loss rates, which could be attributed, in part, to the implementation of conservation measures. Interestingly, despite the decrease in erosion rates, the eroded landscape at the watershed level increased by 20% [18]. However, our model estimates support the positive impact of conservation measures on soil erosion, coherent with research conducted in Tanzania [32]. Other researchers (e.g., [3,105,106]) congruently confirm that sustainable land management solutions provide benefits beyond SWC, such as improving soil nutrients, infiltration, and moisture, enhancing the soil’s physio-chemical properties, and protecting biodiversity.

Meanwhile, in the EHZ, as shown in Figure 8, SWC practices and a slight drop in rainfall erosivity might partially explain the reduction in soil erosion magnitude observed over the past two decades (from 2000 to 2020). However, despite this incidence, estimated soil erosion rates have increased by 18% (equivalent to 1.7 × 10⁶ t) in the last 30 years, from 1990 to 2020 (Tables S6–S9). The trend of increased soil erosion in the study area requires further investigation and intervention to mitigate its detrimental impact. Furthermore, our findings were strongly supported by compelling evidence. Specifically, we identified significant areas of soil erosion where conservation measures in the past were ineffective or not yet fully implemented, as shown in Figure 9a–f.

Failure to implement, enforce, or sustain successful conservation measures in eroded landscapes can be the main factor in increased soil loss during the studied period and could further exacerbate the problem. In connection to this, Guzman et al. [107] reported that past interventions in SWC have been insufficiently effective, contributing to the persistent sediment concentration in rivers of the Ethiopian Highlands, despite extensive conservation efforts. Thus, renewed efforts are required to implement efficient and effective measures to reduce soil loss rates, prevent or minimize further erosion risk, and tailor the required attention to each district’s needs. This will help preserve the ecosystem and ensure the sustainability of agricultural practices and other human activities.

4.3. Evaluating Erosion Susceptibility on Different Slope Gradients

In this subsection, we present the results of our study that examined the rates of soil erosion across different slope gradients. We found that soil loss increased proportionally with increasing rainfall and slope gradients. This could be attributed to the erosive impact of heavy rainfall, which reduces infiltration and increases runoff, leading to soil erosion [108]. Our findings indicate that cropland areas experienced higher soil erosion rates as the slope increased. For instance, croplands with high-intensity rainfall ranging from 693 to 767 mm yr⁻¹ showed steeper slopes (β > 27%) and had higher average soil loss rates. Conversion of land use and its long-term effects on land productivity may have led to the expansion of croplands into upslope regions, resulting in increased soil erosion rates. Unfortunately, upslope areas are more susceptible to erosion due to cropland expansion, lack of vegetation cover, and unsustainable land use management practices. These findings support the conclusions made by Fenta et al. [23]. The results suggest that unsustainable human activities have significantly increased erosion over time. As the population grows and food demands rise, more land is being cleared for crop production, even in unsuitable terrains with steep slopes and protected areas. This results in higher rates of soil erosion and reduced land productivity potentials in sloping landscapes [35,43].

Table 3. Erosion rates per slope gradient and the contribution to annual soil loss, the proportion of cropland areas with average soil loss rates, and mean annual rainfall in 1990, 2000, 2010, and 2020.

Year	Slope Class	Mean Soil Loss (t ha−1 yr−1)	% of Annual Soil Loss	% of Cropland in the Class	Cropland Soil Loss (t ha−1 yr−1)	Mean Rainfall (mm y−1)
1990	I	11	15	62	15	635
	II	22	9	12	40	685
	III	31	12	9	67	714
	IV	43	16	8	110	736
	V	86	49	10	259	764
	Total	17	100	100	5	687
2000	I	30	8	13	40	612
	II	49	12	12	67	644
	III	75	18	11	110	667
	IV	140	52	15	236	693
	V	333	100	100	7	610
	Total	12	10	49	15	549
2010	I	28	8	12	42	677
	II	42	12	11	70	709
	III	62	17	9	113	733
	IV	118	50	12	243	759
	V	25	100	100	6	672
	Total	13	13	55	15	608
2020	I	24	8	12	43	686
	II	38	12	11	76	718
	III	59	18	10	126	741
	IV	106	50	12	258	767
	V	23	100	100	66	660
	Total	10	12	55	15	623

I = 0–7.0%; II = 7.0–11.3%; III = 11.3–17.6%; IV = 17.6–26.8%; V ≥ 26.8%.

illustrates the changes in soil erosion rates across different slope gradients in the four reference periods.

Our comparative analysis revealed that soil erosion rates increased significantly from 1990 to 2000, followed by a slight decrease from 2000 to 2010 and a more significant reduction from 2010 to 2020. For example, the risk of soil erosion in the lower slope grade (β < 7%) decreased by 4% between 1990 and 2000. In particular, there was an increasing trend in average soil erosion for a specific slope class between 2000 and 2010; this trend reversed in the following decade, with a significant reduction in estimated soil erosion between 2010 and 2020. The study also found that areas with steep slopes (β > 26%) and high rainfall intensity (693–767 mm yr⁻¹) had the highest soil erosion rates throughout the study periods. In contrast, lower slope gradients had a low average soil erosion rate. Flat terrains with low-intensity average rainfall ranging from 549 to 634 mm yr⁻¹ accounted for negligible soil erosion during each period (refer to Table 3).

Our research findings are consistent with the conclusions from previous studies conducted in various regions, in which steeper slopes are related to higher soil erosion. In contrast, flat/gentle slopes experience the lowest erosion rates [99,108,109]. For example, a study shows that steep slopes and land use patterns are the most influential factors contributing significantly to susceptibility to higher soil loss rates, which exceed an acceptable threshold [104]. Yesuph and Dagnew [100] also reached the same conclusion in the Ethiopian highlands’ Beshillo Catchment (Blue Nile Basin), where gentle slopes become saturated and have a higher runoff. Another model-based study found that topographies with steep slopes had the most increased soil erosion, and lower average soil loss rates were found on gentle slopes [104]. Rodrigo-Comino et al. [110] similarly determined that sloping Mediterranean agricultural fields with limited protective vegetation and bare soil were associated with the highest rate of runoff and soil erosion. The growth of agriculture and unsustainable human practices, such as deforestation in steep slope areas, significantly increases soil erosion, as highlighted by Karamage et al. [12] and Nambajimana et al. [15].

4.4. A District-Level Investigation of Soil Erosion Risk

Prioritizing districts with the highest risk of soil erosion is a critical first step in guiding future efforts in sustainable land use planning and addressing conservation concerns. By identifying areas with the highest risk, policymakers and land-use planners can focus their efforts on implementing management practices that effectively mitigate erosion and protect the soil while ensuring efficient and effective allocation of resources. This approach can help prevent further soil degradation, protect watersheds, and support sustainable agriculture, thus maximizing the impact of conservation efforts (Figure 10).

Our analysis of erosion risk and severity levels, presented in the supplementary material (Tables S6–S9), reveals that in 1990, no districts had very or extremely high soil erosion risks, and most had low to medium risk levels. Four districts had average erosion rates of 40–80 t ha⁻¹ yr⁻¹, accounting for almost half of the annual soil loss (Figure 10). Nine other districts, (with average soil erosion rates ranging from 22 to 38 t ha⁻¹ yr⁻¹) representing 36% of the study area (namely Kombolcha, Jarso, Haro Maya, Meyu Muleke, Bedeno, Girawa, Kersa, Deder, and Kurfa Chele), had relatively moderate erosion rates. Lastly, districts with low-medium erosion risk had soil loss rates ranging from 10 to 20 t ha⁻¹ yr⁻¹ and were characterized by low rainfall and slope values.

The results of our study indicate an increase in the rate of soil erosion from 1990 to 2000. The Meta district experienced the highest soil loss rate growth, followed by Melka Balo and Deder (Tables S6 and S7). On the contrary, Golo Oda and Meyu Muleke had the highest estimated reduction in soil erosion rate. During the same period (2000), Meta, Deder, and Melka Balo were all exposed to an increased erosion risk, with an average rate of 172, 146, and 120 t ha⁻¹ yr⁻¹, respectively. Our analysis revealed that between 2000 and 2010, soil erosion decreased in fourteen districts, which accounted for approximately 61% of the total area. Notably, Deder showed the most substantial reduction in soil loss, with 6 × 10⁵ t (83 t ha⁻¹ yr⁻¹), followed by Meta with 3.9 × 10⁵ t (54 t ha⁻¹ yr⁻¹), and Goro Gutu with 2.8 × 10⁵ t (61 t ha⁻¹ yr⁻¹). However, some districts recorded higher soil loss rates during the same period. For instance, Golo Oda had an average of 8 t ha⁻¹ yr⁻¹, Meyu Muleke had 3 t ha⁻¹ yr⁻¹, Midega Tola had 4 t ha⁻¹ yr⁻¹ in the south and southwestern parts, and Fedis had 2 t ha⁻¹ yr⁻¹ in the southeast (Figure 10). Between 2010 and 2020, the district of Deder experienced the highest increase in estimated soil erosion, while Golo Oda and Meyu Muleke saw the most significant reductions, shown in Tables S8 and S9.

In 2020, the highest average annual estimated rates of soil loss via water erosion were found in Meta (184 t ha⁻¹ yr⁻¹), Deder (146 t ha⁻¹ yr⁻¹), Goro Gutu (103 t ha⁻¹ yr⁻¹), and Malka Balo (77 t ha⁻¹ yr⁻¹) regions, while the lowest was identified in Chinaksen (16 t ha⁻¹ yr⁻¹), Meyu Muleke (16 t ha⁻¹ yr⁻¹), and Fedis (18 t ha⁻¹ yr⁻¹). The four districts with the highest average soil erosion rates, which comprise 13% of the area and contribute 36% of annual soil loss (as shown in Table S9), are high-priority areas for comprehensive conservation practices due to their higher risk of erosion. The average slope in these districts ranged from 22% in Malka Balo to 30% in Deder and Meta. Focusing on the sections with the highest risk of erosion makes possible sustainable land use planning and prioritizing of the implementation of appropriate SWC practices. Girawa, Kurfa Chele, Jarso, Bedeno, and Kersa, which account for 17% of the annual total soil loss and cover 15% of the study area, have been identified as second-priority areas. The average soil loss rate in these areas ranges from 34 t ha⁻¹ yr⁻¹ in Girawa to 65 t ha⁻¹ yr⁻¹ in Kersa. Aside from the top-priority districts with the highest erosion rates, other neighborhoods also experience water-induced soil erosion, although to a lesser extent. However, their soil erosion rates are still significantly higher than the national soil loss tolerance level [39]. These districts, including Babile, Golo Oda, Kombolcha, Gursum, Midega Tola, Haro Maya, Fedis, Meyu Muleke, and Chinaksen, have average erosion rates ranging from 16–30 t ha⁻¹ yr⁻¹ and account for 72% of the area, making them the third priority for intervention measures.

4.5. Examining the Link between LULC Changes and Soil Erosion

Figure 11a–d displays the spatial extents and proportions of the LULC categories, while Table S10a–d provides the classification accuracy computed for each LULC map. The OA was 93% in 1990 and increased to 95% in 2000, 2010, and 2020, with a K^ coefficient of 0.91, 0.93, 0.93, and 0.94, respectively. When comparing classification accuracy per LULC classes, it was found that PA was 86% for all categories except for built-up (82%) on the 1990 map. Similarly, UA exceeded 89% for all types of LULC except for bare land (84%) on the 1990 map, with the highest accuracy observed for water bodies (UA = 100%; PA = 97–100%). However, higher OE and CE were calculated for the 1990 map. This is partly due to the unavailability of reference field data matching the image acquisition date. On the contrary, there was less confusion between reference points and classified images for 2000, 2010, and 2020, resulting in insignificant classification errors (Table S10). Generally, the LULC classification accuracy for the four-reference periods surpassed the minimum recommended OA of 85% and 70% for each class, as confirmed by our validation analysis for post-classification comparison [79,80].

As shown in Figure 11, the study area was mainly covered by three classes of LULC: bare land, cropland, and shrubland. These classes totaled 89%, 90%, 92%, and 91% for 1990, 2000, 2010, and 2020, respectively. The spatial variations revealed the aerial extent of LULCs attributed to changes between each class during the first (1990–2000), second (2000–2010), and third (2010–2020) periods. For instance, during the first period, croplands increased significantly, going from 12,240.3 km² to 17,018.7 km². However, they significantly decreased to 15,307.9 km² in the second period and 12,518.1 km² in the third. The results suggest an increase in croplands of 278 km² from 1990 to 2020. Similarly, bare land increased from 2148 km² in 1990 to 2377 km² in 2000 and 2534 km² in 2010. However, it decreased to 1474 km² in 2020.

On the other hand, areas covered by built-up land have increased continuously, although their expansion varied between the four reference periods. Additionally, forest areas remained stable during the first period but declined rapidly in the second period and increased slightly in the third. The grassland has undergone some changes in the last three decades (1990–2020), while the waterbody has nearly doubled its initial area. This increase could be attributed to various factors, including the construction of artificial water bodies due to the growing human population and the recovery of Lake Haramaya, which had previously dried up. However, it is also acknowledged that the observed increase could be affected by misclassification due to spectral issues.

The impacts of natural and anthropogenic processes have made LULC change unavoidable, resulting in complex and intertwined effects that alter the characteristics of landscapes [6,43]. Over the past few decades, numerous researchers have reported increasing evidence that population growth and rapid urban development have caused changes in LULC. These changes have prompted the deforestation of natural and semi-natural vegetation covers, transforming these areas into porous lands and increasing the likelihood of soil erosion. Studies have further indicated that water erosion, leading to the displacement of topsoil, is the primary factor that exacerbates land degradation in upland ecosystems [10,11,12]. Tables S11–S13 detail the conversion between the LULC classes from 1990–2000, 2000–2010, and 2010–2020, while Figure 12a–d presents estimated soil loss rates for each LULC type. Figure 13 and Table S14 provide information on changes in soil loss rates attributed to conversions among different LULC classes over the same periods.

We examined soil erosion in different LULC classes in 1990. The results revealed that cropland had the highest average soil erosion (53 t ha⁻¹ yr⁻¹), while bare land, grassland, and shrubland had moderate soil loss rates ranging from 7–39 t ha⁻¹ yr⁻¹. Conversely, forest land had the lowest average soil loss rate (2 t ha⁻¹ yr⁻¹), and built-up areas exhibited resistance to soil loss with an average soil loss rate of 2.4 t ha⁻¹ yr⁻¹. From 1990 to 2000, the average soil erosion rates increased for cropland and grassland areas, while the average soil erosion rates decreased for bare land, built-up areas, forests, and shrublands. LULC change occurred in 49% of the total landscape area, with 17% experiencing an increase in average soil erosion of 80 t ha⁻¹ yr⁻¹ and 15% showing a reduction of 79 t ha⁻¹ yr⁻¹ (Figure 12a).

Similarly, the estimated soil erosion for LULC classes in 2010 and 2020 demonstrate that cropland had the highest average rates, followed by bare land (Figure 11c,d; Figure 12c,d). Results have shown that forested areas, which account for a small share of the overall size, have lower soil erosion rates, suggesting their significant role in regulating soil loss via water erosion [6,47,104]. However, canopy cover can influence the amount of rainfall intercepted before it reaches the soil surface, which affects soil erosion.

We found the highest rates in croplands established by clearing forest areas, compared to areas concurrently changed from cropland areas to forestland. In contrast, converting croplands to built-up areas, forests, shrublands, and grasslands reduced soil erosion by water. These results suggest that converting non-croplands to croplands can increase soil erosion by water (Figure 13a,b). The results of a study conducted in Kenya supported this finding, showing that soil loss was highest when forests and grasslands were altered into farmland [6]. Deforestation and converting natural forests to croplands have also led to soil degradation in Northern Iran [111]. On the contrary, shade-tolerant herb species have effectively moderated soil erosion rates in forest grassland ecosystems in China’s Loess Plateau Yangling, Shaanxi Province [112]. Thus, it is paramount to restrict deforestation practices, implement improved agricultural management, and enact appropriate policy measures to restore degraded forests to mitigate the resulting detrimental ecological impacts.

Consistent with the studies above, our findings have confirmed that cropland expansion and deforestation are the primary drivers of LULC changes, leading to worsening erosion risk. In connection with this, Abera et al. [113] reported that population growth has driven cropland expansion and fuelwood consumption, leading to a change in wood cover in Eastern Africa. They also asserted that rapid vegetation changes could put enormous pressure on the woody cover, highlighting the need for effective management [113]. Moreover, agricultural expansion and climate change have exacerbated the risk of soil erosion. A more recent meta-study has shown that agricultural expansion under climate change causes a significant increase in soil loss rates compared to a reference rate of erosion without any climate or land-use modification scenarios [2]. On the other hand, agricultural abandonment and reforestation have reduced the estimated rates of soil loss [2].

Our findings indicate a reduction in the magnitude of changes in LULC during the second and third periods (Tables S12 and S13), resulting in a decrease in soil loss of approximately 13% and 13.06%, respectively (Figure 13b,c, Table S14). These changes were associated with a decrease in the average erosion rate of landscapes that underwent LULC alteration during the second (4281 km²) and third (3871 km²) periods, with values of 17 t ha⁻¹ yr⁻¹ and 16 t ha⁻¹ yr⁻¹, respectively. However, an increase in average soil loss rate was observed in the remaining areas that underwent LULC change during the second (5655 km²; 114 t ha⁻¹ yr⁻¹) and third (7454 km²; 56 t ha⁻¹ yr⁻¹) periods, mainly due to the expansion of croplands at the expense of other LULC classes. These results contrast with those of Nambajimana et al. [15], who reported the highest soil erosion in landscape areas that changed from grasslands to forestlands. Furthermore, we have shown a decrease in average soil erosion rates for the same LULC class during the periods of 2000–2010 and 2010–2020, with reductions of 0.7 t ha⁻¹ yr⁻¹ and 0.8 t ha⁻¹ yr⁻¹, respectively.

Following the methodology of Ustaoglu and Collier [47], we defined land abandonment as the percentage loss of cropland areas that were neither converted into urban lands nor afforested during the first, second, and third ten-year periods. The repercussions of land abandonment on biodiversity, ecosystem services, and agricultural soil could have either positive or negative outcomes, depending on local or regional factors such as the degree of vegetation recovery, conservation efforts, and topographical and agroclimatic characteristics [110,114]. Based on the earlier approach, we extracted and examined the soil loss rates accompanying different transformations, including cropland expansion, abandoned cropland surfaces, afforested areas, and deforested areas. We investigated the extent to which agricultural development and deforestation have resulted in increased soil loss. In contrast, agricultural abandonment and reforestation have reduced soil erosion in the study area (Figure 14a,b).

Analysis of land use change shows that during the first, second, and third periods, a total of 2399 km² (10%), 4363 km² (17%), and 5604 km² (22%) of cropland were abandoned, respectively. This reduced the average soil loss rates of 23, 21 t ha⁻¹ yr⁻¹, and 46 t ha⁻¹ yr⁻¹ (Figure 14). Conversely, during the same periods, a total of 7694 km² (31%), 3271 km² (13%), and 3299 km² (13%) of cropland expansion occurred, respectively, resulting in an increase in the average soil loss rates by 36 t ha⁻¹ yr⁻¹, 70 t ha⁻¹ yr⁻¹, and 36 t ha⁻¹ yr⁻¹. Furthermore, the results indicate a moderate increase in soil erosion in deforested areas, highlighting the effectiveness of afforestation measures in reducing soil loss over time.

The forested land area decreased from 1244 km² of the total landscape area between 1990 and 2000 to 952 km² between 2000 and 2010. It then rose again to 1090 km² from 2010 to 2020. As a result, there was a corresponding reduction in average soil erosion estimates of 26 t ha⁻¹ yr⁻¹, 23 t ha⁻¹ yr⁻¹, and 26 t ha⁻¹ yr⁻¹, respectively. On the other hand, during the first, second, and third evaluation periods, 1228 km² (5%), 1604 km² (6%), and 880 km² (4%) of the total land area experienced deforestation, which increased the average soil loss rates by 23.7 t ha⁻¹ yr⁻¹, 23.6 t ha⁻¹ yr⁻¹, and 15.9 t ha⁻¹ yr⁻¹, respectively.

Therefore, it is clear that vegetation growth as a soil erosion control measure is critical for future planning and implementation to mitigate water-caused soil erosion. Although the study reported in this article evaluated the impact of spatiotemporal LULC change on erosion risk, it is essential to conduct future exhaustive investigations to bridge the gap on the adverse consequences of the change process on the ecosystem and ecosystem services. This evidence is critical for forming development policy and improving land management planning.

4.6. Exploring Diverse Approaches to Combat Soil Erosion

Soil erosion caused by water is a primary environmental concern, posing a severe threat to soil and water resources. Various stakeholders have actively participated in SWC programs, including government, development partners, and particularly local communities, through consecutive free labor mobilization events to restore degraded land and improve productivity. Although the effectiveness of these efforts in reducing soil erosion in cropland areas has not been thoroughly investigated in the present study area, studies have shown both challenges and ineffectiveness [94,107,114], as well as positive contributions in some cases [102], related to restoration initiatives and the implementation of SWC practices. This highlights the need for exhaustive local-specific investigations to identify critical erosion-prone areas and prioritize successful management approaches, thereby developing effective measures to address degradation caused by water erosion, ensuring the long-term sustainability of ecosystems and their services [18,26].

This study examines and quantifies the effectiveness of established conservation support in reducing soil erosion rates in cropland areas of Ethiopia’s EHZ, as shown in Table S5. Our objective was to evaluate the efficacy of these support strategies by integrating them into the soil erosion prediction model and answering the following questions: (i) To what extent can the application of the three support practices reduce the average soil erosion in cropland areas? (ii) Which of the conservation support practices is the most effective? (iii) Do the soil loss rates in the current study area align with estimates of cropland erosion rates in other regions when proposed practices are implemented?

Our analysis based on the predicted soil loss rate of the model showed that the cropland areas accounted for 50% of the total LULC (Figure 11d), with an estimated average soil loss rate of 66 t ha⁻¹ yr⁻¹ (Figure 15a), representing 52% of the annual soil loss. The model estimate is significantly higher than the Ethiopian Highland Reclamation Study [51], which modeled the average yearly soil loss predicted for cropland fields at 42.1 t ha⁻¹ yr⁻¹. Our estimate is also nearly twice the national average soil erosion rate for croplands, covering 26 × 10⁶ ha [23]. Other studies in Ethiopia (e.g., [34,100]) and various areas worldwide, such as those studies conducted by Watene et al. [11], Borrelli et al. [20], and Karamage et al. [22], have predicted higher soil loss rates in cropland areas compared to other classes of LULC. Our baseline soil erosion estimate indicates that entire croplands (12,518 km²) are at significant risk of soil erosion, with 25% of them experiencing unsustainable soil loss rates (>1 t ha⁻¹ yr⁻¹). However, nearly 14% of the cropland area showed soil loss rates that exceeded the tolerance limit of 10 t ha⁻¹ yr⁻¹ [39]. These areas were mainly observed on steep slopes (β = 22–44%) with a high rainfall intensity ranging from 741 to 1227 mm yr⁻¹ (Table 3). Thus, it is vital to prioritize soil conservation measures in these high-risk cropland areas to prevent further degradation.

We examined three management support practices, namely terracing, strip-cropping, and contouring, and compared their effectiveness to a baseline scenario. The results showed that the three support practices reduced the estimated soil erosion rates more effectively than the baseline scenario. Nevertheless, the extent of the reduction varied among the conservation scenarios considered in the model. Specifically, terracing was the most effective practice, potentially reducing average soil erosion rates by 87% to 9 t ha⁻¹ yr⁻¹ (Figure 15b). On the other hand, strip farming decreased average soil erosion rates by 65% from baseline to 23 t ha⁻¹ yr⁻¹. Contouring was also effective, reducing the average soil erosion by 29% compared to the baseline, resulting in a corresponding soil erosion rate of 47 t ha⁻¹ yr⁻¹. The estimated rate based on contours was lower than the baseline average soil erosion rate by 19 t ha⁻¹ yr⁻¹. However, it was higher than the estimated terracing and strip-cropping rates by 38 t ha⁻¹ yr⁻¹ (82%) and 23 t ha⁻¹ yr⁻¹ (50%), respectively. In summary, our study underscores the potential of different SWC practices in reducing cropland soil loss.

As we lacked field data, we compared our results to previous studies to verify the consistency of our model in predicting soil erosion in tropical African croplands under similar management practices. Our findings align with recent work, which demonstrated that implementing terraces, strip-cropping, and contouring could reduce average erosion rates compared to the baseline scenario by 84.4%, 60%, and 24.44%, respectively [11]. In another study conducted in the Nyabarongo River Catchment in Rwanda, Karamage et al. [12] consistently predicted soil erosion using the RUSLE model. They showed that implementing terraces, strip-cropping, and contouring could reduce average annual soil erosion by 7 t ha⁻¹ yr⁻¹, 18 t ha⁻¹ yr⁻¹, and 35 t ha⁻¹ yr⁻¹, respectively. However, terracing has the potential to reduce the average cropland soil loss by up to 20%, while contouring and strip cropping can intensify soil loss by two to four times [12]. In contrast, our study found that these practices could effectively reduce soil erosion in cropland areas (as shown in Figure 15b).

Furthermore, our analysis supported a study conducted in Uganda’s cropland areas. It showed that terraces and strip-cropping could reduce predicted soil erosion by 80% and 47% below the sustainable soil loss tolerance value for SWC [22]. However, contour-based conservation practices could only maintain an average soil erosion rate similar to the baseline estimate of 1.5 t ha⁻¹ yr⁻¹ in cropland areas, with a moderate average soil erosion estimate of 1.6 t ha⁻¹ yr⁻¹ [22]. Likewise, a study from the Finchaa catchment of the Upper Blue Nile Basin confirmed that conservation scenarios help reduce sediment yields [26]. The research showed that contour strips and soil bunds could decrease sediment yield by 64% and 65%, respectively, bringing it below the acceptable level. However, practices such as contouring and terracing that increased sediment yield beyond the tolerable point did not significantly decrease the amount of suspended sediment in the study catchment [26]. For their part, Ligonja and Shrestha [32] found that land management practices (e.g., tree plantation, cutoff drains, drainage channels, contour ridges and ridges, mixed cropping, strip cropping, crop diversification, livestock prohibition, and zero grazing) could help reduce soil erosion in the eroded Kondoa area in Tanzania.

To effectively manage SWC, it is essential to have clear goals that provide direction, purpose, and a framework for decision-making and strategic management. They define conservation objectives, guide resource prioritization and allocation, inform planning and implementation of successful conservation practices, facilitate monitoring and evaluation of conservation initiatives, foster stakeholder engagement, and support advocacy for policies promoting sustainable practices. Common conservation management goals include reducing the risk of soil erosion, preserving soil fertility, and improving water quality while maintaining ecological sustainability and biodiversity, which is critical for human well-being and providing ecosystem services. With such goals in position, land managers and conservation planners can develop focused strategies, measure progress, and achieve meaningful improvements in soil and water resources. These measures can also help adapt to climate change [60,106].

Many countries have implemented conservation practices that can serve as models for other regions. For example, the United States Department of Agriculture (USDA) Environmental Quality Incentives Program financially and technically assists farmers and ranchers in adopting conservation practices such as no-till agriculture, cover crops, and conservation tillage to reduce soil erosion and improve soil health [115]. A recent study by Deng et al. [116] found that China’s Sloping Land Conversion Program has decreased soil erosion rates and increased forest coverage on steep slopes. Nyberg et al.’s study found that adopting sustainable agricultural land management (SALM) practices, including agroforestry, led to improved maize yield, higher food self-sufficiency, and significant monetary savings for project farmers compared to control farms in Kenya [117]. These practices showcase the effectiveness of targeted conservation efforts and can serve as models for achieving conservation goals.

However, the most effective conservation practices can vary depending on local conditions [8,39]. Accordingly, conducting a detailed cost–benefit analysis is necessary to establish the best land management practices (BLMP) or “best bet” strategies in a specific region to optimize investment returns on SWC [8,26,118]. In this context, the results of our study can contribute to a broader understanding of land susceptibility to water erosion and its severity level, particularly in croplands. These findings can help inform policy decisions on appropriate land use planning and coordinated intervention measures to restore degraded land and maintain the sustainability of the ecosystem.

5. Conclusions and Perspectives

In this study, we estimated water-induced soil erosion rates, assessed the severity of erosion risk, and identified potential management practices to mitigate cropland soil loss in the Ethiopian EHZ. We used multi-source remote sensing data analysis and a GIS-based RUSLE model. The average annual soil erosion rates were estimated to be 33 t ha⁻¹ yr⁻¹ in 1990, 50 t ha⁻¹ yr⁻¹ in 2000, 44 t ha⁻¹ yr⁻¹ in 2010, and 39 t ha⁻¹ yr⁻¹ in 2020. The results suggest an 18% increase in overall soil loss over the last 30 years (1990–2020), particularly in districts with high erosion factors. The analysis showed the significance of land use, topography, and rainfall in influencing soil erosion. The findings highlight the need for targeted SWC planning and management strategies to curtail topsoil loss and rigorous erosion control intervention measures in highly eroded districts.

Our findings showed that applying well-planned structural and agroeconomic management support practices, such as terracing, stripping, and contouring, could reduce soil loss in cropland areas. Our study’s contribution to the knowledge base on soil erosion and facilitating operational conservation planning is significant. Policymakers, land managers, and conservation practitioners can use the erosion risk maps in this study as evidence-based spatial decision support for designing and implementing more effective measures, thus mitigating the adverse impacts of water-induced soil loss, enhancing agricultural productivity, and promoting ecosystem resilience. We, however, bring to light the following shortcomings:

• First, it is critical to note that independent sets of actual field measurements and soil loss monitoring at specific sites should be used to validate RUSLE-based estimates to ensure the model’s accuracy and better understand its performance, replicability, and relevance to policy decisions. However, validating model-derived soil erosion with on-site measurements is often challenging due to a lack of field observations. This rings true in the studied landscape, a field-data-scarce zone in eastern Ethiopia. Thus, when applying the model to broad regional studies encompassing unobserved areas, it is decisive to recognize that uncertainties may stem from constraints or the absence of access to high-quality datasets for accurately determining model sub-components.
• Second, we have incorporated values aligned with RUSLE factors, such as C and P, adopted from published coefficients. Interestingly, given the development of ranges of alternative empirical methods for generating the model factors, we propose examining the workability of various model parameterization approaches in the context of the study area. Furthermore, one must consider the uncertainties prevalent in the RUSLE model factors that depend on data quality, expert decisions, and assumptions made during the modeling process, which may introduce errors/biases in sol loss estimates.
• Third, the study does not extensively explore the socioeconomic elements that could potentially impact decisions about land use and SWC methods. Integrating the socioeconomic factors and human interventions that contribute to soil erosion would provide a more comprehensive knowledge of the determinants of soil erosion. This could involve surveying farmers, analyzing land use policies, and exploring the socio-economic dynamics of the region.
• Lastly, it is essential to recall that this study only addressed rill and inter-rill erosion given that, based on the information provided, the RUSLE model applied to estimate soil loss has limitations in considering all forms of erosion, such as gullies. Therefore, future studies should pay more attention to gully erosion, which contributes to land degradation in the landscape studied by causing loss of nutrients and sediment accumulation in downstream areas.