Linking climate change to soil loss estimation in the Kosi river basin, India

Soil erosion has been the most serious problem for sheet and rill erosion, which thereby affects the pastoral production on red soil in India, covering an area of 72 Mha (Dhruvanarayana & Babu 1983), and according to UNEP (1997), 1,964.4 Mha of area globally is affected due to anthropogenic activities. Soil erosion is mainly caused by the rainfall as individual raindrops have a tendency to detach the soil particles and cause them to runoff down the slope if the terrain is not flat (Angulo-Martínez et al. 2012). The high intensity of rainfall may cause erosivity to create saturation and lagooning which leads to enhanced soil particle detachment (Salles et al. 2002). Precipitation in a higher amount was the major driver in the spatial and temporal rainfall–runoff erosivity (Amanambu et al. 2019). There is no linear relationship between erosivity and soil erosion but erosivity is highly dependent on raindrops, the intensity of rainfall and its duration (Salles et al. 2002).

The extremes of precipitation may impact directly or indirectly on human activities and environment (Arab Amiri & Mesgari 2018). The spatio-temporal variability of precipitation was well studied in recent research using various techniques (Martins et al. 2012; Arab Amiri & Mesari 2019). Arab Amiri et al. (2017) investigated the spatial variation of overall precipitation over northwestern Iran using time-series analysis for a period of 1991–2010. Cheng et al. (2019) analyzed the spatio-temporal characteristics of precipitation over the Tianshan Mountains, Central Asia for a period of 55 years (1961–2016). With the application of geographic information system (GIS) and the technique of time-series analysis, the spatio-temporal variation of climatic factors can be evaluated and visualized in a better way (Arab Amiri & Mesari 2016). Hence, in this study, the effect of rainfall characteristics on erosivity has been analyzed considering the spatio-temporal variation of rainfall.

Around 82 models are used for the estimation of soil, and some of them are Universal Soil Loss Equation (USLE) (Wischmeier & Smith 1978), Revised USLE (RUSLE) (Renard et al. 1997) along with the Water Erosion Prediction Project (WEPP) (Nearing et al. 1989) and European Soil Erosion Model (EUROSEM) (Morgan 2005). The RUSLE is usually preferred for the estimation of soil loss and is very much used to find the soil erosion risk and implemented regionally and globally (Borrelli et al. 2016). The RUSLE is a kind of deterministic approach and is based on the empirical model including parameters which might not be suitable as these parameters are based on the erosion in US-based studies despite remotely sensed data (Wischmeier & Smith 1978; Renard et al. 1997). Rainfall–runoff erosivity factor (R) is related to the climatic factors like rainfall (Yang et al. 2003), and cover management factor (C) is the quickest responsive parameter among all parameters (Estrada-carmona et al. 2017) that can be managed by farmers in order to control the erosion (Panagos et al. 2015).

Soil erosion gets affected due to climate change either directly or indirectly, as the direct effects depend on the amount and time of rainfall and, on the other hand, indirect effects depend on both the increased concentration of CO₂ in the atmosphere affecting the growth of crop and climate-driven changes in land use (Boardman 2001). The total precipitation and its distribution are important factors in the climatological studies (Arab Amiri et al. 2017). Renard & Foster (1983) expressed the climate as one of the factors that affect the erosion. The effect of climate change and the spatio-temporal variation of hydro-meteorological parameters are well discussed in the recent studies (Cao & Pan 2014; Arab Amiri & Mesgari 2017; Amanambu et al. 2019). Due to the rainfall erosivity, the rate of soil erosion is expected to change in accordance with the climate change (Nearing 2001). The trends of rainfall–runoff erosivity have been studied due to climatic parameters like rainfall for the past and future periods in the Himalayan region (Gupta & Kumar 2017). Estrada-carmona et al. (2017) argued that the uncertainty in the RUSLE model could be reduced with the help of better parameterization in the cover management factor and topography factor avoiding extreme soil losses by targeting the practices of conservation of soil in the region where both factors enhance soil loss. It is essential to quantify the soil erosion considering past and future variation of rainfall–runoff erosivity and crop management and identifying areas prone to erosion due to climate change in the Kosi river basin.

