1. Introduction
Urbanization is a continuous process that has boldly accelerated with population increase, expansion and spread of built-up structures in a designated urban area [1]. Urbanization in Africa has been growing at alarming rates with an anticipated annual growth rate of approximately 3.9% [2]. Population growth and increasing economic activities have been linked to aggravate Land Use and Land Cover (LULC) changes [3,4,5]. As such, urban development inevitably involves construction and sealing activities that alter natural landscapes [6] resulting in an increase in impervious surfaces, which replace natural vegetation and reduce the capacity for water infiltration.
This in-turn results in surface runoff that substantially threatens soil loss in vulnerable landscapes through erosion processes [7]. While multiple studies on urbanization processes focus on the social and planning aspect and the assessment of the effects of urbanization processes on the (quasi-) natural environment, in the presented study we assess the effects of urbanization on surface processes by water in the strongly urbanizing area of Harare in tropical Africa; as case study we selected the 35 km2 largest growing informal urban district of the Harare Metropolitan Province, which showed urbanization processes in the sense of surface-sealing of 18 km2 since the year 2000 [8,9]. Due to the unprecedented growing rates of urban-built up area in Harare Metropolitan Province [9,10,11], there is profound need for mapping potential erosion, and estimating potential soil erosion risk.
Soil erosion by water is a naturally occurring process, which is accelerated by human interference and activities including agriculture, deforestation, urban expansion and potentially climate change [12,13,14,15]. It is a continuous process that involves the bodily detachment, transport and deposition of particles originating from soils and weathered bedrock [16,17]. This resultant effect of overland flow and surface runoff has degrading effects on soil resources and affects agriculture and infrastructure [18,19]. The mechanisms involved in soil erosion by water underly a spatio-temporal variability and they are affected by factors such as rainfall, soil characteristics, ground cover and terrain [20]. In regions with seasonal rainfall dynamics, vegetation is sparsely distributed during dry season, making surfaces vulnerable to high raindrop energy impacts at the onset of wet season [21,22,23,24]. Accelerated soil erosion is also reported especially on rainfed agriculture during the pre-monsoon due to high intensity rainfall events occurring concurrently with high wind speeds [25,26].
The spatial heterogeneity and dynamic mechanisms of soil erosion can be attributed to on-site and off-site effects that trigger landscape degradation [18,27]. On-site soil erosion impacts encompass soil sealing and reduced crop productivity [28]. Furthermore, change of soil structure emanating from the loss of topsoil and soil aggregates due to high rainfall erosivity and low infiltration rates severely perturb the sustainability of agricultural systems [24,29]. Off-site damages effect widespread damages such as clogging of canals and drainage systems, siltation of reservoirs and destruction of roads thereby impeding transportation of goods and accessibility to properties [23,27,30]. Construction activities such as grading and use of heavy machinery on steep and lengthy slope are reported as major sources of sediment loss due to their effect in accentuating erosion potential [31]. Long-term LULC changes negatively impact river systems through the alteration of channel flow, soil deposition processes, soil textural organization and habitat loss [25,32]. Henceforth, the need is to assess the spatial distribution and extent of soil erosion in the rapidly transforming Epworth district of the Harare Metropolitan Province.
Globally, almost 84% of land loss results from soil erosion processes [30]. The estimated mean rates of soil erosion across the world range between 12 and 15 t ha−1 yr −1 [33]. In Africa, it is estimated that 19% of the total reservoirs are under siltation threat due to soil erosion by water [33]. About 494 million hectares of land in Africa are subjected to different types of degradation with degradation influenced by water estimated at 227 million ha [34]. A countrywide annual soil loss of 1.5 * 109 t was predicted for Ethiopia, with cultivated lands recording the highest soil loss rate at 42 t ha−1 yr −1 [35]. Soil erosion rates in Ethiopia’s cultivated lands and highlands have been highly influenced by increasing demographic pressure, climate change, terrain and depletion of vegetative cover [36,37]. Only a few studies focus on an assessment of soil erosion risk in urbanized areas. For Kinshasa/DR Kongo, approximately 4.3% of the total urban area was predicted under high-risk of soil erosion with over 15 t ha−1 yr −1 annual average loss [38]. A study on Kuala Lumpar/Malaysia Metropolitan city shows that approximately 38.1% of the city area was above the soil loss tolerable rates (> 1 t ha−1 yr −1) due to high soil erodibility and steep slopes in some parts of the city [3]. At all, urban flooding falls among other problems associated with water induced urban soil erosion. This further aggravates risks on downstream water quality from overflowing sewers and street solids transported in concentrated urban runoff [39].
A wide range of empirical, conceptual and physical models have been developed to estimate soil loss risks and these models vary in complexity, data requirements, processes considered and calibration [40,41]. These models include, among others, the statistical model of the Universal Soil Loss Equation [42] and its derivatives the Revised Universal Soil Loss Equation [43], the Soil and Water Assessment Tool [44] or the Water Erosion Prediction Project [45]. In general, model selection is particularly dependent on the availability of data, attributes of a working area and intended use [41,43]. However, [23] we point out that the complexity of urban land use makes standard soil erosion recording and sampling techniques difficult to apply.
The suitability of satellite remote sensing and Geographic Information Systems (GIS) applications to extract, delineate and manipulate land characteristics, and their integration with the Revised Universal Soil Loss Equation (RUSLE), makes them fundamental tools for spatial soil erosion estimation [46,47,48,49]. The RUSLE model is widely used and has been validated over decades, in addition to the fact that its limitations have been well documented [43,50]. Various studies on the prediction and assessment of soil erosion rates using RUSLE have been reported for Africa, with much attention to highlands and river watersheds [51,52,53]. However, there is paucity in research on mapping potential erosion and soil erosion risk modelling in urban built-up areas [23,30,38]. The advances in urban development and soil erosion management have shown that land managers and policy makers consider the spatial distribution of soil erosion risk more than actual soil loss values [3].
The objectives of this paper are: (a) to model the spatial distribution of the soil erosion risk and potential erosion for Epworth district over different time slices: 2000 and 2018 in order to assess the dynamics of soil erosion over time in heterogenous urban landscapes; (b) to assess the influence of rampart land use and land cover changes on soil erosion risk in the Epworth district through the analysis of temporal soil erosion estimations between 2000 and 2018. The temporal investigation makes use of the RUSLE model for estimating quantitative and spatial data on potential erosion and soil erosion risk in the Epworth district of the Harare Metropolitan province [54] in order to enhance and support soil conservation planning [43]. The study implements field observations, direct erosion feature registration and quantification towards establishing comparison and validation tools for Epworth district soil erosion models [55]. Overall, our study results aim to provide scientific advice towards the sustainability of urban growth.
2. Materials and Methods
2.1. Study Area
Harare Metropolitan Province is the capital city of Zimbabwe (Figure 1), located between 17°40′–18°00′ S and 30°55′–31°15′ E. Harare Metropolitan Province is situated in the headwater areas of the Mukuvisi, Marimba and Manyame rivers, which supply Lake Chivero, a reservoir that supplies Harare with water [39]. The Harare Metropolitan Province is composed of four districts: Harare urban, Harare rural, Chitungwiza and Epworth. Epworth district is located about 12 km southeast of the central business district (Figure 1). Epworth district has been purposively selected as a study site because it is the largest informal settlement across urban districts in Zimbabwe [8]. Epworth counted approximately 500 families in 1950 [56] while in 2012 the district had an estimated total population of 167,462 [57].
Settlement activities such as urban agriculture, construction, illegal sand mining, brick molding and effluent discharge from industries pose water quality threats to the Lake Chivero drainage basin [8,58]. Epworth district is largely dominated by high density residential areas. These are characterized by overcrowdings and concentrated housing residential structures due to densification; only little vegetation occurs across the settlements [8].
The general topography of Epworth district is characterized by undulating and slightly rolling terrain in the upland areas, being part of the southern Africa Highveld. Elevations range from 1455 m to 1556 m a.s.l. For Epworth district, clayey Fersiallistic soils (moderately leached soils of the kaolinitic order) occur in contact zones and Paraferrallitic soils (comprises of highly leached soils) are widely distributed across the entire district [59,60]. The sub-tropical climate of Harare Metropolitan Province is dominated by four distinct seasons: the cool-dry season (mid-May to August), hot-dry season (September to mid-November), rain-wet season (mid-November to mid-March) and the post rainy season (mid-March to mid-May) [11]. Harare Metropolitan Province receives 470–1350 mm of precipitation annually; rainfall predominantly occurs during the four months of rainy season. Average temperatures range from 7 °C to 20 °C during dry periods and from 13 °C to 28 °C in hot-dry periods [11].
2.2. Soil Erosion Modelling
2.2.1. Parameter Estimation for Soil Erosion Risk Assessment Using RUSLE
The RUSLE model is simple and the mostly used computerized version of the Universal Soil Loss Equation (USLE), a statistical model developed to estimate the annual soil loss per unit area based on erosion factors [43,61]. The RUSLE model has been widely implemented for the prediction of average annual soil losses caused by sheet and rill erosion and to display the spatial distribution of potential erosion risk [23,55,61,62,63,64,65]. The application of the RUSLE model for soil erosion risk considers the rainfall erosivity factor (R), soil erodibility factor (K), slope length and steepness factor (LS), land cover and management factor (C) and the support practice factor (P) [43]. In the current study, the RUSLE model was adapted for mapping potential erosion using C and P factors as identity elements (C and P = 1) and for the spatial distribution of soil erosion risk.
According to [43], the Revised Universal Soil Loss Equation (RUSLE) states that:
(1)
where: A is the annual average of soil erosion rate factor (t ha−1 yr−1); R is the rainfall erosivity factor (MJ mm ha−1 h−1 yr−1); K is the soil erodibility factor (t hMJ−1 mm−1); LS is the dimensionless slope length and steepness factor; C is the dimensionless crop management factor (ranging between 0 and 1) and P is the dimensionless conservation support practice factor (ranging between 0 and 1). The calculation of the potential erosion is based on the same formula while adjusting factors C and P to one [3,65].Spatial analyses in the study were performed using ArcGIS 10.2 software to assess the dynamics of soil erosion over time in the heterogenous urban landscapes for the years 2000, 2018 and the overall, 1984–2018 long-term rainfall data were used for rainfall erosivity factor computation (Figure 2). The LS factor RUSLE geospatial input factor was computed using a hydrology module (LS-factor field base) in SAGA GIS [66]. The acquired geospatial input parameters for the RUSLE model (Table 1) were used to produce thematic maps for the estimation of potential erosion and soil erosion risk generated within every cell grid [67,68,69]. For data harmonization, we resampled all data sources to determine 30 m × 30 m grid cell size using nearest neighborhood technique so as to retain original pixel value before carrying out grid cell calculations. This was performed to enhance data compatibility from varying data sources used for modelling [70].
2.2.2. Data Integration for Soil Erosion Risk Assessment
Rainfall Erosivity Factor (R)
Rainfall erosivity factor (R) describes the ability of rainfall to trigger soil erosion [69,73,74]. The RUSLE model makes use of this erosivity factor (R [MJ mm ha−1 h−1 yr−1]) in integrating the effects of raindrop impact, rainfall duration and resulting runoff rates, which are coupled with the amount and the energy within each recorded rainfall pattern [43,69]. Rainfall erosivity was calculated using mean annual rainfall data collected from three gauging stations in Harare Metropolitan Province (Harare Belvedere, Harare airport and Harare Kutsaga; Table 2) following Equation (2) [40,64]:
(2)
where R = Rainfall erosivity factor, M = Mean annual rainfall.To analyze possibly changing rainfall erosivity since 1984, annual rainfall averages were calculated for the time steps: 1984–2000, 2000–2018 and overarching, 1984–2018. The mean annual precipitation data were interpolated over entire district by applying the inverse distance weighting (IDW) interpolation technique and converted to rainfall erosivity by applying Equation (2).
Soil Erodibility Factor (K)
The responsive effect of a particular soil in a given location to the erosive power of rainfall and runoff impacts is referred to as the soil erodibility factor (K) [16,75]. Soil erodibility is regarded as a function of the soil texture, organic matter content, soil structure and the degree of permeability [41,76]. Soils being highly susceptible to erosion have soil erodibility values close to 1, whereas corresponding values close to 0 indicate a resistive nature of the soil [51,77]. In the current study, information on soil structure and profile permeability was not available. Therefore, the K factor was estimated using the ISRIC (International Soil Reference Information Centre)-World Soil Information data [78], following the equation by [79].
Topographic Factor (LS)
The LS factor summarizes the effects of topography on soil erosion and combines the influence of slope length and slope angle on soil loss; while the S-factor measures the effect of slope steepness, the L-factor defines the impact of slope length. The slope length L is defined as the distance between the upslope starting point of a slope segment to the downslope point where deposition begins [33,66]. Increasing slope length and steepness per unit area results in increased runoff and flow velocity and consequently in increased soil loss exposure [42,75,80]. The Shuttle Radar Topography Mission (SRTM) digital elevation model (DEM) with a spatial resolution of 30 m was used for the LS factor computation applying SAGA software 2.3/Hydrology module [66,80]. For data pre-processing, fill sink algorithm was employed using the spill elevation method for filling sinks on the DEM [81]. Multiple flow direction tool (MFD) incorporated in the Hydrology module was applied to the DEM to assign flow directions and flow accumulation [82,83].
Land Cover and Management Factor (C)
The Land cover and management factor C’s pivotal role is to capture differences in soil loss in vegetated areas by dissipating raindrop impact on the soil surface compared to bare areas [43,84,85]. The C factor decreases from 1 to 0 depending on vegetation cover and cropping management systems implemented to mitigate soil erosion [41,80]. In order to determine the C factor, land use maps generated from the supervised classification of satellite images were adapted for use [9]. The C factor values for each land use are a result of weighted average of the soil loss ratio deduced from a reference plot (bare) with a C factor of 1 [43,66]. Henceforth, the C factor values of each land use type were evaluated from literature [3,86,87,88]; weighting of the data was performed following field observations and the biophysical characteristics per sampling plot (canopy cover, prior land uses, and surface cover) [43,89]. The weighted C factor evaluation considered plant growth, height and the extent of canopy cover in-situ [89,90], in relation to bare and sealed area. Further, previous urban farming practices and residues of the plant material influence were majorly factored in during surface cover subfactor assessment per sampling plot. For croplands the C factor was estimated by averaging the values of predominant crops within the study area’s croplands [91]. The availability of remote sensing data, especially the supervised classification maps, appropriately aided the evaluation of spatial variability of the C factor [89]. Overall, the weighted C factors were estimated as a result of multiplying the scaled percentages of the evaluated C subfactors (Table 3), reviewed RUSLE C factors according to the literature [86,87,88] and the ratio of sealed area proportion [9] in relation to the reference plot (bare).
Support Practice Factor (P)
The Support practice factor P expresses the effects of surface management practices that are applied to reduce soil loss through erosion processes [30,43]. These practices include among others terracing, strip cropping and contour ploughing [43]. The P factor value ranges between 0 and 1, where 0 shows the highest effectiveness of the conservation practice and 1 indicates that there are no support practices or measures implemented [43,92]. In this study a P factor of 1 has been all-over applied due to the area wide absence of support or management practices.
2.3. Mapping and Surveying Soil Erosion Dynamics
The complexity of the setup of urban built-up areas and the distribution of the different land use require soil erosion field survey mapping to receive reference data [23]. In this regard, a simple snap-shot sampling procedure was implemented during the field survey in December 2019 to estimate the spatial extent of eroded areas on plots 40 m × 80 m in size [23]. In total, 49 sites were randomly surveyed, covering varying land use. Mapping of spatial soil erosion assessment was geocoded using a hand-held GPS (Garmin 60Cx); parallel on-site soil erosion features were measured, considering the erosion features (inter-rills and rills) individually calculating their area and volumes [55,93]. Rills were defined as linear erosion channels not more than 0.5 m deep and with a cross-sectional area <929 cm2 to make them distinct from ephemeral gullies and deep incised gully features [94,95].
Overall total percentages of disturbed surface area data were registered for the sample plots in the Epworth district to estimate the spatial extent of perturbed regions of the district. The spatial extent of eroded areas was estimated as percentage per plot [23] and classified into five severity classes. Mapping results are displayed in a positional diagram map.
3. Results
3.1. Factors Controlling Soil Erosion
The rainfall erosivity factor map for Harare’s Epworth district depicts very small variations over the study periods between 1984 and 2018 (Figure 3). The annual area wide R factor averaged 329 MJ mm ha−1 h−1 yr−1 between 1984 and 2000 (std = 0.55, n=16), 315 MJ mm ha−1 h−1 yr−1 between 2000 and 2018 (std = 0.27, n=18) (Figure S1). As for both observation periods, the R factor did not vary significantly (α >0.05), the area wide averaged R factor for the time period 1984–2018 with 322 MJ mm ha−1 h−1 yr−1 (std = 0.38, n=34) was applied (Figure 3).
Soils in the Epworth district upland’s areas are dominated by sandy to clayey loams and sandy loam along the alluvial plains; also, at sloping positions sandy loam soils frequently appear. Correspondingly, K factor values vary in the Epworth district between 0.06 and 0.09 t h MJ−1 mm−1 and are highest along the slopes flanking the valleys while they are lowest in the plateau area. Relief in the Epworth district is due to its location at the northern extension of the Highveld slightly rolling. Slopes at the plateau vary between 0.74° and 2.0°, while along the valley flanks they increase in steepness up to 4.6°. Correspondingly, the values of the topographic factor (LS) are highest along the valley flanks, increasing from the headwater areas moving downstream of the Epworth district (Figure 3). Overall, the topographic factor values range from 0 in the plateau areas to 21.74 in steep areas.
The land cover management factor weighted values are based on field observations and satellite images on LULC analysis (Table 1 and Table 4, Figure 4). Data for the LULC changes from 1984–2000 and 2000–2018 were adapted from [9].
3.2. Potential Erosion Risk Analysis
The potential erosion risk map was derived from the application of the “natural” RUSLE factors for soil characters, rainfall erosivity and topography (K, R and LS) for the Epworth district. Resulting data are classified into five potential erosion risk classes showing how erosion varies in the Epworth district (Figure S2). The potential erosion risk map indicates the vulnerability of the landscape independent from vegetation cover and crop management. The findings reveal that very high to extreme erosion risk areas occur in areas of steep slopes (Figure S2); the only spatially slight variations of R and K factors cause a strong control of potential erosion risk by topographic factors’ LS. Due to the general orientation of the drainage network to the south and downstream with increasing inclination of the slopes along the valley flanks, the areas with high erosion risk expand from the north to south. Very high to extreme erosion risk areas are also observed in the eastern parts of the Epworth district resulting from locally occurring steep slopes towards the middle course of the Jacha river and tributaries. In contrast, low to moderate potential erosion risk areas occur in the plateau areas of the Highveld. Extreme potential erosion risk zones occur in the immediate vicinity of streams and at steep slopes. The area wide average potential erosion risk rate in the Epworth district was 13.2 t ha−1 yr−1 since 1984, referring to the precipitation data 1984–2018 (Figure S2, Table 5).
3.3. Soil Erosion Risk
The estimated soil erosion risk maps were generated for 2000 and 2018. The estimated soil erosion risk averaged 1.31 t ha−1 yr−1 in 2000 and 1.12 t ha−1 yr−1 in 2018; with highest total soil loss rate for the Epworth district amounting to 92.79 t ha−1 yr−1 in 2000. The spatial patterns of the estimated soil erosion risk indicate areas with high soil erosion loss predominantly along the river courses (Figure 5). Correspondingly, areas of high soil erosion risk can also be found in the southwestern and southeastern parts of the Epworth district with highest soil erosion risk close to the river courses.
The soil erosion damage map shows plot-wise surface area damages in the Epworth district calculated from the soil erosion feature dimensions, expressed in percentages of the total surface area for each plot (Figure 5). According to this soil erosion damage map, the southeastern and southwestern plots of the Epworth district experience higher soil erosion compared to the other areas investigated. The spatial extent of the soil eroded area ranges from 0% to 1.4% for the surveyed plots in the Epworth district. The magnitude of soil erosion observed during field mapping in 2019 indicated that slope, high proportion of sealed and impervious surfaces attributed to increased soil erosion damages in the Epworth district (Figure 6).
The estimated soil erosion risk for the year 2000 highlights that 56.3% of the Epworth district was exposed to low soil erosion risk and 25.9% to moderate soil erosion risk, while 15% of the Epworth district was exposed to high, and 2.8% to very high and extreme soil erosion (Figure 6, Table 6). For 2018 modelling of soil erosion risk displays a slight decline of risk of exposure with 59.5% of the area being exposed to low soil erosion risk and 29.3% to moderate soil erosion risk; the spatial extent of areas exposed to high soil erosion risk declined to 10% and areas exposed to very high to extreme soil erosion risk covered 1.2% of the Epworth district (Table 6).
3.4. Magnitude of Soil Erosion in Epworth District
A spatial extent of about 765 m2 was eroded with an average area of 31 m2 affected by soil erosion as calculated from the 49 randomly selected sample plots in Epworth district during the field survey in 2019 (Table 7). The soil erosion damage measured approximately 0.5% of the total area mapped (15.7 ha). The occurrence of soil erosion features varied in the surveyed plots corresponding to vegetation cover, slope characteristics and human activities.
Model validation was done using the empirical RUSLE model data in comparison with on-site field measurements. The results showed good RUSLE model performance as there was satisfactory moderate positive correlation between field measurements and model results for sample areas (r = 0.76 and R2 = 0.581, p <0.05) (Figure 7). This provides confidence in the application of the model for sustainable land use planning and decision-making processes.
3.5. Land Use and Soil Loss Analysis
The results show that about 50,408 tons of soil were estimated to be lost under 2000 LULC conditions, while an estimated total soil loss of 42,934 tons was calculated for 2018 (Table 8). For the industrial areas of Epworth district, approximately 40 tons of soil loss were estimated for 2000, while an increase of up to 47 tons of soil loss was estimated for the same land use type for 2018. For 2000, for the land use type “less concentrated residential area” (15.5% of the Epworth district in 2000) 6218 tons of soil loss were estimated while, for the land use type “concentrated residential areas” (38% of the Epworth district in 2000) about 14,018 tons total soil loss were estimated. An increase in soil erosion risk for less concentrated and concentrated residential areas were estimated to amount 12,203 tons for the “less concentrated residential areas “(31.5% of the Epworth district in 2018) and 19,858 tons for the “concentrated residential areas” (52.6% of the Epworth district in 2018). A decline in the estimated soil loss was observed for land use types either of agricultural use or covered by green spaces (undifferentiated) between 2000 and 2018, decreasing proportional to the reduction of the areas of these land use types (Table 8).
4. Discussion
For the tropics, studies reported average soil loss rates of 5 t ha−1 yr−1 [96,97], while [98] highlights that a soil loss limit of 11 t ha−1 yr−1 may be accepted as reasonable mean annual loss due to soil erosion. However, [99] argues that for sensitive and fragile land areas average soil loss tolerance rate of 2 t ha−1 yr−1 could be recommended. In contrast, considering the slow rate of soil formation and spatio-temporal effects of soil loss on water quality and productivity, tolerance limits for soil erosion loss are set for the tropics at 1 t ha−1 yr−1 [3,65,100]. The occurrence of soil loss exceeding 1 t ha−1 yr−1 was considered as the critical rate for the Epworth district due to the low rate of soil formation as typical for the tropics [3]. This agrees with the recommendation of [98] that the rate of soil loss through soil erosion should be relatively balanced with the soil formation rate to minimize excessive environmental damage. As such, average soil loss rates exceeding the suggested 1 t ha−1 yr−1 are classified as unsustainable to continue supporting land use [101,102].
The empirical RUSLE model implemented in this study was applied to predict potential erosion risk and soil erosion risk for the Epworth district of the Harare Metropolitan Province (Figure S2 and Figure 5). The calculation of the potential erosion risk is based on the assumption that there is no land use and no land management as well as no support practice; potential erosion is understood as the erosion processes only controlled by physical factors. Consequently, potential erosion risk depicts areas vulnerable to erosion even without considering land use [65]. For the Epworth district the potential erosion risk was averaged at 13.2 t ha−1 yr−1 between 1984 and 2018, significantly exceeding the soil loss tolerance limit of 1 t ha−1 yr−1 for the tropics [3,100].
Estimated soil losses due to soil erosion risk averages amounted to 1.31 t ha−1 yr−1 and 1.12 t ha−1 yr−1 in the years 2000 and 2018, respectively. Correspondingly, the revealed soil loss due to soil erosion risk for Epworth district was slightly above the recommended tolerable limits of 1 t ha−1 yr−1 [3,65,100]. Considering the proposed range of tolerable maximum annual soil loss in tropical regions, it can be deduced that the mean estimated soil loss during all study periods slightly causes irreversible soil erosion. The occurrence of high potential erosion risk compared to low soil erosion risk is due to the assumption that for the calculation of potential erosion risk, the factors land cover and management (C) and support practice (P) are not considered, and consequently, these factors are mathematically handled as identity elements in order to assess the impact of RUSLE “natural factors” on the study area. In the application of the RUSLE this corresponds to dealing with C and P factors as bare ground [3]. Such conditions earmark the impact of soil erosion on cleared land area and also reveal the significance of vegetation cover in dissipating raindrop energy impact on the bare ground. Nevertheless, such scenarios merely occur in urban areas as a result of built-up densification and spread of impervious surfaces unless croplands and disturbed green spaces exist at larger scale.
The decrease of soil loss due soil erosion risk for the Epworth district calculated for the period 2000–2018 probably resulted from the expansion of built-up areas at the expense of green spaces and cropland areas and thus with sealing the underground impeding surface exposure to erosion [5,103]. However, soil disturbance has been triggered by human activities in the district due to high population pressure and demand for shelter [104,105]. This contributes to increased soil loss within the urban built-up areas compared to the previous years (Table 8). The replacement of green spaces with impervious surfaces causes a reduction in surface area for water infiltration [106,107], causing increased overland flow either by sheet flow or by concentrated surface runoff. This substantially threatens soil loss through erosion where the overland flow reaches areas where soils are exposed to the surface making them highly vulnerable landscapes with mostly bare grounds [5,7].
Construction activities create artificial slopes and reduce vegetation canopy cover exposing soil surfaces to raindrop impact thereby exacerbating rates of soil detachment and transportation during rainstorm events [30]. Beyond, construction activities repeated earth movements affect the stability of the soil structure and increase soil erodibility [74,108] and soil compaction resulting in the reduction of infiltration capacity and increase surface runoff generation [107,108]. Therefore, continued LULC changes dominated by the expansion of built-up areas largely perpetuate soil loss within the built-up areas [103]. To ensure that soil erosion risk thresholds remain sustainable, land suitability analysis should be considered to enhance land use allocation and proper land servicing require to be implemented by responsible authorities in order to meet competing demands for land. The findings of this study highlight that zones with high potential erosion risk and soil erosion risk correspond to zones of strong relief and high slope lengths (high LS factors).
The weighted C factor values derived based on the biophysical properties observed and measured in the field [42,66]. The approach integrated LULC maps derived from satellite images. However, there are limitations on the use of LULC data resulting from misclassifications, heterogeneity and spatial distribution of vegetation densities across the entire district. The resulting approach anticipates that the same LULC class poses the same land cover and management factor C value [47,89,109,110]. However, uncertainties and limitations influencing the results could have also emanated from omitting other biophysical characteristics during estimations of the C factor including surface roughness and below-ground biomass [42,111,112].
An analysis of the relationship between LULC and estimated soil erosion risk was performed by overlaying LULC and soil erosion risk maps for the time slices 2000 and 2018 (Table 8). This relationship has been observed as a useful tool monitoring patterns of LULC change against soil erosion risk [3,100]; the analysis reveal spatial patterns and changes for each LULC class majorly as influenced by human activities in relation to soil erosion risk. Soil erosion risk was extremely high in rain fed croplands in the year 2000, highlighting their vulnerability to water induced erosion. This might have been propelled by leaving little to no crop residues as soil cover on the fields, exposing bare soil to rainfall at the onset of rainy season [22,24,113]. In addition, pre-season land preparation exposes fine tilled land to raindrop impact, exacerbating soil loss due to reduced surface and ground cover in dissipating and scattering raindrop energy [114], while concurrently ploughing increases surface roughness and pore volume, both fostering infiltration [115,116].
The estimated soil loss in areas covered by green spaces initially was predicted to be high considering the vegetation’s ability to intercept raindrop impacts [66]. However, it has to be clearly differentiated what character the vegetation cover has and whether the area is disturbed or undisturbed by human impact. Especially areas with sparse disturbed green spaces are exposed to erosion processes by the first rains coming after a dry period due to the hydrophobic character of dried out soils while parallel the soil stabilization by roots is insufficient [22,23]; footpaths spreading across areas increase soil sealing along the paths and thus foster generation of concentrated runoff [24,117]. Other local human activities including harvesting of firewood for domestic use, burning of bricks during their processing and fencing of gardens have negative impacts on vegetation cover. The more intense land use gets, especially transferring fallow land covered by sparse green spaces into built-up areas, the higher the soil erosion risk in this area, which is highly affected by the increasing areas characterized by sealed and impervious surfaces [118]. Soil erosion risk investigations on the Epworth district show that, despite the decreasing of the overall estimated soil loss from 50,408 tons to 42,934 tons between 2000 and 2018, the expansion of urban built-up areas at the expense of croplands and green spaces has locally increased soil loss risk.
The expansion and spread of residential areas have been linked with increases in soil erosion risk within areas of intense human activities (Table 8). The relationship between LULC and soil erosion risk analysis reveals that in concentrated residential areas in 2000 about 14,018 tons of estimated soil loss occurred, while about 19,858 tons of estimated soil loss occurred in 2018. The estimated increase of soil loss in concentrated residential areas corresponds to the massive growth of built-up areas and the coinciding increase of impervious surfaces [119]. The observed changes in built-up areas are significantly attributed to population growth in the Epworth district. The population was estimated to have increased from approximately 114,067 in 1992 to about 167,462 in 2012 [57]. The resulting reduced infiltration rates contributing to high rates of surface runoff [120], result in soil damage down the slope and parallel to roads [117]. The increase in estimated soil loss in the residential areas of the Epworth district are comparable with those in Kinshasa/DRC, where it is observed that highest soil erosion risk rates spatially correlate with steep slopes along river flanks and increasing density of informal settlements [38]. Estimated soil loss in less concentrated residential areas increased from 6218 tons in 2000 to approximately 12,203 tons in 2018. This corresponds to a doubling in size of the land use class “less concentrated residential areas” in the Epworth district from 15.5% in 2000 to 31.5% in 2018 (Table 8). The increase in soil loss for less concentrated residential areas may have resulted from the registered decrease in green spaces (Table 8). This is possibly attributed to high population pressure and demand for shelter, hence propelling landowners to informally construct shacks on their backyards to curtail housing demand. These activities inevitably occur at the expense of green spaces and ground cover resulting from clearing of land for construction [6,121].
Steep areas have high potential erosion risk as well as estimated soil erosion risk, especially on the flanks of the streams [122,123]. This corresponds to the high loading of slope steepness in the RUSLE model [124,125]. The combination of slope inclination with slope length contributes to the cumulative effect of increasing surface flow with increasing the drainage basin resulting in increased soil erosion risk [30,38]. For the Epworth district, wall around homesteads and industrial areas most likely acted as physical barriers for surface runoff, reducing slope length and affecting flow velocity and flow direction of surface runoff.
The spatial pattern of current soil erosion damage as documented in a diagram map is based of field survey in December 2019 (Figure 5 and Figure 6). The field measurements data serve as a validation tool for the estimated soil erosion risk modelling for the Epworth district. The plot-based field measurements recorded reflect the extent of erosion prevalent within the mapped area. The RUSLE model is presumed to predict the amount of soil moved on the field [126,127]; henceforth, the spatial damage data from the field measurements represented the proportions of soil moved. These mapped sites of the spatial extent damage map spatially concur with soil loss estimated using the RUSLE model. This, however, improves the model evaluation despite the lack of sheet erosion assessment. Even though, [128,129] reiterated that empirical modelling requires long-term field measurements and the analysis of sedimentation rates for validation purposes. The utilization of point-like plot-based data from field surveys for validation improves the evaluation of the model outputs and its understanding [55,130,131]. Validation of soil erosion risk modelling by comparison with outcomes of current damage mapping was further coupled with an analysis of the relationship between LULC and estimated soil erosion risk [62,100].
A total area of 765 m2 was subjected to soil erosion in 2019 recorded for the 49 randomly selected sampling plots covering an area of 15.7 ha in the Epworth district. The randomized locations of the plot measurements indicate that eroded areas occurred in high frequency in the southwestern and southeastern parts of the Epworth district. Evidently, the spatial extent soil erosion damage in the diagram map indicates altogether low percentages of disturbance with a range of 0%–1.4% damaged area in the surveyed plots. However, during the field survey, the observed damage resulting from soil erosion (Figure 6) appears to be greater in extent and magnitude than the depicted damage illustrated on the soil erosion damage map (Figure 5). In comparison to the maps displaying the soil erosion risk based on the application of the RUSLE model, the spatial extent of soil erosion damage has been observed to spatially concur majorly in the southwestern and southeastern parts of the district where LS factor is presumably high. Areas of high surface damage could be identified in the southeast of the district and were predominantly observed in croplands and areas undergoing construction.
The field survey positively contributes to the study through the identification and registration of areas vulnerable to soil erosion (Figure 5), which therefore lessens the burden of the resource strained land managers and local boards in developing conservation strategies direct on hotspots rather than concentrating on the entire district [5,132]. Areas of high soil erosion risk which occur in zones with high topographic (LS) factors can be confirmed by strong surface damages (Figure 5). Nevertheless, the exclusion of sheet erosion recording and quantification during conducting damage mapping reduces the usability of the soil erosion damage map for the validation of the soil erosion risk mapping applying the RUSLE model [133]; this is because sheet erosion has a major contribution to erosion damage and is included in the RUSLE model [64,86,125]. However, the conclusions were drawn on the basis of field observations due to the heterogeneity of the urban set-up and the widespread nature of impervious surfaces.
The modelling findings reveal that topographic characteristics (LS factor) significantly influence potential erosion and soil erosion risk in the Epworth district, which concurs with the findings by [23] that model-based distribution patterns of soil erosion risk in the area of Windhoek, Namibia were mainly defined basing on the spatial structure of slope. The high soil erosion risk observed on sloppy areas in Epworth district corresponds to areas with convex to straight profile curvature and to the occurrence of ridges [134,135]. Furthermore, along the channel flanks modelled soil erosion risk was in general high, predominantly controlled by relief [122]. Due to the southward drainage of the stream network in Epworth district and thus, southward increasing incision of the streams, increased soil erosion risk along the river flanks can be predominantly observed and also on the southeast areas (Figure 5). In contrast, the analysis of the effects of LULC change on soil erosion pointed out that increasing distribution of built-up areas as a result of high population pressure and demand for shelter substantially propels soil erosion risk within residential areas.
5. Conclusions
Soil erosion is a global environmental concern impacting negatively on agricultural productivity, accessibility to properties and also posing flooding risk in urban areas. The empirical RUSLE model was implemented to assess vulnerable areas and the computation of soil erosion risk through the integration of GIS and remote sensing techniques. The quantitative assessment of average annual soil loss for the Epworth district using the RUSLE model considers climate, soil, land use and topographic datasets as input parameters. Areas with high soil erosion risk were found to spatially correlate with topographic characteristics, especially slope length and slope steepness. The unrestricted LULC changes resulting from rampart informal settlements growth have accentuated soil erosion risk in the Epworth district. The analysis of LULC and estimated soil erosion risk improves the understanding of the spatial distribution patterns of soil loss for the different land uses in the years 2000 and 2018. The predicted soil erosion loss in the Epworth district amounted to 50,408 tons in 2000 while, 42,934 tons were estimated for 2018. Thus, the findings reveal that estimated soil erosion risk in total decreased over the study period (2000–2018). This is attributed to the reduction of croplands and areas covered by green spaces at the expense of built-up areas. Soil loss massively increased in the residential areas from 20,236 tons in 2000 to 32,061 tons 2018, regardless of the concentration of built-up areas (concentrated and less concentrated residential areas); in total the area covered by residential areas almost doubled between 2000 and 2018. Increasing impervious surfaces, sealed areas and avoidance of paved areas during high traffic flow have been observed as contributing factors towards increased generation of surface runoff and hereby affected soil erosion risk in the growing residential areas of the Epworth district.
The soil erosion damage map generated from the field measurements served as a validation tool for the study as it revealed areas vulnerable to soil erosion within the Epworth district that concur with the results of the application of soil erosion risk models. The area affected by soil erosion in the surveyed plots showed damages of 0%–1.4% in the spatial extent of the mapped plot area. Therefore, field mapping data have been observed as necessary in ascertaining and improving an understanding of quantitative soil erosion modelling. Model validation demonstrated that the RUSLE model performance was good due to positive correlation between field measurements and model results basing on sample areas (r = 0.76, R2 = 0.58, p < 0.05). It can be concluded that the spread of urban built-up areas without implementation of sound conservation practices, such as proper land suitability analysis and the construction of runoff drainage canals, will increase soil erosion damage by water. Although, the research has predicted potential erosion risk and soil erosion risk, it is important to outline that there are uncertainties with the modelled data provided by the RUSLE model arising from the lack of site-specific parameterization. Limited studies on water-induced soil erosion for urban areas within the region reduce options for data comparison. However, the computed soil erosion risk maps may assist environmental managers and land and policymakers on planning mitigatory measures for the study area.
Supplementary Materials
The following are available online at
Author Contributions
All authors: conception of the paper, A.K.M.: data acquisition, data analysis, writing of original draft; B.S.: supervision, writing—review and editing, funding acquisition. All authors have read and agreed to the published version of the manuscript.
Funding
This research was funded by the Deutscher Akademischer Austauschdienst—DAAD.
Acknowledgments
The publication of this article was funded by Freie Universität Berlin. We thank our colleagues from Freie Universität Berlin, who provided valuable insights and expertise that greatly assisted this study. We are grateful to the United States Geological Survey (USGS), for providing one Arc-Second Global data from the Shuttle Radar Topography Mission (SRTM) for providing terrain corrected data and the Department of Meteorology in Harare, Zimbabwe for the provision of precipitation data.
Conflicts of Interest
The authors declare no conflict of interest.
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Figures and Tables
Enlarge this image.Figure 1. Location of the Harare Metropolitan Province composed of the Harare urban, Harare rural, Epworth and Chitungwiza districts. (a) Zimbabwe district boundaries depicting the Harare Metropolitan Province; (b) districts of Harare Metropolitan Province; (c) Epworth district showing areas of recorded current soil erosion damages.
Enlarge this image.Figure 2. Flowchart for modelling potential erosion and soil erosion risk in Epworth district. LULC: Land Use and Land Cover, ISRIC: International Soil Reference Information Centre.
Enlarge this image.Figure 3. Input factor maps for modelling potential soil erosion for the Epworth district. (a) shows the overall rainfall erosivity factor between 1984 and 2018; (b) topographic factor (LS) and (c) soil erodibility factor (K).
Enlarge this image.Figure 4. Land Use and Land Cover maps for the Epworth district over the years (a) 2000; (b) 2018 and crop management factor maps for the years (c) 2000; (d) 2018.
Enlarge this image.Figure 5. Estimated soil erosion risk maps for the Epworth district for the years 2000 (a) and 2018 (b); map (c) shows the spatial extent soil erosion damage in percentages based on field mapping in December 2019 for selected plots (plot size measured 40 m × 80 m; in total 49 sites were randomly surveyed, covering varying land use).
Enlarge this image.Figure 6. The extent of soil erosion observed in Epworth district. (a) paved roadside erosion feature; (b) erosion occurring in an unpaved road.
Principal data used for soil erosion modelling showing the resolution and the source.
| Data Type | Resolution | Source |
|---|---|---|
| K Factor | 250 m | Global Soil map and attributes in raster (TIF format) from ISRIC (International Soil Reference Information |
| R Factor | - | Average monthly rainfall data from the Zimbabwe Department of Meteorological Services (Harare) database. |
| LS Factor | 30 m | Digital Elevation Model (DEM) from the United States Geological Survey (USGS) website [72] |
| C Factor | 30 m | Obtained by assigning weighted C factor values to the LULC maps adopted from Marondedze and Schütt, 2019 |
| P Factor | 30 m | The value of 1 was assigned to all-over the study area. |
Location of gauging stations and mean annual precipitation for the study periods.
| Rain Stations | Coordinates | Mean Annual Precipitation (mm) | ||
|---|---|---|---|---|
| 1984–2000 | 2000–2018 | 1984–2018 | ||
| Belvedere | 17°50′ S, 31°01′ E | 880.2 | 851.1 | 865.7 |
| Airport | 17°55′ S, 31°06′ E | 834.3 | 774.2 | 804.3 |
| Kutsaga | 17°55′ S, 31°08′ E | 804.6 | 792.8 | 798.7 |
C factor values and relative proportion of LULC classes for 1984, 2000 and 2018 [9].
| Land Area (%) | ||||
|---|---|---|---|---|
| Land Use Class | Weighted C Factor Value | 1984 | 2000 | 2018 |
| CBD/Industrial areas | 0.017 | 0.1 | 0.4 | 0.5 |
| LMD (Less concentrated) | 0.066 | 2.8 | 15.5 | 31.5 |
| HD (Concentrated) | 0.083 | 16.6 | 38 | 52.6 |
| Irrigated cropland | 0.166 | 1 | 7 | 0.4 |
| Rainfed cropland | 0.239 | 17.9 | 17 | 9.1 |
| Green spaces | 0.03 | 61.5 | 22 | 5.9 |
| Water | 0 | 0.1 | 0.1 | 0.04 |
Note: *CBD: Central Business District; *LMD: Low-Medium Density; *HD: High Density.
Table 4Description of study area Land Use and Land Cover (LULC) classes (data source: [9]).
| ID | Land Use and Land Cover Class | Description |
|---|---|---|
| 1 | CBD/Industries | Industries and central business district defined with high fraction of impervious surfaces mainly buildings and little proportion of vegetation |
| 2 | LMD residential (less concentrated) | Leafy and well established low and medium density suburbs surrounded with high vegetation |
| 3 | HD residential (concentrated) | High density residential areas with low vegetation cover or clustered settlements with areas undergoing developments and bare exposed land |
| 4 | Irrigated cropland | Cultivated land under irrigation schemes |
| 5 | Rainfed cropland | Cultivated land or land with crop residues after harvesting |
| 6 | Green spaces | All wooded areas, shrubs and bushes, riverine vegetation and grass covered areas |
| 7 | Water | Areas occupied by water, rivers, wetlands, reservoirs and dams |
*CBD: Central Business District; *LMD: Low-Medium Density; *HD: High Density.
Table 5Potential erosion risk classes with erosion rate and area covered proportion.
| Soil Erosion Risk | Soil Loss |
1984–2018 | |
|---|---|---|---|
| Area (km2) | Area (%) | ||
| Low | 0–1 | 0.5 | 1.4 |
| Moderate | 1–2 | 2.3 | 6.6 |
| High | 2–5 | 0.4 | 1.1 |
| Very high | 5–10 | 11.3 | 32.3 |
| Extreme | >10 | 20.5 | 58.6 |
Estimated soil erosion risk in Epworth district for 2000 and 2018.
| Soil Loss |
Soil Erosion Risk | 2000 | 2018 | ||
|---|---|---|---|---|---|
| Area (km2) | Area (%) | Area (km2) | Area (%) | ||
| 0–1 | Low | 19.6 | 56.3 | 20.7 | 59.5 |
| 1–2 | Moderate | 9.0 | 25.9 | 10.2 | 29.3 |
| 2–5 | High | 5.2 | 15.0 | 3.5 | 10.0 |
| 5–10 | Very high | 0.9 | 2.5 | 0.4 | 1.1 |
| >10 | Extreme | 0.1 | 0.3 | 0.04 | 0.1 |
The extent of soil erosion in Epworth district summarized for 2019 field survey.
| Spatial Eroded Area (m2) | |
|---|---|
| Number of mapped sites | 49 |
| Total extent of erosion | 765 |
| Mean | 31 |
| Standard Deviation | 10.5 |
| Standard Mean Error | 1.5 |
Estimated soil loss for the different Land Use and Land Cover (LULC) classes in Epworth district based on the assessment of soil erosion risk and LULC analysis for the years 2000 and 2018.
| LULC Class. | 2000 | 2018 | ||||
|---|---|---|---|---|---|---|
| Soil Loss (tons) | Area (km2) | Percentage (%) | Soil Loss (tons) | Area (km2) | Percentage (%) | |
| CBD/Industrial area | 40 | 0.12 | 0.4 | 47 | 0.19 | 0.5 |
| LMD (less concentrated) | 6218 | 5.41 | 15.5 | 12,203 | 10.96 | 31.5 |
| HD (concentrated) | 14,018 | 13.17 | 38 | 19,858 | 18.32 | 52.6 |
| Irrigated cropland | 6970 | 2.45 | 7 | 733 | 0.13 | 0.38 |
| Rainfed cropland | 19,228 | 5.85 | 17 | 9239 | 3.16 | 9.1 |
| Green spaces | 3934 | 7.78 | 22 | 854 | 2.06 | 5.9 |
| Water | 0 | 0.05 | 0.1 | 0 | 0.01 | 0.04 |
| Total | 50,408 | 34.83 | 100 | 42,934 | 34.83 | 100 |
*CBD: Central Business District; *LMD: Low-Medium Density; *HD: High Density.
© 2020 by the authors.