Quantiﬁcation of Soil Losses along the Coastal Protected Areas in Kenya

: Monitoring of improper soil erosion empowered by water is constantly adding more risk to the natural resource mitigation scenarios, especially in developing countries. The demographical pattern and the rate of growth, in addition to the impairments of the rainfall pattern, are consequently disposed to adverse environmental disturbances. The current research goal is to evaluate soil erosion triggered by water in the coastal area of Kenya on the district level, and also in protected areas. The Revised Universal Soil Loss Equation (RUSLE) model was exercised to estimate the soil loss in the designated study area. RUSLE input parameters were functionally realized in terms of rainfall and runo ﬀ erosivity factor (R), soil erodibility factor (K), slope length and gradient factor (LS), land cover management factor (C) and slope factor (P). The realization of RUSLE input parameters was carried out using di ﬀ erent dataset sources, including meteorological data, soil / geology maps, the Digital Elevation Model (DEM) and processing of satellite imagery. Out of 26 districts in coastal area, eight districts were projected to have mean annual soil loss rates of > 10 t · ha − 1 · y − 1 : Kololenli (19.709 t · ha − 1 · y − 1 ), Kubo (14.36 t · ha − 1 · y − 1 ), Matuga (19.32 t · ha − 1 · y − 1 ), Changamwe (26.7 t · ha − 1 · y − 1 ), Kisauni (16.23 t · ha − 1 · y − 1 ), Likoni (27.9 t · ha − 1 · y − 1 ), Mwatate (15.9 t · ha − 1 · y − 1 ) and Wundanyi (26.51 t · ha − 1 · y − 1 ). Out of 34 protected areas at the coastal areas, only four were projected to have high soil loss estimation rates > 10 t · ha − 1 · y − 1 : Taita Hills (11.12 t · ha − 1 · y − 1 ), Gonja (18.52 t · ha − 1 · y − 1 ), Mailuganji (13.75.74 t · ha − 1 · y − 1 ), and Shimba Hills (15.06 t · ha − 1 · y − 1 ). In order to mitigate soil erosion in Kenya’s coastal areas, it is crucial to regulate the anthropogenic disturbances embedded mainly in deforestation of the timberlands, in addition to the natural deforestation process caused by the wildﬁres.