According to the IPCC (AR5, 2013) report, the percentage change of mean precipitation trend varies from +10 to +30% by 2100, which is very much uncertain in the region. Kumar et al. (2011) projected the Indian summer monsoon rainfall and found that there will be an increase of 20% by the end of this century, but the percentage rate of change varied according to the input from various general circulation models (GCMs) as well as their corresponding scenarios. Kumar et al. (2020) used three representative concentration pathways (RCPs: RCP2.6, RCP4.5 and RCP8.5) to find the trend of rainfall in the Sone command area in Bihar, India.

Various studies looked over a change in erosivity in many ways (Almagro et al. 2017; Amanambu et al. 2019) as well as the combination of erosivity and soil erosion using GCMs in different regions. Various GCMs have shown the changing pattern of rainfall which affects the soil erosion as well as rainfall–runoff erosivity (Plangoen et al. 2013). Precipitation data obtained from GCMs can be used to project the rainfall erosivity change for future periods (Hoomehr et al. 2016) under various RCPs as per the IPCC report (AR5, 2013).

Trend analysis guides the possibility of trends whether it is present or absent. The nonparametric method of trend analysis is one of the most common approaches for trend analysis since the parametric method of approach has the limitation that the dataset must be normally distributed and independent (Tamaddun et al. 2019). Kundu et al. (2014) has applied the Mann–Kendall (MK) trend analysis approach to obtain the trend of rainfall for the whole of India with the data of 141 years (from 1871 to 2011) with seven major regions. Bezak et al. (2020) utilized the MK trend test to find the trend of station-wise rainfall erosivity. Arab Amiri & Mesgari (2018) investigated trends spatially along the latitude and longitude in the time slice of 10 years having realized a spatial analysis of precipitation extremes in Iran.

Land use is one of the factors which affect the soil erosion due to anthropogenic activities. There may be an increase in soil erosion due to less vegetation cover, and this relationship is the main reason to be a part of soil erosion studies (Cotler & Ortega-larrocea 2006; Solaimani et al. 2009). Hence, it is essential to model change in land use which ultimately can be used for soil erosion prediction.

Due to high sediment yield (Mishra & Sinha 2020), it is essential to compute soil erosion due to the impacts of climate change in the Kosi river basin. Climate change is one of the major contributors of soil erosion in the lower Himalayan region as well as other basins in the world which are prone to frequent floods. Due to the accelerated and intense erosion of the Kosi river basin, various studies (Ganasri & Ramesh 2016; Uddin et al. 2016) have focused only on the severity and prioritize the soil erosion level. Amanambu et al. (2019) quantified the nature of spatial and temporal variations of rainfall erosivity using various GCM data; however, the land-cover changes are lacking in their study, which play a vital role in the soil loss estimation. As yet, no studies have been well documented for soil loss for the future period considering the change in rainfall erosivity as well as land use in the Kosi river basin. Hence, in this paper, we aim to quantify the spatial and temporal rainfall erosivity and their nature under RCP2.6, RCP4.5 and RCP8.5 scenarios for four different GCMs (CanESM, CCMC, MPI and CSIRO) for the projected time period (2020–2100). We quantified annual soil loss utilizing land-cover values and rainfall erosivity derived from average rainfall data of different GCMs for the projected time period in the three time slices (2020–2039, 2040–2069 and 2070–2100). We also estimated the soil loss for each RCP (RCP2.6, RCP4.5 and RCP 8.5) and compared it with that of the soil loss obtained from average rainfall data of different GCMs for the period of 2020–2100. These efforts contribute to the development of suitable and effective strategies to reduce the soil loss and related environmental problems.

