Spatial-Temporal Variability of Future Rainfall Erosivity and Its Impact on Soil Loss Risk in Kenya

: Ongoing climate change poses a major threat to the soil resources of many African countries that mainly rely on an agricultural economy. While arid and semi-arid lands (ASALs) take up most of Kenya’s land mass, approximately 64% of its total croplands lie within mountainous areas with high rainfall, hence, areas highly vulnerable to water erosion. Flooding of the Great Lakes and increasing desertiﬁcation of the ASALs are illustrative cases of the implications of recent precipitation dynamics in Kenya. This study applied the Revised Universal Soil Loss Equation (RUSLE) to estimate future soil erosion rates at the national level based on four Coupled Model Intercomparison Project v5 (CMIP5) models under two Representative Concentration Pathway (RCP) scenarios. Results showed the current soil loss rate to be at 4.76 t ha − 1 yr − 1 and projected an increase in average rainfall erosivity under the two scenarios, except for RCP-2.6 (2030s) and (2080s) for the MIROC-5 model. Future projections revealed an incremental change in rainfall erosivity from the baseline climate by a cumulative average of 39.9% and 61.1% for all scenarios by the 2030s and 2080s, respectively, while soil loss is likely to increase concomitantly by 29% and 60%, respectively. The CCCMA_CANESM2 model under the RCP 8.5 (2080s) scenario projected the highest erosion rate of 15 t ha − 1 yr − 1 over Kenya, which is a maximum increase of above 200%, with the Rift Valley region recording an increase of up to 100% from 7.05 to 14.66 t ha − 1 yr − 1 . As a ﬁrst countrywide future soil erosion study, this assessment provides a useful reference for preventing water erosion and improving ecosystem service security.


Introduction
The global soil erosion rate (36 billion metric tons of soil per annum) is projected to increase by as much as 30-66% over the next half century due to climate change, with the highest rise projected over Sub-Saharan Africa (SSA) [1,2]. Accelerating global warming severity is expected to have adverse effects on crop production, adversely affecting food security, especially in the drought prone SSA countries that mainly rely on rain-fed agriculture [3][4][5][6][7][8]. Future climate variability will most likely intensify rainfall runoff patterns, leading to increased soil water erosion, thus negatively impacting agricultural production in this region [9][10][11][12][13]. This endangers more than half of the human population in SSA countries that are predominantly over-dependent on subsistence farming [14,15]. Within Numerous studies have incorporated projected precipitation data from General Circulation Models (GCMs) to predict rainfall erosivity climatic factor within RUSLE [25,41,[50][51][52]. In Africa, Amanambu et al. [53] observed soil loss increases ranging between 12.2-20.6% within the Niger basin, while in Asia, Doulabian et al. [54] reported a maximum erosion increase of >135% over Iran. Closer to Kenya, a recent study by Moges et al. [10] showed an increase of 23% in soil loss vulnerability in Ethiopia by the year 2050. Nevertheless, only a small number of researchers have estimated soil loss risk in SSA countries using future GCMs [55]. The aim of this study has been to (i) assess the current soil erosion rates at the national scale by physiographic regions and topography and (ii) evaluate the spatial and temporal variability of future rainfall erosivity and soil loss based on GCMs of CMIP 5 under RCP 2.6 (optimistic) and 8.5 (pessimistic) scenarios.

Description of the Study Area
Kenya is an equatorial country that geographically extends between latitudes 4°63′ N and 4°68′ S and longitudes 33°9′ E and 41°9′ E (Figure 1). The country has a moderate tropical climate with large regional variations influenced by multiple factors [13]. It has a high undulating topography and heterogeneous landscape (Figures 2 and 3     Approximately 21.6% of the country experienced a browning trend as per the 2015 land degradation-neutrality (LDN) national baseline [30]) due to climate change and increasing human population pressure. This has endangered sustainable agriculture that forms the country's mainstay economy. In recent years, rainfall-related risks have plagued the country, including floods, landslides and droughts. Currently, most of the Great Lakes have drastically bulged (Figure 4) flooding large tracts on adjacent farmlands and settlements. This has been attributed to a 50-year period of climatic phenomena [39], as well as the recent enhanced OND East African short rains linked to the Indian Ocean Dipole (IOD) [37]. Rampant deforestation in the country's water towers and adverse farming practices within highland areas could also have altered surface water runoff patterns into these lakes [58][59][60].
Appl. Sci. 2021, 11, 9903 6 of 34 Approximately 21.6% of the country experienced a browning trend as per the 2015 land degradation-neutrality (LDN) national baseline [30]) due to climate change and increasing human population pressure. This has endangered sustainable agriculture that forms the country's mainstay economy. In recent years, rainfall-related risks have plagued the country, including floods, landslides and droughts. Currently, most of the Great Lakes have drastically bulged (Figure 4) flooding large tracts on adjacent farmlands and settlements. This has been attributed to a 50-year period of climatic phenomena [39], as well as the recent enhanced OND East African short rains linked to the Indian Ocean Dipole (IOD) [37]. Rampant deforestation in the country's water towers and adverse farming practices within highland areas could also have altered surface water runoff patterns into these lakes [58][59][60].