Introduction
Climate change and its corresponding influences on natural resources are one of the greatest challenges that the scientific community has yet to answer The dynamicity and the complexity of climate change enforce the drastic impacts on natural ecosystems, especially the fragile ecosystems [1]. Climate prediction models suggest the alteration of the future climate to be shifted toward higher frequencies and intensities of extreme events on a global scale, with a substantial lack of mitigation plans on the local scale, founded on the Intergovernmental Panel for Climate Change findings In addition, many of the changes in land use and land cover throughout current years have been compelled by the population growth burden, subsequent environmental pressure due to deforestation, the tourism industry, and agriculture [35,36].
It should be noted that the RUSLE model was chosen because the required data is relatively limited, and because of its simplicity, as pointed out by other authors [37][38][39]. The RUSLE model has been studied in depth, and has been widely used on different scales to estimate the loss from soil erosion and to plan the control of erosion imbalances for different land cover categories (e.g., agricultural land, pasture and forests) [23,40]. The RUSLE model and other methods of simulating soil erosion (such as the Coordination of Information on the Environment (CORINE model) can provide the most detailed information about soil erosion risks [27,41,42]. The RUSLE model usually uses auxiliary data provided free of charge in a geographic information system as an alternative method, because it measures soil, and so measuring erosion is expensive and takes time [43,44]. Although RUSLE is considered the main model for evaluating soil erosion, the availability of data to obtain certain RUSLE parameters is the biggest limitation, and it is impossible to maximize accuracy and unify the RUSLE processing method [45]. Model-based methods involve uncertainty in the calculation of each factor [46,47]. This disadvantage is common to all model-based methods [37].
The objectives of the current study are: (a) to assess soil erosion by water in coastal areas of Kenya; (b) to estimate soil erosion for different land cover and land use types, and for protected areas; and (c) to assess soil erosion risks for protected coastal areas of Kenya using the Revised Universal Soil Loss Equation (RUSLE) model [48].

Study Area Description
The coastal areas in Kenya are located between 2 15 • N and 4 0 • S, and encompass around 400 km along the Umba River headed North to the Somalian border, as shown in Figure 1. The climatic condition is defined to be a tropical climate with no all-year-round humidity. Daily temperature ranges between 30 • C and 22 • C. The mean annual precipitation is recorded to be around 1000 mm, falling in two rainy seasons. The long rainy season is from April until July, with average monthly rainfall of 134 mm, and the short rainy season is from October until December, with average monthly rainfall of 100 mm. the average monthly rainfall in the dry season is recoded to be 37 mm [49].  Kenyan protected areas play a crucial role in nature conservation strategies and in maintaining significant biodiversity in the country. The current natural ecosystems assist the Gross National Income (GNI) to stay functional, in terms of serving foraging, woods and clean water, as well as the Mijikenda provenance's spiritual and cultural services [50,51]. Nevertheless, the Kenyan protected areas are under continuous pressure from the surrounding anthropogenic activities, especially when they are exaggerated by the climatic variabilities [52,53].
The total number of protected areas in Kenya is 411, expanding over 72,544.49 Km 2 , 34 of which are listed as coastal protected areas (Figure 2), consistent with the United Nations Environment Program (UNEP) and the World Conservation Monitoring Center (WCMC) [54].

Methodological Framework
The current methodological framework is based on the implementation of the RUSLE empirical model (Table S1). The model shows significant results at different scales of applications [22,55,56]. The estimation of soil erosion based on the RUSLE model practices several natural and human-prompted basics [57].
where A = Annual soil loss (t −1 ·ha −1 ·y −1 ), R = Rainfall-runoff erosivity factor (MJ −1 ·mm −1 ·ha −1 ·h −1 ·y −1 ), K = Soil erodibility factor (t −1 ·ha −1 ·h −1 ·ha −1 ·MJ −1 ·mm −1 ), LS = Slope length and slope steepness factors, C = Cover management factor, P = Conservation support practice factor equal to 1 in our study [58]. Basically, rill and sheet erosion are induced by surface water runoff, therefore the estimation of the rainfall erosivity factor (R) is essential. The seasonal variabilities and intensities of the rainfall (Figure 3) must be taken into consideration, specifically when the raindrop size is questioned. Larger raindrops contain enough kinetic energy to detach soil particles (Merlotto et al., 2016). The average monthly rainfall amounts for the period 1985-2015 were downloaded from http://www.cru.uea.ac.uk/data, considering 10 stations within 22 km of the Kenyan coastline and distributed along the entire coastline [59]. The ratio of soil loss due to the practices of land use and land cover type is called the cover management factor (C) [23,31]. For the coastal area of Kenya under investigation, C factor is estimated based on the Normalized Difference Vegetation Index (NDVI). We used Landsat 8 with temporal resolution of August 24, 2015, and a spatial resolution of 30 m, to extract NDVI information to use in the equation (7). It was acquired from the US geological survey online portal (https://earthexplorer.usgs.gov/). The following equation: Where = 2 and = 1 is the determined parameters of the NDVI. The erosion control practices (P) factor expresses the effects of conservation practices that reduce the amount and rate of water runoff, which reduces erosion due to agricultural management practices such as contour tillage and planting, strip-cropping, terracing and subsurface drainage [60,61]. Unfortunately, due to the lack of financial means, these conservation practices are still undeveloped in developing countries, including the area of interest in this study [62]. Furthermore, establishment of a P factor map at the large watershed scale with complex land use systems is nearly impossible using the Wischmeier and Smith (1978) method, in which P factor is estimable based on the slope gradient and different support practices, such as terracing, contour tillage, etc. [63,64]. Therefore, the P factor value ranges from zero, as an index of good conservation practice, to one, as an index of poor conservation practice or no protective management in our study area [45,65].
The soil susceptibility to erosion is called the soil erodibility factor (K). The present study utilized the K factor [66], estimated using soil properties (sand, clay, silt and organic carbon fractions) compiled by the Africa Soil Information Service (AfSIS) with high spatial resolution of 250 m (http://www.isric.org/data/afsoilgrids250m, Table 1) [67], and Equation 8 proposed by Williams (1995) [68]. For validation reasons the results dataset from Fenta et al. (2020) [69] was kindly provided to compare our own K factor result.
) is a factor that gives low soil erodibility factors for soils with high coarse- According to Wischmeier and Smith (1978), Rainfall erosivity factor can be estimated as follows: 1.735 × 10 1.5 log 10( Pi 2 P ) − 0.08188 (2) where, R factor was calculated and interpolated over the whole coastal area using a geostatistical model, P i is the monthly rainfall in mm, P is the annual rainfall in mm. The influence of the topographical factors in soil erosion estimation cannot be neglected. For that reason, the LS factor is estimated following Panagos et al. (2015a): where S i.j = is the slope factor for the grid cell with coordinates (i.j), D = the grid cell size (m), x i. j = (sin a i.j + cos a i.j ), a i.j = aspect direction for the grid cell with the coordinates (i. j), A i.j−in = Flow accumulation with (i.j) (m 2 ), β = the ratio of inter-rill erosion, θ = the slope.
The ratio of soil loss due to the practices of land use and land cover type is called the cover management factor (C) [23,31]. For the coastal area of Kenya under investigation, C factor is estimated based on the Normalized Difference Vegetation Index (NDVI). We used Landsat 8 with temporal resolution of August 24, 2015, and a spatial resolution of 30 m, to extract NDVI information to use in the equation (7). It was acquired from the US geological survey online portal (https: //earthexplorer.usgs.gov/). The following equation: where a = 2 and β = 1 is the determined parameters of the NDVI. The erosion control practices (P) factor expresses the effects of conservation practices that reduce the amount and rate of water runoff, which reduces erosion due to agricultural management practices such as contour tillage and planting, strip-cropping, terracing and subsurface drainage [60,61]. Unfortunately, due to the lack of financial means, these conservation practices are still undeveloped in developing countries, including the area of interest in this study [62]. Furthermore, establishment of a P factor map at the large watershed scale with complex land use systems is nearly impossible using the Wischmeier and Smith (1978) method, in which P factor is estimable based on the slope gradient and different support practices, such as terracing, contour tillage, etc. [63,64]. Therefore, the P factor value ranges from zero, as an index of good conservation practice, to one, as an index of poor conservation practice or no protective management in our study area [45,65].
The soil susceptibility to erosion is called the soil erodibility factor (K). The present study utilized the K factor [66], estimated using soil properties (sand, clay, silt and organic carbon fractions) compiled by the Africa Soil Information Service (AfSIS) with high spatial resolution of 250 m (http://www.isric.org/data/afsoilgrids250m, Table 1) [67], and Equation 8 proposed by Williams (1995) [68]. For validation reasons the results dataset from Fenta et al. (2020) [69] was kindly provided to compare our own K factor result.
where f csand (Equation (9)) is a factor that gives low soil erodibility factors for soils with high coarse-sand contents and high values for soils with little sand, f cl−si (Equation (10)) is a factor that gives low soil erodibility factors for soils with high clay to silt ratios, f orgc (Equation (11)) is a factor that reduces soil erodibility for soils with high organic carbon content, and f hisand (Equation (12)) is a factor that reduces soil erodibility for soils with extremely high sand contents [70].
where m s is the percent sand content (0.05-2.00 mm diameter particles), m silt is the percent silt content (0.002-0.05 mm diameter particles), m c is the percent clay content (<0.002 mm diameter particles), and orgC is the percent organic carbon of the layer (%).

Results and Discussion
The estimation of soil loss using the RUSLE model necessitates a GIS environment where the multiplication RUSLE parameters take place. Study area slope is shown in Figure 4, and rainfall erosivity factor, soil erodibility factor, cover management factor, slope length and steepness factor, conservation support practice factor, and the land use and land cover map are demonstrated in Figure 5, respectively. Meanwhile, the estimation of soil erosion for the entire coastal area of Kenya is represented in Figure 6. The spatial configuration of the land use and land cover (LULC) types and the corresponding soil erosion estimation is proven in Table 2, except for tree cover areas, shrub cover areas and grassland; all other land cover types had unsustainable soil erosion.       The study area of Fenta et al. (2020) [69] also covers the coastal area of Kenya. Therefore, it could be used for validation of our results. A spatial analysis of the estimated soil loss risk at a district level revealed that 8 of the 26 districts of the coastal area (Table 3) are unprotected to an overall mean erosion risk rate >10 t·ha −1 ·y −1 , and 8 districts, including Kololenli (19.709 t·ha −1 ·y −1 ), Kubo (14.36  The study area of Fenta et al. (2020) [69] also covers the coastal area of Kenya. Therefore, it could be used for validation of our results. A spatial analysis of the estimated soil loss risk at a district level revealed that 8 of the 26 districts of the coastal area (Table 3) are unprotected to an overall mean erosion risk rate >10 t·ha −1 ·y −1 , and 8 districts, including Kololenli (19.709 t·ha −1 ·y −1 ), Kubo (14.36 The study area of Fenta et al. (2020) [69] also covers the coastal area of Kenya. Therefore, it could be used for validation of our results. A spatial analysis of the estimated soil loss risk at a district level revealed that 8 of the 26 districts of the coastal area (Table 3)    The soil erosion acceptance value assists as a foundation for discerning whether the soil poses a prospective threat for drastic soil loss, or mostly for soil poverty [72]. As 8 districts of the coastal area have an average soil erosion of 10 t·ha −1 ·y −1 , this is a prone land as an unsustainable average soil erosion rate is greater than 10 t·ha −1 ·y −1 (Table 3) [73]. The coast's total soil loss, more than 70%, was in areas with a high soil risk, characterized by a steep slope range of 40% up to more than 100% (Figure 4), and high mean rainfall intensity of 1024 mm·y −1 , as well as a monthly rainfall of 103 mm (Figure 3).
The highest average soil loss is found at the LULC types of built-up areas along the coast, with Mombasa County, Malindi, Kilifi and Lamu representing major urbanization infrastructure at the coastal shoreline. This has key effects on soil erosion acceleration, usually by baring the vegetative cover of the land surface, changing the drainage schemes, and compressing the soil particles throughout the construction process, and consequently covering the land with an impermeable layer of asphalts or concreate. This leads to the surface runoff becoming risky [40,74], hence poor planning and a lack of inter-sectoral coordination, and a lack of adequate enforcement items of management, [75] results in soil erosion in the coastal region.
The total soil losses of the 34 coastal protected areas under investigation in the current research are elaborated in Table 4. Taita Hills (11.12 t·ha −1 ·y −1 ), Gonja (18.52 t·ha −1 ·y −1 ), Mailuganji (13.75 t·ha −1 ·y −1 ) and Shimba Hills (15.06 t·ha −1 ·y −1 ) were the most endangered protected areas, where the annual mean soil erosion is higher than 10 t·ha −1 ·y −1 [76]. The higher soil erosion rate recorded in those protected areas could be based on the extensive anthropogenic activities, as well as the expansion of the deforestation process, either naturally or human-induced [42]. For validation reasons (Table S2), in tropical ecosystems, the threshold of soil erosion lenience was projected by Panagos et al. (2015aPanagos et al. ( , 2017 to be 10 t·ha −1 ·y −1 . The tropical threshold is considered to be higher than soil erosion impacts on crop production in Europe and North America, where the threshold is conducted to be 1 t·ha −1 ·y −1 due to the implemented higher standards of environmental protection [77]. Other examples of soil erosion risk include experiments conducted within the Sierra de Manantlán Biosphere Reserve in Mexico, which is recorded to have a wide range of soil erosion risks from 0 to 100 t·ha −1 ·y −1 , which basically depends on the mountainous topography [78]. Meanwhile, the soil erosion risk in the Chemoga watershed in Ethiopia is projected to be a slope resilience risk, with 0 risks at the downstream and over 80 t·ha −1 ·y −1 at the watershed upstream [79].
The driving forces in coastal areas that make the soil erosion magnitude adverse are mostly climatic and topographic factors, assembled in the form of higher rainfall intensities and rugged terrain. The soil's physical and chemical properties have a particular role in soil loss rate variabilities.

Conclusions
Remote Sensing, incorporated with the Geographical Information System, is exercised in the current study as a cost-effective tool to spatiotemporally quantify soil erosion risks on a local scale of a district-level along the shoreline of the Kenyan protected areas. Accordingly, the Revised Universal Soil Loss Equation model parameters were successfully realized and practiced to project potential soil erosion on a watershed scale, and derive robust results to be considered for soil restoration master plans. The slope is a very effective parameter that must be pointed out to limit soil degradation. Therefore, the transformation of gradients into vegetated terraces is the key element to deaccelerate the soil erosion process. This approach drops the strength of surface runoff erosivity and declines the scale influence of soil erosion. The findings of the present study suggest that areas with an annual average soil erosion potentiality of 10 t·ha −1 ·y −1 and higher must be prioritized in soil conservation plans. Above and beyond this is the providing of essential assistance to the decision-makers to implement effective anti-erosive practices to restrict the soil erosion process.
Author Contributions: Y.H. and M.A. were involved in conceptualizing the problem and identifying data while N.Y. and Q.M. mostly conducted data integration, while E.M. assisted in data analysis. All authors designed the research methods and worked on the entire framework of this research, carried out the discussion and analysis as well as the review of the manuscript. All authors have read and agreed to the published version of the manuscript.