The Kosi river is one of the major tributaries of the Ganga river. This river is known as the ‘Sorrow of Bihar’, as it causes very frequent floods in Bihar. The main cause behind the selection of the Kosi basin is that it has been a substantial issue of rapid bank erosion and sediment yield due to which shifting of channels takes place within the basin (Mishra & Sinha 2020). Kosi river flows in the North Bihar plain and covers an area of around 11,410 km². The Kosi river basin is located between latitude 25°19′18″–26°43′30″N and longitude 87°4′35″–87°12′32″E. It originates from the Himalayan region in Nepal and occupies a large area in Tibet and Nepal. The river enters the Bihar (India) region near Bhimnagar and joins the Ganga river near Kursela, Katihar district, Bihar. The major tributaries of the Kosi river are Bagmati, BhutiBalan, and Kamlabalan. The Kosi basin comprises 27 sub-watersheds. Ten districts of Bihar are covered in the drainage basin.

SRTM (Shuttle Radar Topography Mission)-based digital elevation model (DEM) has 90 m resolution for the year 2019, and Landsat 8, Landsat 7 and Landsat 4-5 were obtained from the earth explorer portal (https://earthexplorer.usgs.gov/). The daily rainfall data for a period of 33 years, i.e., from 1985 to 2017, were obtained from the India Meteorological Department (IMD), Pune for seven available rain gauge stations (Bhimnagar, Nirmali, Bhadurganj, Birpur, Murliganj, Galgalia and Kursela) as shown in Figure 1. The highest average annual rainfall is 2,711.95 mm obtained in Galgalia, and the least is 1,249.13 mm in Bhimnagar. The contour plots for average monthly observed rainfall for the baseline (1985–2017) are shown in Figure 2 for all rain gauge stations in the Kosi river basin. It is observed from Figure 2 that most of the concentration of rainfall occurs in May–October at all stations. At all stations, the least rainfall is found in the year 2009–2010, of which Bihar faced drought during this tenure (Kishore et al. 2014). The GCM (CanESM, CCMC, MPI and CSIRO) data for three RCPs, such as RCP2.6, RCP4.5 and RCP8.5, were obtained from the World Climate Research Programme (WCRP) based on the fifth phase of Coupled Model Intercomparison Project (CMIP5) mentioned in Table 1.

Soil erosion has been a major problem, as it aggravates the natural condition of the field. The problems resulted from the soil erosion like sediment deposition in the rivers, lakes and their estuaries have been ongoing and long-term issues throughout the world (Ganasri & Ramesh 2016). There are many models used for the estimation of soil loss, including USLE (Wischmeier & Smith 1978), RUSLE (Renard et al. 1997), WEPP (Nearing et al. 1989) and EUROSEM (Morgan 2005). Among these models, the RUSLE model was preferred for the estimation of soil erosion, as it needs remote sensing data for generating land-use and land-cover maps, management practices and types of soil and their properties. In the present study, the RUSLE model has been used for the estimation of annual soil erosion for the past as well as projected period. The parameters included in the RUSLE model are DEM, types of soil, rainfall data and land cover. Figure 3 shows the exaggerated view of the study.

As the value of S changes, the trend of the analysis also varies. An increasing trend means to say that there is a high positive value of S, and vice-versa. However, the significance of the trend can also be computed with the help of probability with a sample size n. The P-value is the probability value corresponding to the MK statistical parameter. The P-value actually tests whether the hypothesis is correct or not. An incorrect hypothesis means that there is a significant trend and it is declined. If the P-value is less than 0.05, it signifies that there is a significant trend; if the corresponding tau-value is positive, then it has an increasing trend and if the tau-value is negative, then it has a decreasing trend.

In this study, the MK test was adopted using the R-Studio software to find the trend of observed rainfall data (1985–2017). It was found that four stations have an increasing trend, while three stations have no sign of a significant trend. Figure 4 shows trends and spatial variation of average annual rainfall which varies from 1,249.14 mm in Bhimnagar to 2,711.89 mm in Galgalia.

An inverse distance weighted (IDW) interpolation technique was used for the spatial interpolation in this study using ArcGIS10.3. It is an assumption in the IDW interpolation technique that the quantification at the unsampled point can be done with the help of the distance weighted average of the neighboring sampled point values (Xu et al. 2019). The method of the IDW interpolation technique is simple, efficient and easy to understand (Azpurua & dos Ramos 2010). Arab Amiri & Mesgari (2017) utilized the IDW interpolation technique to map the detected trends in annual rainfall. Many researchers used this interpolation technique to quantify the rainfall erosivity in order to evaluate soil erosion (Chen & Liu 2012; Huang et al. 2013; Xu et al. 2019).

In this study, the average annual rainfall for a period of 1985–2017 (baseline period) as well as GCM precipitation data for the projected period (2020–2100) was found above 850 mm, and hence Equation (5) was used for the calculation of R. The rate of soil erosion is more sensitive to the occurrence of rainfall in the catchment (Jain et al. 2001; Dabral et al. 2008). The daily rainfall works as a better indicator for the sediment yield, and an additional advantage is that it can be used for the characterization of seasonal distribution of sediment yield. The advantages of using annual rainfall include its ready availability, ease of computation and greater regional consistency of the exponent (Shinde et al. 2010). The R factor has been calculated for the baseline using the IDW interpolation technique (Figure 5) and found that the value of R factor ranges from 5,470.35 to 27,471.6 MJ mm ha⁻¹ h⁻¹ y⁻¹. This can be visualized from Figure 5 that the R factor is less in the Murliganj, Bhimnagar and Kursela rain gauge stations and gradually increases to the Birpur rain gauge station but drastically changes in the Galgalia rain gauge station since Galgalia has the highest observed average annual rainfall.

The length and slope steepness factor is related to topography or length and slope factor, and it is defined as the ratio of soil loss under the given condition to that of a fallow land with 9% steepness and length of 22.6 m (Dabral et al. 2008). For a given site, the LS factor reflects the susceptibility of topographic erosion and is computed from the DEM (30 m) using the Arc Hydro tool. The slope length is defined as the stretch from where the surface flow starts to the point where the deposition starts taking place or where runoff water starts entering into the well-defined channel (Ganasri & Ramesh 2016). Hence, there is an increment in the soil loss as the slope increases. The steepness of slope has a higher effect on the loss of soil as compared to that of the length of the slope. The steeper the slope, the higher the erosion.

The inputs used for the calculation of the LS factor are flow accumulation grid and percentage slope (Olorunfemi et al. 2020).

The C factor is calculated for the year 2019 and is shown in Figure 8. The C factor ranges from 0.289469 to 1.19271. The C factor indicates the effect of cropping and its management practices on soil erosion. The higher the value of the C factor, the poorer will be the management of crops, and consequently, the higher will be the soil erosion, and vice-versa.

The support practice (P) factor is defined as the ratio of soil loss having specific support practice to the corresponding soil loss with inclination and declination in the ground (Renard et al. 1997). The value of the P factor varies from 0 to 1, where 0 reflects the good conservation practice and 1 reflects the poor conservation practice. A supervised classification was done for the preparation of the LU/LC map using the Landsat 8 image (October to December 2019) and verified by ground truthiness (between June and September 2019) and Google Earth software for proper classification. The value of the P factor (Table 2) was then assigned according to Yang et al. (2003).

Effects of climate change on rainfall and therefore rainfall erosivity factor have been discussed by many researchers (Gupta & Kumar 2017; Panagos et al. 2017). The rainfall erosivity (R) factor is a function of average annual rainfall and is evaluated for both the baseline period and the future period. For GCMs, the fifth phase of CMIP5 is a significant model improvement over CMIP3 which utilizes a new set of emission scenarios referred to as RCPs (Taylor et al. 2012). Kumar et al. (2019) utilized GCMs (CanEMS, CMCC and MPI) to downscale precipitation data in the Bagmati river basin, Bihar. In this study, global climate data from 4 GCMs (CanESM, CCMC, MPI and CSIRO) were used under RCP2.6, RCP4.5 and RCP8.5 scenarios (Kumar et al. 2020). Using these four GCMs, the projected precipitation data from 2020 to 2100 were obtained from the WCRP CMIP5 portal. These data were then used to find out the precipitation data according to the coordinates of rain gauge stations using the bilinear interpolation method. In this study, the linear bias correction technique (Kumar et al. 2020) was adopted for the precipitation data obtained from different GCMs (Equation (10)).

Climatic parameters are very well related to the vegetation index, and the C factor is associated with the variation in the NDVI. Recent studies (Liu et al. 2018, 2020) showed the variation of NDVI with respect to the climate change. In this study, the linear variation of NDVI was taken into consideration. The evaluation of NDVI for the projected period (at the end of 2039, 2069 and 2100) was done using the regression equation obtained from the past NDVI (1999, 2005, 2006 and 2009) data.

The average annual soil loss was evaluated using the RUSLE model (Equation (1)) for the baseline (1985–2017), and the spatial distribution of annual soil loss map was prepared as shown in Figure 9. For the evaluation of annual soil loss, the C factor for the year 2019 (October–December) and the rainfall data for the baseline (1985–2017) were utilized. The annual soil erosion (t ha⁻¹ y⁻¹) was categorized into five different levels as shown in Figure 9. The first level of soil erosion is less than 1.5 t ha⁻¹ y⁻¹, the second level is from 1.51 to 5.15 t ha⁻¹ y⁻¹, the third level is from 5.16 to 10.30 t ha⁻¹ y⁻¹, the fourth level is from 10.31 to 20.00 t ha⁻¹ y⁻¹ and the last and fifth level is greater than 20 t ha⁻¹ y⁻¹. It is found from Figure 9 that the eastern portion of the study area is severely affected by soil erosion, where the rain gauge station Galgalia is located (Figure 4) which has the highest observed rainfall as well as R factor. This may indicate that annual soil erosion is highly dependent on rainfall–runoff erosivity. At the location of the Birpur rain gauge station, annual soil erosion was found to be around in the area of third level. Nearby the Kursela rain gauge station, the annual soil erosion lies in the range of first level.

The percentage variations of average annual rainfall (three time slices) of all GCMs were calculated with respect to the baseline and plotted from 2020 to 2100. For better analysis, three time slices are considered, 2020–2039, 2040–2069 and 2070–2100, and are shown in Figure 10. It is found from Figure 10 that most of the GCMs have low variation apart from CanESM-8.5 (2070–2100). The maximum percentage change in precipitation was found to be 26.2% in the CanESM-8.5 model for the time slice 2070–2100. The least percentage change was found in MPI-2.6 in all time slices. The percentage change in precipitation for the time slice 2020–2039 was found to be negative in CSIRO for all scenarios. For the time slice 2040–2069, all models were having a positive percentage change in precipitation except for CanESM (RCP 2.6 and 4.5). The maximum variations for the time slice 2020–2039 and 2040–2069 were found to be −8.16 and 9.4% in the CanESM-4.5 and CanESM-8.5 models, respectively.

The MK trend test is generally recommended by the World Meteorological Organization (WMO) for the trend analysis of meteorological/hydrological data. In the recent studies, they have utilized the MK test for trend analysis (Xu et al. 2019; Kumar et al. 2020). The average of precipitation data was calculated for four GCMs with three time slices, and a rainfall map (Figure 11(a)) was prepared using the IDW interpolation technique. The MK trend test was done for finding the trend of future precipitation in the three time slices. It was found that all stations have no trend, except the Murliganj rain gauge station in the time slice of 2040–2069 (Figure 11(a)). The rainfall deficit was calculated, which is the difference in future rainfall obtained from the average of all GCMs and observed data (1985–2017), and a map was prepared as shown in Figure 11(b). In the figure, the rainfall deficit is the highest at the Murliganj station for all three time slices (Figure 11(b)). The rainfall deficit ranges from −129.359 to 142.917 mm. The minimum rainfall deficit was found at the Nirmali station.

Due to climate change, the rainfall pattern and its amount affect the erosivity. The impact of climate change has been visualized with the help of rainfall–runoff erosivity and is shown in Figure 12 for the projected period. Using GCM data from 2020 to 2100, the rainfall–runoff erosivity map was prepared for the three time slices (2020–2039, 2040–2069 and 2070–2100). The range of the rainfall–runoff erosivity variation was found to be from 1,137 to 3,441 MJ mm ha⁻¹ h⁻¹ y⁻¹. In the CanESM-2.6 model, rainfall erosivity decreases in the time slice of 2040–2069 as compared to 2020–2039, but again it increases in the time slice of 2070–2100. In the other two RCPs (RCP4.5 and RCP8.5) of CanESM, rainfall erosivity increases, and becomes especially severe in CanESM-8.5 under the time slice of 2070–2100, and the extreme value of all models was in the same model. In GCMs like CCMC-2.6 and CCMC-4.5, rainfall erosivity first increases in the time slice 2040–2069 but again it decreases in 2070–2100, but in CCMC-8.5, rainfall erosivity looks similar in the first two time slices; however, it decreases in the last time slice. In the time slices of 2020–2039 and 2040–2069, the rainfall erosivity map looks almost the same but it decreases in RCP2.6; however, there is an increase in the rainfall erosivity. Finally, rainfall erosivity increases in the first two time slices, and then again decreases in the last time slice, i.e., 2070–2100.

From the above results, it may be concluded that rainfall erosivity increases in the time slice 2040–2069 as compared to 2020–2039, and it is almost the same in the last two time slices in most of the models. The current study may help to find out rainfall erosivity for the projected period in the subtropical belt and the humid region like the Kosi river basin on the basis of rainfall occurrence.

The NDVI for the projected period was calculated using the linear regression method (Zhu et al. 2020) and is plotted in Figure 13. The C factor was then calculated at the end of three time slices, i.e., 2039, 2069 and 2100 for the projected period using Equation (9) and is tabulated in Table 3. From Table 3, the value of C factor goes on increasing for the projected years.

In Figure 13, the year 1999 was taken as a reference year or zero year for NDVI calculation. The predicted value of NDVI was used for the evaluation of the C factor. The calculated C factors are 0.7436, 0.7830 and 0.8622 at the end of the years 2039, 2069 and 2100, respectively. These C factor values were used for the calculation of soil erosion for all three projected time slices. The higher value of C factor resulting from browning vegetation will lead to a higher rate of soil erosion (Liu et al. 2020). A rise in C factor could potentially be affected by deforestation or the burning of agricultural fields.

The consequences of erosive change as a result of past and future climate changes are demonstrated in Table 4, illustrating differences in possible annual soil loss for past climate and future scenarios. The percentage variations of soil loss have been evaluated taking the combined average of annual rainfall data of all stations for the three time slices (Amanambu et al. 2019). In the present study for the estimation of soil loss due to climate change, the variations of the R factor and the C factor were taken into account while keeping the other factors (LS, K and P factors) constant. It is shown in Table 4 that some scenarios show an increase in soil erosion, whereas other scenarios show a decrease. All changes were taken with respect to the baseline. In all scenarios, the percentage change in annual soil erosion was found to be negative in the time slice 2020–2039. The least percentage change in soil loss was found in CanESM-4.5, −24.00% in the time slice 2020–2039, and the corresponding percentage change in R factor was also the least. The highest percentage change in the R factor was found to be 8.93%, and the corresponding percentage change in soil loss was also highest in the time slice 2020–2039. In the time slice 2040–2069, the percentage change in R factor was found positive in all scenarios except apart from CanESM-2.6 and CanESM-4.5. The least percentage change in soil loss was found to be −12.00% in the model CanESM-2.6, and the maximum for CanESM-8.5 (15.00%). The percentage change in the R factor was also found to be highest for the same model, at 20.48%. In all scenarios of the last time slice, i.e., 2070–2100, except for CCMC-8.5, MPI-4.5 and MPI-8.5, showed a decrease in percentage soil erosion with values of −2.00, −4.00 and −1.00%, respectively. The highest percentage increase in soil erosion (71.00%) was found in CanESM-8.5 for the time slice 2070–2100, which is the highest of all scenarios, and the corresponding percentage change in R factor is also the highest, which was also found in the projected R factor (Figure 11). It was found that there was an increment in the percentage of soil erosion change despite the increase or decrease in their percentage change in rainfall–runoff erosivity factor. All the scenarios in the time slice 2040–2069 have an increasing percentage change in rainfall–runoff erosivity factor, except the two scenarios CanESM RCP2.6 and CanESM RCP4.5. The average of percentage change in R factor varies from −4.15 in the time slice 2020–2039 to 4.45 in the time slice 2070–2100. This change in R factor may be due to change in rainfall pattern for the future period. The average of all models with three time slices was seen to be increasing from −13.03 to 10.39%, which can be justified by the work conducted by Yang et al. (2003) that there will be an approximate increase in average soil erosion of 9% globally by 2090. It is evident that all models and their scenarios do not have the same direction of change and magnitude of erosivity, and this can be attributed to a degree of uncertainty in GCMs.

Factors such as rainfall–runoff erosivity and crop management were taken as variables in the RUSLE. The annual soil loss map was prepared for three time slices and is shown in Figure 14. The variation of annual soil loss can be seen in all the time slices, and it was found that soil loss ranges from 0 to 250 t ha⁻¹ y⁻¹. The highest annual soil loss was found in the time slice 2070–2100. It was found that there is a significant increase in annual soil loss between the time slice of 2020–2039 and 2040–2069. A slight difference in the annual soil loss between the time slices 2040–2069 and 2070–2100 is noticed. It was found that the maximum soil loss was in the eastern region of the Kosi basin, which is quite similar to the observed annual soil loss. As R factor was more in the eastern region, hence the soil loss was also found more in the eastern region which reflects that the rainfall erosivity has a tremendous role in annual soil loss.

In this study, there is an increasing trend of estimated annual soil erosion in the Kosi river basin in the three time slices taken. The results observed for the annual soil erosion from Table 4 and Figure 14 do not have a uniform indication of precipitation for individual GCMs. It has been observed from the results in Table 4 for different GCMs that rainfall is not linearly related to the soil erosion in the lower Himalayan region. Thus, cover management is one of the major factors which have to be considered for the estimation of annual soil loss for projected periods.

Annual soil erosion does not only depend on the higher rainfall but also the duration of rainfall, its intensity, size of the raindrops and the management practices of the field. The Kosi river basin is a highly sensitive area with reference to the flood frequency and the sediment deposition within the basin. Soil erosion is expected to increase in the coming years due to climate change and lack of control over erosion control practices. Anthropogenic activities are one of the main reasons for not controlling the soil erosion and this is unavoidable.

The R factor and the corresponding soil losses based on each scenario (RCP2.6, RCP4.5 and RCP8.5) have also been estimated for the future time period (2020–2100). The average of annual precipitation data for RCP2.6, RCP4.5 and RCP8.5 from four GCMs (CanESM, CCMC, MPI and CSIRO) has been utilized in order to estimate the R factor. The C factor has been taken for the year 2019 for the estimation of annual soil loss. Figure 15(a) and 15(b) shows the spatial variation of R factor and the corresponding annual soil loss in the three RCPs. It is visualized from Figure 15(a) that the most severe effect of the R factor is in the RCP2.6 scenario, whereas the other two scenarios (RCP4.5 and RCP8.5) have an exiguous difference. However, RCP8.5 has a slightly higher effect of R factor in the eastern region of the Kosi river basin as compared to the RCP4.5 scenario. Despite the effect of R factor, the soil losses in the three scenarios are apprehended to be similar and it was found that in all the cases, the eastern region of the Kosi river basin is severely affected.

Soil loss based on the average of all GCMs and each RCP is tabulated in Table 5 with five classes based on the histogram and standard deviation of images (Guo et al. 2018) obtained in ArcGIS10.3. Soil losses are classified as 0–10, 11–32, 33–64, 65–126 and 127–250 t ha⁻¹ y⁻¹ for the future time period. The soil loss area in all cases (GCMs and RCPs) is the maximum in the 0 –10 t ha⁻¹ y⁻¹ range, whereas it is minimum in the 127–250 t ha⁻¹ y⁻¹ range. In all ranges, there are meager changes in the area of soil erosion.

The Kosi river basin is considered as a sensitive basin in terms of frequent floods in Bihar. This research has been taken into account for long-term variations in annual soil loss, the trend of projected precipitation and its variation. Variations in major factors like rainfall–runoff erosivity factor as well as crop management factor are unavoidable factors in the basin parameters.

The rainfall–runoff erosivity factor is calculated for three time slices, 2020–2039, 2040–2069 and 2070–2100, using projected rainfall data downscaled from four GCMs: CanESM, CCMC, MPI and CSIRO for RCP2.6, RCP4.5 and RCP8.5 scenarios. The MK trend test shows that except for the Nirmali station (increasing trend) in the time slice of 2070–2100, no station is showing any significant trend in rainfall. The rainfall–runoff erosivity factor and thereby the soil loss are computed using the average rainfall value of all GCMs for all time slices. The results show that there is an increasing trend in the rainfall–runoff erosivity factor as well as the annual soil loss, especially in CanESM-8.5 in the time slice of 2070–2100, which is the last time slice in this study. The average percentage changes of R factor are found as −4.15, 4.19 and 4.45 for the time slice of 2020–2039, 2040–2069 and 2070–2100, respectively, which is found to be increasing. This change may reflect the change of rainfall pattern for the projected period. Beside the R factor, the C factor was also taken as a variable for the computation of annual soil erosion for the projected period with the help of linear regression. The computed C factor values are 0.7436, 0.7830 and 0.8622 for the years 2039, 2069 and 2100, respectively. The increasing C factor values are due to poor supporting factor in the projected period, which may be due to deforestation. The corresponding percentage change in annual soil erosion is also calculated and an increasing trend for the selected time slices was found. The average value of soil loss from all GCMs was seen to be increasing from −13.04 to 9.88% for three time slices.

Besides estimating the soil loss for the three time slices, the effects of rainfall–runoff erosivity and the corresponding soil loss based on each RCP (RCP2.6, RCP4.5 and RCP 8.5) have also been visualized. It is observed that there are meager differences among soil erosion maps in the three RCPs. However, the effect of the R factor may be observed more in the RCP2.6 scenario but less in the other two scenarios (RCP4.5 and RCP8.5). The area of soil erosion for all GCMs and each RCP have been estimated and compared and it was found that in all ranges of soil loss, there are meager changes in the area of soil erosion.

The present research evaluated the long-term variations of annual soil erosion resulting in projected rainfall erosivity as well as cover management. The consequences would help to better control measures in the highly affected regions in areas like the Kosi river basin. These observations are essential for further studies since this research gives an outlook on the relationships of soil loss with rainfall, soil characteristics and land-cover change.

The authors thank the IMD, Pune for providing data for the analysis of the study.

The authors confirm that the content of this article has no conflicts of interest.

All relevant data are included in the paper or its Supplementary Information.