Observed Precipitation
Observed monthly precipitation data for 29 land-based meteorological stations spanning over a period of 30 years (1970-2000) were sourced from the Kenya Meteorological Department (KMD).The stations covered all the Agroecological Zones and Agroclimatic Zones, as well as low to high altitude regions of the investigated area. Figure 5 gives the mean monthly rainfall (seasonal variability and intensities) over the study area.

Observed Precipitation
Observed monthly precipitation data for 29 land-based meteorological stations spanning over a period of 30 years (1970-2000) were sourced from the Kenya Meteorological Department (KMD).The stations covered all the Agroecological Zones and Agroclimatic Zones, as well as low to high altitude regions of the investigated area. Figure 5 gives the mean monthly rainfall (seasonal variability and intensities) over the study area. Approximately 21.6% of the country experienced a browning trend as per the 2015 land degradation-neutrality (LDN) national baseline [30]) due to climate change and increasing human population pressure. This has endangered sustainable agriculture that forms the country's mainstay economy. In recent years, rainfall-related risks have plagued the country, including floods, landslides and droughts. Currently, most of the Great Lakes have drastically bulged ( Figure 4) flooding large tracts on adjacent farmlands and settlements. This has been attributed to a 50-year period of climatic phenomena [39], as well as the recent enhanced OND East African short rains linked to the Indian Ocean Dipole (IOD) [37]. Rampant deforestation in the country's water towers and adverse farming practices within highland areas could also have altered surface water runoff patterns into these lakes [58][59][60].

Observed Precipitation
Observed monthly precipitation data for 29 land-based meteorological stations spanning over a period of 30 years (1970-2000) were sourced from the Kenya Meteorological Department (KMD).The stations covered all the Agroecological Zones and Agroclimatic Zones, as well as low to high altitude regions of the investigated area. Figure 5 gives the mean monthly rainfall (seasonal variability and intensities) over the study area.

Genereal Circulation Models
Studies by Ongoma et al. [12,61] noted that CCCMA_CANESM2, CESM1-CAM5, CSIRO_MK3.6.0 and MIROC-5 CMPI5 GCMs (Table 1) are among the top eight models that well simulate rainfall patterns for the EA region. The projected precipitation data of these four GCMs were downloaded from the CMPI5 archive of the Climate Change Agriculture and Food Security (CCAFS) database (http://www.ccafs-climate.org/ accessed on 3 January 2021) for the 020-2049 and 2060-2089 periods under two Representative Concentration Pathways: RCP 2.6 (low scenario) and RCP 8.5 (high baseline emission scenario). The RCP 2.6 and 8.5 represent radiactive forcing of 3.0 and 8.5 W/m 2 , respectively, by 2100. The baseline climatic data, at 1 km 2 spatial resouliton, was sourced from the WorldClim database [62]. The WorldClima data form an intergral input in the formulation of future GCMs; thus, they are seen as a good source of baseline climatic data [52]. In this study, future climatic data were statistically bias corrected and downscaled to unify their horizaontal resolution with that of the baseline data using the delta method [63]. Projected rainfall erosivity for the 2030s and 2080s periods were then quantified based on the R factor methods from the RUSLE model ( Figure 6).

Genereal Circulation Models
Studies by Ongoma et al. [12,61] noted that CCCMA_CANESM2, CESM1-CAM5, CSIRO_MK3.6.0 and MIROC-5 CMPI5 GCMs (Table 1) are among the top eight models that well simulate rainfall patterns for the EA region. The projected precipitation data of these four GCMs were downloaded from the CMPI5 archive of the Climate Change Agriculture and Food Security (CCAFS) database (http://www.ccafs-climate.org/ accessed on 3 January 2021) for the 020-2049 and 2060-2089 periods under two Representative Concentration Pathways: RCP 2.6 (low scenario) and RCP 8.5 (high baseline emission scenario). The RCP 2.6 and 8.5 represent radiactive forcing of 3.0 and 8.5 W/m 2 , respectively, by 2100. The baseline climatic data, at 1 km 2 spatial resouliton, was sourced from the WorldClim database [62]. The WorldClima data form an intergral input in the formulation of future GCMs; thus, they are seen as a good source of baseline climatic data [52]. In this study, future climatic data were statistically bias corrected and downscaled to unify their horizaontal resolution with that of the baseline data using the delta method [63]. Projected rainfall erosivity for the 2030s and 2080s periods were then quantified based on the R factor methods from the RUSLE model ( Figure 6).  = * * * * (1) Figure 6. Integration of the climate change scenarios with the R factor used in this study.

RUSLE Model Application
Recent advancements in GIS and remote sensing technologies have made it relatively easy to estimate and monitor soil loss, both at local and global scales (with varying climates), using the empirical RUSLE soil erosion model. As an updated version of the Universal Soil Loss Equation (USLE), the RUSLE equation forecasts soil erosion risk using five environmental and anthropogenic variables: rainfall erosivity, topography, vegetation cover and land support practices. The average annual soil loss (A) in the RUSLE model is expressed by Equation (2) [64] where: A is expressed in MJ mm ha −1 h −1 yr −1 , R is the rainfall erosivity factor (MJ mm ha −1 h −1 yr −1 ), K is the soil erodibility factor (t ha h ha −1 MJ −1 mm −1 ), LS is the slope length and steepness factor (dimensionless), C is the cover management factor (dimensionless, ranges from zero to one) and P is the conservation practice factor (dimensionless) The various geospatial datasets, with their respective sources used in this study, are as summarized in Table 2. The RUSLE model was run to estimate the current soil erosion rate based on baseline rainfall erosivity, while future soil losses were forecasted based on respective projected rainfall erosivity scenarios.

Rainfall Erosivity (R) Factor
Rainfall erosivity is the most significant factor in the RUSLE model ( [23,67]); it influences approximately 80% of the total soil loss. The rainfall erosivity factor is the product of raindrop kinetic energy (E) and the maximum 30 min rainfall intensity (I 30 ) [64]. Acquiring complete pluviographic data (with a minimum of 20 years in order to capture cyclical rainfall dynamics) is a huge setback in African countries due to insufficient gauged meteorological stations [49]. Different approaches have been used to derive the R factor within Kenya. Maeda et al. [9] and Schürz et al. [44] utilized the Fourier index (Equations (2)-(4)) and Modified Fourier Index (Equations (5) and (6)), respectively, while Yves [49] applied the model suggested by reference [64] (Equation (7)) to compute rainfall erosivity for the entire coast of Kenya. With access to some rainfall intensity data, Angima [48] and Akali [68] used (Equations (8)-(10)) to determine estimates for the western and central regions, respectively. Various regression realizations have also been exploited, including Renard and Freimund (Equations (11) and (12)) [44], Moore (Equations (13) and (14) where P i is the average monthly rainfall (mm) for month i, P is the mean annual precipitation (mm) and r i is the average monthly erosivity (MJ mm ha −1 h −1 month −1 ).
where α and β are regression coefficients (taken as 50.7 and −1405, respectively, for the African continent [69]), Pi is the average monthly rainfall (mm) for month i and P is the mean annual precipitation (mm).
where Pi is the average monthly rainfall (mm) for month i, P is the mean annual precipitation (mm) and R is the rainfall erosivity (MJ mm ha −1 h −1 y −1 ).
where R is the mean annual rainfall erosivity (MJ mm ha −1 h −1 y −1 ), n is the number of years of data, m j is the number of erosive events in the j year, EI 30 is the rainfall erosivity index of a storm k, e r is the unit rainfall energy (MJ h −1 ), v r is the rainfall depth (mm) during a time period r, I 30 is the maximum rainfall intensity during a 30 min period of the rainfall event (mm h −1 ) and i r is the rainfall intensity during the period (mm h −1 ).
where P ≥ 850 mm (12) where P is the average mean annual rainfall (mm) and R is the rainfall erosivity where KE is the kinetic energy, R is the rainfall erosivity (MJ mm ha −1 h −1 y −1 ) and value 17.02 is a conversion factor from imperial to International (SI) units.
where P is the average mean annual rainfall (mm) and R is the rainfall erosivity (MJ mm ha −1 h −1 y −1 ). In this study, the rainfall erosivity factor was derived using (Equations (11) and (12)), since it has been widely applied in previous similar studies.

Erosivity Density Ratio
The erosivity density (ED) value is the ratio of rainfall runoff erosivity to precipitation and is given by the expression (Equation (16)): where ED (MJ ha −1 h −1 ) [52] is the erosivity density, P annual rainfall (mm) and R is the rainfall erosivity.

Rainfall Erosivity (R) Parameter Evaluation
Monthly precipitation data observed from 29 meteorological synoptic stations were used to evaluate the validity of R factor values derived from the baseline data. Three statistical techniques were used to assess the performance of the R factor of the observation data with that of the baseline data: the coefficient of determination (R 2 ), root mean squared error (RMSE) and Nash-Sutcliff Efficiency (NSE) [70], as shown in Equations (17), (18) and (19), respectively; where Y model i is the baseline rainfall erosivity, Y obs i is the observed rainfall erosivity and Y mean i is the mean of observed and baseline rainfall erosivity.

Soil Erodibility Factor
The soil erodibility (K) factor takes into account the inherent soil properties (including soil texture, organic matter content and the permeability) to quantify soil erodibility or susceptibility to disintegration due to surface water runoff action. This present study utilized the EPIC (erosion-productivity impact calculator) model, as contended by Williams [71], to derive the K factor using the sand, organic, silt and sand soil fractions of the area, as compiled by the Africa Soil Information Service (AfSIS) [66].
where SAN, SIL and CLA are the percentage of sand, silt and clay content, respectively; C is the organic carbon content and SN1 is the sand content subtracted from 1 and divided by 100. F csand (Equation (4)) gives a low soil erodibility factor for soil with coarse sand and a high value for soil with little sand content. F si-cl (Equation (5)) gives a low soil erodibility factor with a high clay to silt ratio; F orgc (Equation (6)) is the factor that reduces soil erodibility for soil with high organic contents. F hisand (Equation (7)) is the factor that reduces soil erodibility for soil with extremely high sand content. Multiplication by the constant value 0.1317 converts the values to International (SI) units.

Slope Length and Slope Steepness (LS) Factor
The dimensionless LS geophysical factor is a product of two terrain derivatives (slope length (L) and slope steepness factor (S)) that expresses the influence of the landscape on soil loss. Kim [72] defined slope length as the distance between the point of origin of overland flow to the point where the surface runoff waters flow into a well defined channel. High slope lengths along steep slopes result in increased surface runoff and, thus, high erosion rates. A 30 m void-filled SRTM DEM was first pre-processed to correct for sink errors before estimating slope length and steepness attributes using the Spatial Analyst extension of ArcMap 10.2 (Environment Systems Research Institute (ESRI) Inc., Redlands, CA, USA). Equation (25) was used to compute the slope length factor, while the S factor was generated using the McCool et al. [73] (1987) method (Equation (28)).
where L i.j = slope length factor for the grid cell with coordinates (i.j); D = the grid cell size (m); X i.j = sina i + cosa i.j ; a i.j = aspect direction for the grid cell with coordinates (i.j); A i.j-in = flow accumulation or contributing area at the inlet of a grid cell with coordinates (i.j) (m 2 ), β = the ratio of inter-rill erosion and θ = the slope in degrees 2.4.7. Cover Management Factor (C) Factor The C factor varies from 0 to 1 and models the conservative property of ground vegetation cover against surface water runoff. This implies that holding other parameters constant, thick canopies, e.g., forests, are more likely to prevent soil erosion better than bare or sparse cover; thus, they are assigned a low C factor value [74]. Phinzi et al. [75] evaluated the use of the Durigon algorithm [76] in deriving C factor values based on NDVI for tropical areas, while Watene [47] recommended the use of LULC maps to assign C parameters in Kenya. Similar to other regional studies in other parts of the world focusing on future rainfall erosivity [3,24,50,77], the nationwide 2016 LULC map was used to obtain C coefficients to compute the present soil erosion status, as well as for the 2030s and 2080s periods. Table 3 shows C factor values used in this research, as sourced from past literature within the East Africa region.

Support Practice (P) Factor
The P factor reflects the effects of man-made soil conservation measures with reference to how they impact the intensity and flow of surface runoff. Renard [23] expressed it as a ratio of the amount of soil erosion in a given area with a unique conservation measure to the corresponding erosion following upward and downward tillage. Various approaches have been employed to determine the P factor, e.g., applying DEM, high resolution imagery and using LULC maps in combination with field inspections. However, only a little has been done to determine support practices across the EA region [80]. A P factor value of one (1), indicating no major support practice in place, was thus adopted for the entire study area [14,49].

Estimated Baseline Soil Erosion Rates in Kenya
The baseline and observed precipitation were statistically compared to evaluate the quality of the simulation produced by the WorldClim data ( Figure 7). A correlation coefficient of approximately 0.91 was observed, indicating the strong relationship between these two variables. The mean annual R factor generated by observed rainfall data varied from 267-31,667 MJ mm ha −1 h −1 yr −1 , with an average of 1714.31 MJ mm ha −1 h −1 yr −1 . In comparison, the WorldClim data presented an R factor (baseline) with an average of 1665.9 MJ mm ha −1 h −1 yr −1 (Figure 8a). The observed and baseline R factors showed 0.93, 273 MJ mm ha −1 h −1 yr −1 and 0.74 for R 2 , RMSE and NSE, respectively, demonstrating a satisfactory model performance. Highland regions with high mean annual rainfall presented R factor values >5000 MJ mm ha −1 h −1 yr −1 .
The Soil erodibility Factor (K) had a mean value of 0.019 t ha h ha −1 MJ −1 mm −1 and ranged from 0.014 to 0.028 t ha h ha −1 MJ −1 mm −1 (Figure 8b), while the cover management Factor (C) had a mean value of 0.089 (Figure 8d). The northern parts of Kenya with bare and arid surfaces recorded the highest C values. Almost 78% of the country had an LS factor between 0 and 10, while 4% had values of more than 10 (even though steep gradients may not necessarily show high erosion susceptibility [50]). The mean baseline soil erosion rate was observed to be 4.76 t ha −1 yr −1 which had a spatial distribution that varied equally with that of the annual precipitation of the country (Figure 9). The baseline erosivity density ED had an average value of 2.77 MJ ha −1 h −1 and a standard deviation 1.06 MJ ha −1 h −1 (Figure 9b). The P factor reflects the effects of man-made soil conservation measures with reference to how they impact the intensity and flow of surface runoff. Renard [23] expressed it as a ratio of the amount of soil erosion in a given area with a unique conservation measure to the corresponding erosion following upward and downward tillage. Various approaches have been employed to determine the P factor, e.g., applying DEM, high resolution imagery and using LULC maps in combination with field inspections. However, only a little has been done to determine support practices across the EA region [80]. A P factor value of one (1), indicating no major support practice in place, was thus adopted for the entire study area [14,49].

Estimated Baseline Soil Erosion Rates in Kenya
The baseline and observed precipitation were statistically compared to evaluate the quality of the simulation produced by the WorldClim data ( Figure 7). A correlation coefficient of approximately 0.91 was observed, indicating the strong relationship between these two variables. The mean annual R factor generated by observed rainfall data varied from 267-31,667 MJ mm ha −1 h −1 yr −1 , with an average of 1714.31 MJ mm ha −1 h −1 yr −1 . In comparison, the WorldClim data presented an R factor (baseline) with an average of 1665.9 MJ mm ha −1 h −1 yr −1 (Figure 8a). The observed and baseline R factors showed 0.93, 273 MJ mm ha −1 h −1 yr −1 and 0.74 for R 2 , RMSE and NSE, respectively, demonstrating a satisfactory model performance. Highland regions with high mean annual rainfall presented R factor values >5000 MJ mm ha −1 h −1 yr −1 .
The Soil erodibility Factor (K) had a mean value of 0.019 t ha h ha −1 MJ −1 mm −1 and ranged from 0.014 to 0.028 t ha h ha −1 MJ −1 mm −1 (Figure 8b), while the cover management Factor (C) had a mean value of 0.089 (Figure 8d). The northern parts of Kenya with bare and arid surfaces recorded the highest C values. Almost 78% of the country had an LS factor between 0 and 10, while 4% had values of more than 10 (even though steep gradients may not necessarily show high erosion susceptibility [50]). The mean baseline soil erosion rate was observed to be 4.76 t ha −1 yr −1 which had a spatial distribution that varied equally with that of the annual precipitation of the country (Figure 9). The baseline erosivity density ED had an average value of 2.77 MJ ha −1 h −1 and a standard deviation 1.06 MJ ha −1 h −1 (Figure 9b).

Classification of Baseline Soil Erosion by Severity
Similar to the study by reference [81], the mean annual soil erosion rate was grouped into six levels of severity: slight, moderate, high, very high, severe and very severe (Table  4) Most of the investigated area had a slight severity erosion rate, while approximately 3.7% of the total area fell within the severe and very severe categories.

Baseline Soil Erosion by Land Use and Land Cover
From the study, croplands which mostly occupy areas with a high mean annual precipitation (central highlands and Lake Victoria parts of the country), had the highest mean erosion rates of 17.86 t ha −1 yr −1 (Table 5). This was followed by Barelands and sparsely vegetated lands in the northern ASALs, which had mean erosion rates of 3.39 t ha −1 yr −1 and 2.62 t ha −1 yr −1 , respectively.

Classification of Baseline Soil Erosion by Severity
Similar to the study by reference [81], the mean annual soil erosion rate was grouped into six levels of severity: slight, moderate, high, very high, severe and very severe (Table 4) Most of the investigated area had a slight severity erosion rate, while approximately 3.7% of the total area fell within the severe and very severe categories.

Baseline Soil Erosion by Land Use and Land Cover
From the study, croplands which mostly occupy areas with a high mean annual precipitation (central highlands and Lake Victoria parts of the country), had the highest mean erosion rates of 17.86 t ha −1 yr −1 (Table 5). This was followed by Barelands and sparsely vegetated lands in the northern ASALs, which had mean erosion rates of 3.39 t ha −1 yr −1 and 2.62 t ha −1 yr −1 , respectively.

Baseline Soil Erosion by ACZ and Climate Zones
Humid areas with an average annual rainfall erosivity of 9100.7 MJ mm ha −1 h −1 yr −1 had the highest mean erosion rate, followed by sub-humid and semi-humid areas, as represented in Table 6. Regions under the tropical rainforest climate around Lake Victoria recorded the highest soil erosion rate (27.84 t ha −1 yr −1 ), while the arid desert climate zones had the lowest mean annual soil loss rate (0.38 t ha −1 yr −1 ) ( Table 7).

Baseline Soil Erosion by Kenya Regions and Basins
The Coastal and Rift Valley regions had mean annual erosion rates below 10 t ha −1 yr −1 as demonstrated by the recent regional works by Yves et al. [49] and Watene et al. [47] ( Table 8). The Lake Victoria region, which is characterized by a high mean annual precipitation (tropical rainforest climate zone) and mostly under croplands, recorded the highest soil loss rate of 25.25 t ha −1 yr −1 ( Table 9). Out of the five basins that form the drainage landscape of Kenya, the Lake Victoria basin was shown to have the highest soil erosion rate, followed by the Great Rift Valley and Tana River basins, respectively.

Baseline Soil Erosion by Slope and Elevation
Lands with slope (>26.8%) had the highest mean erosion rates of about 21.41 t ha −1 yr −1 , while flatter regions (0-7%) recorded the lowest rates of 1.88 t ha −1 yr −1 (Table 10). Highland areas, which are the wettest among all the five elevation classes of the country, also represented the highest mean annual erosion rates (Table 11).

Baseline Soil Erosion by Soil and Landform Types
Soils in areas with a mean annual rainfall erosivity greater than 4500 presented high mean annual soil loss rates above 10 t ha −1 yr −1 (Table 12). Nitosols had the highest erosion rate of 22.69 t ha −1 yr −1 , followed by Andosols and Planosols. Yermosols, which are predominant in the north-eastern parts of the country, recorded the lowest mean soil loss rates. The steep escarpment landforms of the Great Rift Valley presented the highest annual erosion rates, followed by mountainous landforms (Table 13).  Tables 14 and 15 give the statistical details and comparisons of the baseline rainfall erosivity, with projected future values computed from the CMIP5 models under the two RCPs. For the 2080s period, all the GCMs, except MIROC-5 RCP 2.6, revealed an increase in precipitation, with CCCMA_CANESM2 RCP 8.5 recording the highest rise of 85.2% for the entire country. Similarly, only MIROC-5 RCP 2.6, MIROC-5 RCP 8.5 and CSIRO_MK3.6.0 RCP 2.6 GCMs resulted in a precipitation decrease in all of the 2030s projections. Nationwide, there was a general upward trend in rainfall erosivity in all the GCMs in both periods, except for MIROC-5 RCP 2.6 (2030s) and MIROC-5 RCP 8.5 (2080s), which showed −3.9% and −8.2% decreases, respectively (Figures 10 and 11). CCCMA_CANESM2 RCP 2.6 (2080s) produced the largest percentage change of erosivity density, followed by CCCMA_CANESM2 RCP 8.5 (2030s) (Figure 12). The projected mean annual erosion rates for the CCCMA_CANESM2 GCMs were notably high amongst all the other GCMs under all scenarios, with CCCMA_CANESM2 RCP 2.6 recording an increment of 99.16% in the 2030s and a 217.23% increase in the 2080s under RCP 8.5, respectively ( Figure 13). The overall average soil erosion rate of all the GCMs ensembles in the 2030s time-slice demonstrated a 29% increase compared to the baseline period, while for the 2080s period, a +59.87% change was noted ( Figure 14). All the GCMs revealed increments in future erosivity values across all Kenyan regions in the years 2030s and 2080s under RCP 8.5, except CSIRO_MK3.6.0 RCP 8.5 (2030s) (−3.9% for the coastal region), MIROC-5 RCP 8.5 (2030s) (-3.2% and -15.7% for the northern and Rift Valley regions, respectively) and MIROC-5 RCP 8.5 (2080s) (-5% for rift valley region) ( Table 15). The MIROC-5 GCM produced the most varied changes, including the lowest change of −25.1% for the Rift Valley Region in the 2080s under RCP 2.6. The North Eastern Region (NER) had the lowest rainfall erosivity among all the Kenya regions, with a mean rainfall erosivity of 644.8 MJ mm ha −1 h −1 yr −1 in the baseline period. The CSIRO_MK3.6.0 RCP-2.6 (2030s) forecasted the lowest change of −16.3% in this region, followed by MIROC-5, which predicted a decrease ranging from −3.2% to −11.8% during the two time periods for the two emission scenarios. In contrast, the CESM1-CAM5 RCP-8.5 (2080s) predicted the highest increase in rainfall erosivity of approximately 97.3%. The baseline average rainfall erosivity of the Lake Victoria region that falls within the tropical rainforest climate had the highest rise (98.0%) under the CCCMA_CANESM2 RCP-8.5 (2080s) projection. Most of the GCMs produced positive changes in the mean annual erosion rates, with the rift valley region indicating the highest change (107.94% under the CCCMA_CANESM2 RCP 8.5 (2080s)) ( Table 16).  (Table 16).

Discussion
The current study emphasized the merits of the well-tested RUSLE model for monitoring soil erosion rates, which include fast analyses of the impacts of future climate variability on soil erosion (climatic factor), easy adaptability in many GIS platforms (which allows simple computation, even for nationwide erosion risk assessments) and applicability of easily accessible data (both at global and local scales) [82]. The baseline results of the present study were in agreement with other similar works performed within Kenya and around the EA region. The baseline mean annual soil erosion rate was within the range of erosion rates for the EA region, as assessed by references [14,44,83], as well as for the African continent [84]. The baseline soil erosion rate results concurred with the results reported previously in the multiple regional studies [14,16,44,46,49]. The current soil loss rate estimation also corroborated the range of mean erosion rate values for Kenya (5.67 t ha −1 yr −1 ) reported in the Global Soil Erosion Modelling (GloSEM) [1] However, our baseline mean rainfall erosivity was lower than the value derived from the Global Rainfall Erosivity Dataset (GloREDa) [85] for Kenya (2588.1 MJ mm ha −1 h −1 yr −1 ). Schürz et al. [44]    The exception here was the MIROC-5 model, which indicated negative soil erosion rate changes in all the projections under the RCP-2.6 scenario. It was notable that only the MIROC-5 model for the Rift Valley Region gave a reduction in soil erosion rates in all of the future RCP-8.5 (2080s) scenarios (Table 16).

Discussion
The current study emphasized the merits of the well-tested RUSLE model for monitoring soil erosion rates, which include fast analyses of the impacts of future climate variability on soil erosion (climatic factor), easy adaptability in many GIS platforms (which allows simple computation, even for nationwide erosion risk assessments) and applicability of easily accessible data (both at global and local scales) [82]. The baseline results of the present study were in agreement with other similar works performed within Kenya and around the EA region. The baseline mean annual soil erosion rate was within the range of erosion rates for the EA region, as assessed by references [14,44,83], as well as for the African continent [84]. The baseline soil erosion rate results concurred with the results reported previously in the multiple regional studies [14,16,44,46,49]. The current soil loss rate estimation also corroborated the range of mean erosion rate values for Kenya (5.67 t ha −1 yr −1 ) reported in the Global Soil Erosion Modelling (GloSEM) [1] However, our baseline mean rainfall erosivity was lower than the value derived from the Global Rainfall Erosivity Dataset (Glo-REDa) [85] for Kenya (2588.1 MJ mm ha −1 h −1 yr −1 ). Schürz et al. [44] also noted such a discrepancy in the RUSLE model ensemble for Kenya and Uganda. Such visible differences were also noted after making comparisons between the high-resolution rainfall erosivity (GloREDa) and other similar studies. For instance, reference [41] computed a mean rainfall erosivity (8470 MJ mm ha −1 h −1 yr −1 ) that was higher than the GloREDa value for Brazil (7560.8 MJ mm ha −1 h −1 yr −1 ). Similarly, Duulatov et al. [52] observed an erosivity of 1447.7 MJ mm ha −1 h −1 yr −1 against a GloREDa value of 772.7 MJ mm ha −1 h −1 yr −1 over Tajikistan. Table 17 presents a comparison of rainfall erosivity results of this study with past literature R factor values based on multi-years EI 30 calculations. The baseline results of the present study were also compared to previous soil erosion plot studies in Kenya for validation purposes (Table 18). The entire croplands in the country were estimated to have the highest baseline average annual soil loss rate (17.86 t ha −1 yr −1 ), which is comparable to the results reported by Fenta et al. [14]. Most of the croplands are found within areas of steep terrain and highly erodible soils. Okoba and Sterk [86] further attributed these high rates to poor farming practices, where farmers leave the soil bare for extended periods, thus enhancing soil erosion, especially during extreme weather conditions. The total baseline soil erosion rate estimated in Kenya (4.76 t ha −1 yr −1 ) was slightly less than in most regions of the globe: 39.2 t ha −1 yr −1 for Rwanda [87], 0.7-17.9 t ha −1 yr −1 for Europe, 10.8-146 t ha −1 yr −1 for Africa [84] and 0-273 t ha −1 yr −1 for Nepal [81]. In addition, Kenya's erosion rate estimates were below the tolerable limit of 25 t ha −1 yr −1 for environments with young mountains topography, as suggested by Stocking [84]. The spatial pattern of the mean annual rainfall figure had a high degree of similarity with other rainfall distribution maps reported at both the continental and regional scales [88]. Nicholson (2017) [89] largely attributed the aridity of the northeastern parts of Kenya to the divergence of the Somali Jet that drives the Indian Monsoon. Both global and regional phenomena contribute to the heterogeneous nature of Kenya's rainfall pattern. These include the seasonal migration of the ITCZ that accounts for the EA's seasonal rainfall variability [97], solar flux fluctuations and global sea surface temperatures (SST) in the Atlantic Ocean, as well as the Indian and Pacific Ocean anomalies that are linked to the Indian Ocean Dipole (IOD) and the El Niño Southern Oscillation (ENSO), respectively [98]. Other factors, including proximity to large inland lakes and orographic effects, greatly influence moisture exchanges, thus contributing to variable local rainfall patterns [89]. Similar to our results, reference [99] also reported rainfall erosivity values greater than 5000 MJ mm ha −1 h −1 yr −1 in the highland areas of Kenya as well as over the Lake Victoria region [18]. Furthermore, the rainfall erosivity distribution map showed a consistent pattern when compared to previous erosivity maps reported for Kenya.
As presented in other previous regional research [50], this study predicted the soil erosion in Kenya with respect to changes in the climatic factor (rainfall erosivity) within the RUSLE model, while holding the other factors constant. Our results indicate that the overall future soil loss rate is projected to increase significantly for both time slices, due to the increased erosive power of precipitation. There was an overall dramatic positive change in rainfall erosivity in the investigated area across all the models, except for model MIROC-5 under the RCP 2.6 scenario for the 2030s and MIROC-5 under RCP 8.5 for the 2080s. The projected future erosion rates were within the range of estimates given other works in EA [10,24,55], as with the continental predictions [1] (Table 19). For instance, references [25,100,101] observed increases of the mean soil loss rate ranging between +17.8%-28.3% in China, Europe and Uzbekistan (Table 19). Similar to observations made by Panagos et al. [100] across Europe, Kenya having a lower erosion rate than Ethiopia (16.9 t ha −1 yr −1 ) [14] is expected to record larger increments due climate change. High amounts of rainfall can increase the erosive power of runoff and significantly magnify soil loss rates, especially in steep terrains. It is also possible that increases in precipitation can promote vegetation growth and, in turn, enhance plant canopy protection against soil erosion [99]. The highest change in rainfall erosivity was noted in the rift valley, Lake Victoria and central and eastern regions, which consequently resulted in high soil erosion rates. This may lead to other hazardous off-site effects, including sedimentation in water bodies, landslides and flooding in these regions. For most cases, the highest variations in soil loss occurred in the central highlands and western regions of Kenya, which represent areas of intensive agricultural production including the Kenya's grain basket region [43]. Agricultural output in these regions will most likely be negatively affected by climate variability, due to possible increases in soil erosion rates and soil degradation. It is important for the country to adapt to these expected changes and devise soil and water conservation policies that will enhance agricultural resilience (ensure a balance between food security and acceptable soil loss rates). Effective soil conservation methods, as eco-sustainable agricultural technologies, should be continuously implemented in such areas. On the other hand, public policies should encourage more agricultural expansions in regions projected to have reductions in average rainfall erosivity [41] (and consequent low soil loss rates), e.g., the North Eastern region, Coastal region and northern parts of the Rift Valley. The recently launched Galana Kalalu irrigation scheme within the NER and CR can be seen as a profitable shift from the country's overdependence on rain-fed agriculture aimed at improving food security [6,8]. Nevertheless, caution should be taken when converting such natural landscapes into croplands via irrigation, as increased farming in arid lands with no accompanying conservation practices can greatly enhance soil loss rates.
The inability to reflect major rainfall erosivity variability arising from extreme rainfall phenomena is a major setback of the RUSLE empirical method applied. The monthly precipitation data used in the study would not be sufficient to account for the size of raindrops and rainfall duration, which bear a direct correlation to surface erosion. In addition, uncertainties inherent to GCMs, including bias correction methods (Delta method), spatial resolution (resampling to 1 km resolution) and other model assumptions, also pose a major limitation to the present study [53]. Different researchers across the globe have applied multiple climate models with variable future scenarios with an aim to mitigate such uncertainties. Likewise, four GCMs that address differences or inconsistencies arising from climate change signals, e.g., MIROC-5 and CCCMA_CANESM2, were employed based on previous climatology research in Kenya and other similar tropical regions. Different studies have utilized multiple GCMs in soil erosion predictions in order to reduce the uncertainties inherent to climate models. In this study, we selected four models among the top eight models (out of 22 CMPI5 models), which Ongoma et al. [61] demonstrated to well produce rainfall estimations over East Africa. The present research only focused on the impact of the climatic factor (K) and overlooked both anthropogenic (P) and vegetation factors (C), which have been shown to have great significance on soil erosion rates. Holding such pertinent factors at a constant certainly raises questions about the reliability of the future RUSLE predictions [74]. Colman et al. [31] utilized LULC projections based on the contribution of various land uses to greenhouse gas emissions to study future soil erosion rates within the tropical Brazilian Pantanal region and revealed an increasing trend ranging between 40-100%. It is thus important to integrate projected vegetation cover coefficients (based on projections of LULC) in order to improve on the integrity of the results or simulations. The ESRI very recently released a high-resolution Global Land Cover map (10 m) using state-of-the-art artificial intelligence (https://livingatlas.arcgis.com/landcover/ accessed on 3 August 2021), as part of its contribution to climate change action, and it plans to continue this initiative on a yearly basis. This database will be of great importance to environmental modeling studies, e.g., climate change and food security, that require largescale vegetation cover information (especially for developing countries). Xiong et al. [80] noted that there is limited knowledge of EA's current support practices; thus P estimations are often neglected in large-scale research or set to a value of "1", which represents no soil conservation measures as being the worst case scenario. Moreover, Taye et al. [102] reported that the impact of soil conservation measures rapidly decline over time, especially in EA's ASAL regions.

Conclusions
Soil conservation is an integral environmental goal in Kenya as it embarks on its ambitious Kenya Vision 2030 of attaining food security and increasing agricultural production. Kenya also strives to align itself with the United Nations Sustainable Development Goals (SDGs) by enhancing soil protection in order to attain land degradation neutrality by 2030. Recent IPCC reports projected increases in the frequency and/or the intensity of heavy precipitation and pluvial flooding over EA. This increases the risk of soil loss by water in the next century. The present study examined the magnitude of current soil erosion rates (baseline climate) for Kenya, as well as the potential impacts of climate change on rainfall erosivity using multiple GCMs (CCCMA_CANESM2, CESM1-CAM5, CSIRO_MK3.6.0 & MIROC-5) under two greenhouse emissions (RCP2.6 and RCP8.5) for the 2030s and 2080s periods. The baseline soil loss rate was estimated at 4.76 t ha −1 yr −1 , with agricultural land being the most susceptible to erosion. There was a steady increase in mean soil erosion rates across most of the model ensembles from the 2030s to 2080s, with only four of all the 16 future scenarios recording a decline in erosion rates. All of the GCMs showed a positive change in average annual rainfall erosivity under all scenarios, except for MIROC-5-RCP2.6 in both the 2030s and 2080s time slices. The cumulative average of rainfall erosivity and soil erosion rates for all of the GCMs in the 2030s (both RCP2.6 and RCP8.5 combined) showed an increase of 39.9% and 29%, respectively, compared with the baseline climate. Similarly, the aggregate average of rainfall erosivity and the soil loss rate for all the GCMs projected increases of 61.1% and 59.87% respectively, in the 2080s. The highest positive change (217.23%) of the soil loss rate over Kenya belonged to CCCMA_CANESM2 RCP 8.5 in the 2080s, while the minimum (−14.29%) was projected under the MIROC-5 RCP 2.6 in the 2030s period. Future climatic projections indicate that major expected increases in rainfall erosivity and soil loss rates are concentrated in the CER and RVR, which represent the main crop production regions of Kenya. In the Rift Valley Region, the relative variation of the mean soil loss rate rose up to 107.94% (CCCMA_CANESM2 RCP 8.5 (2080s). This study provides a good reference to policy makers to enable them devise adaptive soil erosion conservation measures.