Sub-basin prioritization for assessment of soil erosion susceptibility in Kangsabati, a plateau basin: A comparison between MCDM and SWAT models.

Kangsabati basin located in tropical plateau region faces multiple problems of soil erosion susceptibility (SES), soil fertility deterioration, and sedimentation in reservoirs. Hence, identification of SES zones in thirty-eight sub-basins (SB) for basin prioritization is necessary. The present research addressed the issue by using four multi-criteria decision-making (MCDM) models: VlseKriterijumska optimizacija I Kompromisno Resenje (VIKOR), technique for order preference by similarity to ideal solution (TOPSIS), simple additive weighing (SAW), compound factor (CF). To determine the best fitted method from MCDM for erosion susceptibility (ES), a comparison has been made with Soil and Water Assessment Tool (SWAT), where fifteen morphometric parameters were considered for MCDM, and meteorological data, soil, slope and land use land cover (LULC) were considered for SWAT model. Two validation indices of percentage change and intensity change were used for evaluation and comparison of MCDM results. With SWAT model performance, SWAT calibration and uncertainty analysis programs (CUP) was used for sensitive analysis of SWAT parameters on flow discharge and sediment load simulation. The results showed that 23, 16, 18 SB have high ES; therefore they were given 1 to 3 ranks, whereas 31, 37, 21SB have low ES, hence given 38 to 36 rank as predicted by MCDM methods and SWAT. MCDM validation results depict that VIKOR and CF methods are more acceptable than TOPSIS and SAW. Calibration (flow discharge R2 0.86, NSE 0.75; sediment load R2 0.87, NSE 0.69) and validation (flow discharge R2 0.79, NSE 0.55; sediment load R2 0.79, NSE 0.76) of SWAT model indicated that simulated results are well fitted with observed data. Therefore, VIKOR reflects the significant role of morphometric parameters on ES, whereas SWAT reflects the significant role of LULC, slope, and soil on ES. However, it could be concluded that VIKOR is more effective MCDM method in comparison to SWAT prediction.

Introduction

Soil is the most sustaining natural resource for all living organism, and it determines geomorphic processes widely (Keesstra et al., 2016; Ameri et al., 2018; Hembram and Saha, 2018). Sustainable agricultural development and natural resource utilization both are largely determined by soil erosion; furthermore, this erosion is controlled by several hydrological functions in the watersheds worldwide (Molla and Sisheber, 2017; Arabameri et al., 2018). Weathering process leads to disintegration of soil particles from allied rocks and minerals, and eroded materials are entrained, transported and deposited to another place by flowing water (Masselink et al., 2017; Hembram and Saha, 2018). Water is an important agent for soil erosion, land degradation and deterioration of soil fertility from the topmost layer, and thus soil erosion hampers the plant growth, sedimentation in river valley and reservoir, delta generation in estuarine sites, which gives rationality for sustainable human societies (Sharma et al., 2017; Ameri et al., 2018). Indeed, ES at regional and global level, ended with several types of environmental and economic consequences throughout the world (Gayen et al., 2019); where such direct measurement techniques like erosion pins, silt fences, and rainfall simulation are applied to determine the critical erosion sites (Nasre et al., 2013; Ameri et al., 2018). In this context, analysis of drainage basin or watershed become more relevant, and its several geo-environmental elements like lithologic, geomorphic, morphometric, land use land cover (LULC), and soil play significant role in sustainable watershed management (Chauhan et al., 2016; Balasubramanian et al., 2017; Bhattacharya et al., 2019b). Moreover, drainage basin play dominant role to determine the ES, surface runoff, and sedimentation in channel, all of which are caused by various erosional agents of water and wind, and thus are expressed various landscape processes that signifies the effective morphometric characterization (Patel et al., 2012, Patel et al., 2013; Ameri et al., 2018; Hembram and Saha, 2018). Mathematical measurement of morphometric parameters including three major aspects i.e. linear, areal and relief, are prime requirements for prioritization of sub-watersheds or hydrological units that can have considerable effects on land and water conservation in a drainage basin (Aher et al., 2014; Rahaman et al., 2015; Ahmed et al., 2018); where effective morphometric parameters are well analysed the runoff volume, geomorphic landform, soil physical properties, along with basin ES (Keesstra et al., 2014; Ameri et al., 2018; Masselink et al., 2017).

Another three major dominant parameters like geology, topography, and climate control the spatial distribution of flow system geometry, drainage density and drainage patterns in a watershed (Mesa, 2006; Ameri et al., 2018). In recent times, many research works showed that advance effective remote sensing and GIS data like satellite imagery, digital elevation model, plays wide role in extraction of morphometric variables for watershed prioritization, as well as soil and water conservation (Biswas et al., 1999; Chatterjee et al., 2014; Magesh et al., 2013; Rai et al., 2017). Indeed, sub-basin prioritization signifies the geo-environmental characteristics in each hydrological units following the key driving variables: morphometric parameters (Aher et al., 2014; Balasubramanian et al., 2017), lithological setup (Chauhan et al., 2016), LULC (Altaf et al., 2014), groundwater potentiality (Deepika et al., 2013), ES etc. (Ameri et al., 2018; Hembram and Saha, 2018).

For quantification of ES, many researchers used multi-criteria decision making (MCDM) techniques, which are applied for resolving decision-making with high accuracy decision that are widely accepted for accurate decision, where decision rank depends not only on single criterion, but also taken some other significant criteria (Georgiou et al., 2015; Govindan and Jepsen, 2016; Mulliner et al., 2016). Several MCDM models including AHP (Analytical Hierarchy Process), ANP (Analytic Network Process) and VIKOR are also used for accurate decisions (Saha, 2017); which are now combined as more effective tool for prioritization of sub-basin along with the integration of fuzzy logic techniques (Aher et al., 2014; Kharat et al., 2016). Spatial distribution of soil erosion vulnerability entirely depends on fundamental parameters that are estimated with the link of systematic knowledge from alternative or effective parameters, as a result, high complexity forms are generated, which also created such hindrance on effective estimation of potential soil erosion and accurate sub-basin priority (Abdul Rahaman et al., 2015; Hembram and Saha, 2018). In contrast, for ES, application of MCDM methods of SAW, VIKOR, TOPSIS, and CF helps to determine the sub-basin prioritization using several morphometric parameters in a drainage basin (Ameri et al., 2018). In India, most effective MCDM methods such as fuzzy logic and analytical hierarchy process (FAHP) are applied to determine the soil erosion sensitive zones for sub-basin prioritization, where adopted methods considered such morphometric variables for effective control strategies in order to help soil conservation (Mekonnen et al., 2017; Hembram and Saha, 2018). However, the previous studies show that MCDM methods give different results in the same watershed following the previous studies, and create confusion among the researchers to take rational decision, on the best method for their study area.

On the other hand, several hydrological parameters in a drainage basin, including runoff volumes, sediment delivery, sediment yield, potential evapo-transpiration and ES are assessed by various hydrological models (Kamble, 2001; Sridhar et al., 2018; Bhattacharya et al., 2020a, Bhattacharya et al., 2020b); Soil and Water Assessment Tool (SWAT) as one of them which automatically extracted sub-basins boundary delineation (Tripathi et al., 2003, Tripathi et al., 2006), and also estimate soil erosion (Ameri et al., 2018), sediment delivery (Gassman et al., 2007), and all other hydrological parameters at sub-basin level (Vigiak et al., 2017; Anand et al., 2018; Markhi et al., 2019).

Therefore, this study has made an effort to assess the accurate sub-basins prioritization of a tropical plateau basin in order to ES in Kangsabati basin, India using both models i.e. MCDM models (SAW, VIKOR, TOPSIS, and CF) that is prepared from morphometric parameters, and integrated hydrological model as SWAT that is prepared from climatic characteristics, LULC, slope and soil parameters, and then make a comparison among MCDM methods furthermore, effective MCDM methods are selected by applying of SWAT model for estimation of soil erosion susceptibility in each sub-basin units. Additionally, this research offers two kinds of the validation process for better accuracy of MCDM model performance, where first kind of validation has been made through the intensity change, and second kind of validation is percentage change of basin rank of these methods, and best MCDM method is validated by SWAT model in order to similarity percentage of erosion priority rank. This study also attempts to identify the critical sub-basins, where effective conservation strategy is needed to resist soil erosion or land degradation, using complex decision approach.

Study area

Kangsabati basin is a tropical plateau basin with a total area of 9658km2, which is extended from 21°45′ N to 23°30′ N latitudes and 85°45′ E to 88°15′ E longitudes (Fig. 1a, b). The basin spreads over four districts, namely, Purulia, Bankura, Paschim Midnapore and Purba Midnapore. In order to lithological set up, oldest Archaean rock formation mainly granite-gneiss in upper basin makes elevated plateau proper in Chota Nagpur (641 m), and laterite dominance in middle basin makes undulating plateau fringe during pre-Cambrian era, whereas laterite track of Tertiary-Quaternary formation in lower floodplain (5 m) allow huge alluvial deposition during Pleistocene to Holocene era (Ghosh and Guchhait, 2015). Characteristic of climate in this basin, is characterized by southwest tropical monsoon climate, where mean rainfall ranges from 300 mm (drought prone) to 1650 mm (humid or very high rainfall zone). Drainage pattern in this study area is mainly dendritic to sub dendritic nature, where 1st order (5321), 2nd order (1254), 3rd order (298), 4th order (68), 5th order (22) and 6th order (8) are created hierarchical stream orders set-up following Strahler's method. Morphometric status of Kangsabati basin demonstrated that there has significant correlation with drainage, geology, soil, geomorphic landscape and land cover in order to ES (Bhattacharya et al., 2019b). Weaken soil profile, low vegetable cover and prominent rill gullies formation makes high ES in upper basin, whereas dominant laterite tract, steep slope, and agricultural practices creates mild susceptibility in undulating intermediate basin, but dense vegetation, presence of conservation practices resisted ES in the lower basin (Bhattacharya et al., 2020b).

Fig. 1

Study area: a sub basin boundary demarcation, b SWAT drainage outlet.

Materials and methods

Data sources

Geospatial data sources i.e. Survey of India (SOI) 1:50000-scale topographical maps (73I/10,11,12,13,73J/9,15,73N/2,3,4,5,7,8), 30 m spatial resolution Shuttle Radar Topography Mission (SRTM) Digital Elevation Model (DEM), 80 Ground Control Points (GCP) from Google Earth, were taken for extraction of morphometric parameters at sub-basin level. In SWAT model, various weather datasets including temperature, humidity, wind speed, solar radiation at 12 grid points in Kangsabati basin were automatically collected from climate forecasting system reanalysis (CFSR) world climatic database (1970–2014) (http://globalweather.tamu.edu/). While rainfall database has been modified with thirty five years average dataset (1980–2015) of sixteen local rain gauge stations in this basin as provided by Indian Meteorological Department (IMD), for getting better accuracy of model. Average rainfall datasets were further fed into SWAT database through the incorporation with nearest CFSN grid points. In term of spatial resolution scale, temperature data was only available with coarser resolution scale of 1″ × 1″ (approximately 110 km × 110 km), but other weather databases were available with finer resolution scale of 0.25″ × 0.25″ (approximately 27.5 km × 27.5 km), respectively (Himanshu et al., 2019). Furthermore, DEM, Thematic Mapper (TM), Landsat Operational Land Imager (OLI) image (2016), and soil sheet (Indian Council of Agriculture Research) were used to classify the topographical expression, spatial pattern of LULC and soil distribution map in Kangsabati basin. Stream hydrological datasets like sediment concentration, sediment load, sediment transport rate etc. are not available in this basin except flow discharge is available in Mohanpur and Mukutmonipur station, however, to adjust SWAT-CUP calibration and validation, flow discharge (record data) and sediment load (field estimation) at monthly time-scale, were collected at Mohanpur station during 2010–2015. Fig. 2 represents methodological flow chart to obtain the rank of ES in each sub-basin following different MCDM methods and SWAT model.

Fig. 2

Methodological flow chart.

Drainage estimation

Data processing phase involved with two major steps, where first step comprises on scanning georeferences and image rectification. After that, rectified images were converted into resembled sheet, passes with masking and mosaic processes, derived from Universal Transverse Mercator Projection of WGS 1984 (45°N) using Arc GIS 10.2. Second step comprises of boundary demarcation for extraction of stream number and its contributing area of specific outlets, flow direction, and flow accumulation in each DEM pixel, using ArcSWAT tools. In order to attain high accuracy value, SRTM DEM was used for comparison and clarification of sub-basin boundary with topographical maps, and extracted three dominant morphometric aspects i.e. basin shape, linear and landscape. Table 1 presents the measuring techniques of nineteen major morphometric parameters following such literature. After extracting morphometric parameters, four MCDM methods of TOPSIS, VIKOR, SAW, and CF were applied to predict the sub-basins priority for estimating ES in Kangsabati basin at sub-basin level, retrieved from their validated results of percentage change and intensity change. In contrary, SWAT model was applied to estimate the flow discharge and sediment load, as well sediment yield for the assignment of actual sub-basin priority rank, and then make a comparison between effective MCDM methods and SWAT model.

Table 1

Computation method of linear, areal and relief morphometric parameters.

Aspects	Formula	Where	Sources
Linear aspect
Basin length (L)	L = 1.32A0.363	L = basin length (km), A = area of the basin (km2)	Nooka et al. (2005)
Stream order (u)	Hierarchical rank		Strahler (1964)
Mean stream length (Lsm)	L_sm=LuNu	Lsm = mean stream length, Lu = total stream length of order ‘u’ Nu = total no. of stream segments of order ‘u’	Strahler (1964)
Stream length ratio (Rl)	R_L=LuLu−1	RL = stream length ratio, Lu = total stream length of order ‘u’, Lu-l = the total stream length of its next lower order	Horton (1945)
Mean bifurcation ratio (Rbm)	Rbm = average of bifurcation ratios of all orders		Strahler (1952)
Drainage density (D)	D=L_uA	D = drainage density, Lu = total stream length of all orders, A = area of the basin (km2)	Horton, 1932, Horton, 1945
Stream frequency (Fs)	F_s=∑N_uA	Fs = stream-frequency, ∑Nu = total no. of streams of all orders, A = area of the basin (km2)	Horton (1945,1932)
Drainage texture (T)	TS = D × Fs	T = drainage-texture, D = drainage-density, Fs = stream frequency	Horton (1945)
Constant of channel maintenance (C)	C=1D	C = constant of channel maintenance, Dd = drainage density	Schumm (1956)
Length of overland flow (Lo)	L_o=12D	Lo = length of overland flow, Dd = drainage density	Horton (1945)
Infiltration number (If)	If = Fs × D	Fs = stream frequency; D = drainage density	Faniran (1968)
Shape aspect
Form factor (Ff)	F_f=AL^2	Ff = form factor, A = area of the basin (km2), L = basin length (km)	Horton, 1932, Horton, 1945
Circularity ratio (Rc)	RC=4πAP^2	Rc = circularity ratio, π = 3.14, A = area of the basin (km2), P = perimeter (km)	Miller (1953); Strahler (1964)
Shape factor (Bs)	BS=L^2A	Bs = shape factor, L = basin length (km), A = area of the basin (km2)	Horton (1932)
Compactness co-efficient (Cc)	C_c=0.282IPA^0.5	Cc = compactness coefficient, P = perimeter (km), A = area of the basin (km2)	Gravelius (1914)
Relief aspect
Basin relief (R)	R = H − h	R = basin relief, H = maximum elevation in meter, h = minimum elevation in meter	Schumm and Hadley (1961)
Relief ratio (Rr)	Rr = R/L	Rr = relief ratio, R = basin relief, L = longest axis in kilometre	Schumm (1956)
Ruggedness number (Rn)	R_n=R×DK	Rn = ruggedness number, R = basin relief, Dd = drainage densityK = a conversion constant 1000 when relative relief is expressed in meter and drainage density in kilometre/square kilometre.	Schumm (1956)
Mean slope (S)	s=RA×100	R = basin relief; A = basin area (km2)	Nautiyal (1994)

Multi-criteria decision-making models (MCDM)

VIKOR model

Vlse Kriterijumska Optimizacija I Kompromisno Resenje (VIKOR) is used for comprising of alternative conflicting criteria based on ranking sets under optimizing complex systems, introduced by Opricovic and Tzeng (2004). Performance of VIKOR method as closeness of ideal is entirely depended on well organization following nine steps (Huang et al., 2009; Ameri et al., 2018). where Pij means element of normalized decision matrix and Rij mean i-th alternative sets in j-th criteria. where Tij means element of the weighted normalized decision matrix, Pij means element of normalized decision matrix and Wj mean computed weight criteria given by AHP model. where wi mean jth criterion weight that means relative importance of criteria, V⁎ j mean maximum Xij, and V⁎ j mean minimum V− j. where S⁎ mean minimum Sj, S mean maximum Sj, R⁎ mean minimum of Ri, R means minimum of Ri, and V as introductory weight strategy of maximum group utility (Sj) and maximum criteria (Rj). Q value generally ranges from 0 to 1; in terms of range value, v > 0.5 denotes that trend value of Q index reaches in majority rule.

First step mainly involved with comparison of decision matrix to summarize as follows:

Second step involved with computation of normalized decision matrix using Linear method as Eq. (1) (Ameri et al., 2018):

Third stage comprises of weight assignment for each requiring criteria using Analytical Hierarchical Process (AHP). VIKOR method also helps to re-evaluate all criteria, makes pair wise comparison, using Saaty's (1978) method in Expert Choice software.

After that, weighted normalized matrix has been computed through the multiplying normal matrix in each criterion weight, given by Huang et al. (2009) and Sanayei et al. (2010) as Eq. (2):

Fifth stage mainly involved with the determination of best value (X-j) and worst value (X– j) of all standard criteria functions; if j = 1, 2…n: where j-th criterion as most benefited criterion for the causes of maximum j-th value. This criterion becomes more relevant for this purpose. Then X- = max ƒ , ƒ̶  = min ƒ .

Purpose of this stage is to compute maximum group utility index (Si) and minimum individual regret index of opponent group values (Ri) following respective Eqs. (3), (4):

Computation the value of Qi following i = 1, 2…m as Eq. (7) (El-Santawy, 2012):

Alternative rank order has considered as sorting measurement from S, R and Q values, where suitable alternative values are assigned from least value of three major parameters.

Two conditions are used to determine the best alternative value from highest rank order of Q values as follow Eq. (6) (El-Santawy, 2012; Ameri et al., 2018):

C1 or Acceptable advantage:where A1, A2 means first and second alternative position in the ranking list; N is the number of alternative criteria.

C2 or Acceptable stability in decision making: In terms of satisfaction level, alternative rank values are assigned as rank 1 in Q, S and R or both of them. On the contrary or unsatisfactory level, ranking order of alternative value would be assigned as A1, A2 …Am where Am measured by following Eq. (7):

According to Opricovic and Tzeng, 2004, Opricovic and Tzeng, 2007, A1 and A2 must be selected for best solution, where outcome of C1 does not satisfy.

TOPSIS model

Technique for order preference by similarity to an ideal solution (TOPSIS) is a distance-based method, presented by Hwang and Yoon (2012). Principle of TOPSIS indicates that measurement of alternative choice has been computed by Euclidean distance from shortest distance of positive ideal solution (PIS) and largest distance of negative ideal solution (NIS). In this context, closeness coefficient is used to express the results of PIS and NIS distances. In order to closeness coefficient value, preferred alternative value is taken as a higher coefficient value from a set of alternatives (Liou and Wang, 1992; Kannan et al., 2009). This model has detected rating option following some calculation steps (Hwang and Yoon, 1981): where rij presents normalized decision matrix element and aij mean alternative performance of i-th under j-th criteria. Where rij mean normalized matrix elements and Wj mean weight of j-th criteria. where j mean benefitted criteria and J mean cost criteria, respectively. where Cli+ mean closeness coefficient, di+ mean positive ideal solution (PIS) and di- mean negative ideal solution (NIS).

This step estimated normalized decision matrix (rij) as Eq. (8):

This step computed as weighted normalized decision matrix (Vij) following Eq. (9):

This step computed the PIS and NIS using Eqs. (10), (11) (Ameri et al., 2018)

This step computed as separation measures from ideal solution using Euclidean distance (n-dimension). Separation from ideal distance of PIS and NIS measured as following Eqs. (12), (13)

Purpose of the final step, is detect closeness coefficient of alternatives from ideal solution ranges 0 to 1 using Eq. (14). Superior values are estimated from higher relative closeness alternative values.

SAW model

This is another leading MCDM method, which assigned score in each option, and obtained by aggregating values in different criteria. Relative weights are taken for assigning each criteria score using decision maker (Hwang and Yoon, 1981; Sargaonkar et al., 2011). Simple Additive Weighting (SAW) model estimated rating option following consecutive steps as. where Xij mean initial weights and M means criteria number.

Decision matrix is used for the determination of normalized initial matrix (Rij) of j-th criterion as Eq. (15) (Ameri et al., 2018):

This step involved mainly determined the relation as Rij = Rij min or Rij = Rij max

Two different relation can be obtained from normalized weight, where minimum value in efficiency index represent as Rij = Rij min, and maximum value in efficiency index represent as Rij = Rij max. where Wj mean criteria weight of Analytical Hierarchical Process (AHP) model.

Consistency index helps to determine as normalized value following Eq. (16) (maximum efficiency index) and Eq. (17) (minimum efficiency index):

The purpose in this step is determining the weighted normalized decision matrix (Vij) following Eq. (18):

Final step of SAW model involved mainly data integration, acquired from final score value in each option (Ai) following Eq. (19) (Ma et al., 1999; Ameri et al., 2018):

CF model

Compound factor (CF) is an important MCDM methods, which is used for subject base conversion phenomena obtaining from two estimation of logical learning estimation and knowledge-driven framework codes based numerical estimation (Todorovski and Džeroski, 2006; Ameri et al., 2018; Hembram and Saha, 2018). In this study, CF method was applied for assignment of rank order in each estimated morphometric parameters that are incorporated with ES. Therefore, collective impacts of all parameters on SES are represented by compound factor, estimated from average rank values of all parameters (Altaf et al., 2014). CF is determined as Eq. (20):where R denote parameters rank value, pn denote watershed rank from each parameter.

CF is estimated rank order or erosion priority scale of sub-basin level from fifteen morphometric parameters in thirty-eight sub-basins of Kangsabati basin. According to Patel et al. (2012) and Balasubramanian et al. (2017), high ES has corresponded with maximum compound value of areal and linear aspect, which is given rank 1 and second highest compound value is given rank 2 and others in same assigning manner. While high ES has corresponded with minimum compound value of shape aspect, which is given rank 1 and so on.

Validation of MCDM models

Percentage change and intensity change are uses for the assessment of validation performance and comparison among four MCDM methods results that are predicted SES at sub-basin level (Badri, 2003; Ameri et al., 2018).

Percentage of changes is expressed as Eq. (21):where ΔP denote percentage of changes, N denote alternative numbers, and Nconstant denote alternative numbers of same rank.

Intensity of changes is expressed as Eq. (22):where ΔI denote intensity of changes, r1 denote first method of alternative rank, and r2 denote second method of alternative rank.

SWAT parameter estimation

Several predominant parameters like weather condition, topographic characteristic, soil properties, and land use land cover (LULC) were required for successful set-up and run the SWAT model on a monthly time-scale (1.1.2010–31.12.2015) using the ArcSWAT interface (Neitsch et al., 2011; Markhi et al., 2019). Database of slope, soil and LULC distribution map over the thirty-eight sub-basins including 12,387 HRUs are presented in Fig. 3a, b, c. The unique combinations of LULC and soil type were separately represented by two different HRU. 5% of slope, 5% of soil and 5% of LULC were assigned as user threshold values for the determination of in HRU tool box. Other important data sets like weather databases including rainfall, temperature, and evaporation etc. were extracted from IMD and CFSR World weather gridded database. On the other hand, hydrological data sets of monthly flow discharge or outflow data and sediment load were collected from Mohanpur station, only available monthly flow discharge data set was and collected from Irrigation Office of Paschim Medinipur, but there is no record of any observed data on sediment load. In this context, total sediment load (QT) was computed from field observation database of nine fluvial and sediment hydraulic variables (average flow velocity, water depth, weight of specific water, kinematic viscosity of water, gravity acceleration, sediment diameter, weight of specific sediment, shear velocity, bed shear stress) in Mohanpur station. Total sediment load was calculated in ton unit using fluid weight of Ackers and White (1973) method (Bhattacharya et al., 2019a); furthermore, those results were comparing with the SWAT simulation results. Moreover, sediment load data sets helped in SWAT sensitivity analysis, including calibration and validation test on observed and simulated sediment load. Total Sediment load (QT) has measured from Eq. (23):where Q denote water discharge, X denote sediment concentration.

Fig. 3

SWAT parameter estimation: a slope, b soil and c land use and cover.

Model sensitive and calibration

Sensitive analysis of SWAT parameters must be needed for the determination of significant contribution on model outputs that are measured by ranking method (Holvoet et al., 2005; Markhi et al., 2019). In this context, SUFI-2 (sequential uncertainty fitting approach) is used for large calibration performance with the presence of minimum number of iterations and high calibration capacity (Yang et al., 2008). SWAT-CUP software is adopted for the methodological set-up of SUFI-2 (Abbaspour, 2013; Arnold et al., 2012). Flow discharge and sediment load parameters rank are assigned from SWAT parameters in global sensitive analysis, where statistical indices i.e. co-efficient of determination (R2), Nash-Sutcliffe efficiency (NSE), percent bias (PBIS) are used for SWAT model performance evaluation as goodness of fit of simulated results with observed values (Abbaspour et al., 2004; Rice, 2006). According to Neitsch et al. (2011), SWAT manual provided such default values that are fixed at every parameter range. In this research, ten parameters were considered for monthly flow discharge simulation, and ten parameters were considered for sediment load simulation at monthly time-scale accurately based on literature review on these parameters in similar hydrological characteristics in drainage basin (Bokan, 2015; Mittal et al., 2016; Welde, 2016; Romagnoli et al., 2017; Asl-Rousta and Mousavi, 2019; Sridhar et al., 2018; Himanshu et al., 2019; Markhi et al., 2019). Based on sensitive analysis, best sensitive parameters are selected for the analysis of calibration and validation performance in SWAT model using p-value (Mittal et al., 2016). P-value means fraction values of observed discharge data under 95 Percent Prediction Uncertainty (95PPU) band, which ranges from 0 to 1 (Vaghefi et al., 2013; Mittal et al., 2016). In this study, P-value of <0.05 was considered for selection of sensitive parameters in calibration and validation analysis. To determine the calibration and validation of results, iteration steps were required to find out the optimum values for two different times (Yang et al., 2008). Therefore, the calibration (1.1.2010–31.12.2013) and validation (1.1.2014–31.12.2015) tests between observation and simulation data were applied to determine the model performance. In this context, Mohanpur gauge station was taken for this purpose to make a comparison between observed data and SWAT simulated results for both flow discharge and sediment load estimation.

Model evolution and validation

Nash-Sutcliffe efficiency criterion (NSE) was used for the determination of calibration result in observation data (Nash and Sutcliffe, 1970). This well-known statistical criterion provided the coefficient determination (R2), Nash-Sutcliffe efficiency criterion (NSE), and percent bias (PBIAS) values that were obtained from Eqs. (24), (25), (26), respectively. where Ysin denote simulated value, and Yobs denote observed value; Zobs refer to mean of observed value as n; and Zsin refer to mean of simulated values as n.

Result

Morphometric parameters and its role on ES

Basin morphometric parameters play significant role to determine the earth system processes in surface, incorporating geomorphology, geology and hydrological sets-up (Ifabiyi and Eniolorunda, 2012; Balasubramanian et al., 2017). Moreover, relief and drainage system in the entire basin, have great influence on soil erosion vulnerability, whereas three major morphometric aspects like linear, relief and shape, have control the runoff volume and infiltration capacity (Sharma et al., 1986; Hembram and Saha, 2018). This basin area is extended to 6480 km2 with perimeter of 3306 km2, and texture with 6971 stream number, including 833 km2 stream lengths. 5321 first stream order were identified in high erosion prone sites in the entire catchment area. According to erosion rate, fifteen morphometric parameters were taken from thirty-eight sub-basins to assign the priority rank of ES. Detailed analysis of morphometric parameters in Kangsabati basin can be seen in the supplementary result of ‘morphometric analysis’ section.

Sub-basin prioritization has been done following such relation as morphometric parameters, geological set-up, and others geo-environmental parameters with ES, which also successfully done by several researchers: Biswas et al. (1999), Nooka et al. (2005), Patel et al., 2012, Patel et al., 2013, Ameri et al. (2018), Bhattacharya et al. (2019b). Three different morphometric aspects i.e. linear, shape and relief, including fifteen parameters, have significant relation with ES. According to Biswas et al. (1999) and Nooka et al. (2005), ES is positively correlated with all linear and relief morphometric parameters, whereas Patel et al. (2012) and Patel and Dholakia (2010) denoted that all shape parameters are inversely correlated with ES. In order of erosion priority rank, this study considered maximum value of mean bifurcation ratio and relief ratio that were assigned as highest rank and minimum value of those parameters assigned as lowest rank.

Weight assigning in each morphometric parameters using AHP method

In this research, weight assignment in each criterion was computed from final matrix table that was prepared from fifteen morphometric parameters using AHP method. According to Saaty (1978), weight assignment becomes accepted; if incompatibility rate in a final matrix is lower than 0.1. In this study, consistency index (CI) and randomness index (RI) in comparison matrix table is 0.08 and 1.59, respectively; however, incompatibility rate of final matrix become 0.054; therefore, weight assignment of effective morphometric parameters on ES are within acceptable range (Table 2 ). In addition, weight assignment of fifteen morphometric parameter in matrix table showed that mean bifurcation ratio (0.191), mean slope (0.151), infiltration number (0.129) and stream frequency (0.107) have high significant effects on ES. In contrast, minimum value of shape factor (0.010), form factor (0.014), and length of overland flow (0.015) have least effects on ES (Fig. 4 ). Moreover, drainage density (0.087), relief ratio (0.069), ruggedness number (0.059), compact coefficient (0.048), constant channel maintenance (0.043), circular ratio (0.03), drainage texture (0.029), and basin relief (0.02) have moderate effects on ES, respectively. Therefore, it is a point that all morphometric parameters have not same significant impacts on ES.

Table 2

Computation of composition matrix table (C) following normalized criteria.

	Rmb	S	If	Fs	D	Rr	Rn	C	Cc	Rc	T	R	Lo	Ff	Bs		Criteria weights	%	Remarks
Rmb	1	2	2	3	3	4	4	5	5	6	6	7	8	9	9	1	0.191	19.10%	Consistency (accepted)5%
S	1/2	1	2	2	3	3	4	4	2	6	5	7	8	8	9	2	0.151	15.10%
If	1/2	1/2	1	2	2	3	3	4	4	5	5	6	7	7	8	3	0.129	12.90%	Lambda16.21949815
Fs	1/3	1/2	1/2	1	2	2	3	3	4	5	4	6	6	7	8	4	0.107	10.70%
D	1/3	1/3	1/2	1/2	1	2	2	3	3	4	4	5	6	6	7	5	0.087	8.70%	CI0.087107011
Rr	1/4	1/3	1/3	1/2	1/2	1	2	2	2	3	3	5	5	6	7	6	0.069	6.90%
Rn	1/4	1/4	1/3	1/3	1/2	1/2	1	2	2	3	3	4	5	5	6	7	0.059	5.90%
C	1/5	1/4	1/4	1/3	1/3	1/2	1/2	1	2	2	2	4	4	5	6	8	0.048	4.80%	Randomness Index, RI1.59
Cc	1/5	1/2	1/4	1/4	1/3	1/2	1/2	1/2	1	2	2	3	4	4	5	9	0.043	4.30%
Rc	1/6	1/6	1/5	1/5	1/4	1/3	1/3	1/2	1/2	1	1	3	3	3	5	10	0.03	3.00%
T	1/6	1/5	1/5	1/4	1/4	1/3	1/3	1/2	1/2	1	1	2	3	3	4	11	0.029	2.90%	Consistency ratioCI/RI=0.054784284
R	1/7	1/7	1/6	1/6	1/5	1/5	1/4	1/4	1/3	1/3	1/2	1	2	2	4	12	0.02	2.00%
Lo	1/8	1/8	1/7	1/6	1/6	1/5	1/5	1/4	1/4	1/3	1/3	1/2	1	1	3	13	0.015	1.50%
Ff	1/9	1/8	1/7	1/7	1/6	1/6	1/5	1/5	1/4	1/3	1/3	1/2	1	1	2	14	0.014	1.40%
Bs	1/9	1/9	1/8	1/8	1/7	1/7	1/6	1/6	1/5	1/5	1/4	1/4	1/3	1/2	1	15	0.01	1.00%

Rbm-mean bifurcation ratio, S-mean slope, If-infiltration number, Fs-stream frequency, D-drainage density, Rr-relief ratio, Rn-ruggedness number, C-constant of channel maintenance, Cc-compactness coefficient, Rc-circularity ratio, T-drainage texture, R-basin Relief, Lo-length of overland flow, Ff-form factor, Bs-shape factor.

Fig. 4

The weight of each of the parameters with Analytical Hierarchy Process (AHP).

Erosion priority at sub-basin level using MCDM model

After the weight assignment, all morphometric parameters were calculated for the preparation of decision matrix and data normalization that are essential for priority rank of VIKOR, TOPSIS, SAW, and CF models. In this study, fifteen morphometric parameters were used for the computation of priority rank in each sub-basin, where all data sets were normalized to prepare the rank value in four different MCDM methods. Linear vector method was used for the computation of data normalization in TOPSIS model (Eqs. (8), (9), (10), (11), (12), (13), (14)), whereas data normalization of VIKOR and SAW models were computed following Eqs. (1)–(7), (15), (16), (17), (18), (19) (Ameri et al., 2018). In VIKOR model, best and worst values in each criterion helps to normalization of regret levels; however, based on perspective of regret theory, only best values in each criterion play dominant role to determine the regret levels, and worst values in each criterion play dominant role to determine the effective role of normalized S and R values on regret levels (Huang et al., 2009). In this study, best and worst values are calculated in Table 3 using Eqs. (1), (2), whereas utility index, regret index both were determined the priority rank as ascending order in the thirty-eight sub-basins using Eqs. (3), (4), (5) (Table 3). In TOPSIS model, positive value, negative value and Euclidean distance between them were calculated using Eqs. (8), (9), (10), (11), (12), (13) and is presented in Table 3. While closeness coefficient options of ideal solution of TOPSIS model was done using Eq. (14). On the other hand, final weights of normalized matrix row in each sub-basin were calculated for SAW model according to Eq. (19). In contrast, CF values in each sub-basin were assigned through the calculation of rank value from fifteen parameters; finally summation of all values was divided by such parameter numbers according to Eq. (20). Sub-basin priority has been assigned as descending manner, that means, lowest CF given as first rank for high ES, and highest CF given as last rank for low ES in Table 4 (Patel et al., 2012; Altaf et al., 2014; Bhattacharya et al., 2019b).

Table 3

Computation of the best and worst values, PIS, and NIS indicators in VIKOR and TOPSIS models, respectively.

Morphometric parameters		Rbm	D	S	Fs	T	Lo	If	C	Ff	Bs	Cc	Rc	R	Rn	Rr
TOPSIS	V+	0.053817	0.026172	0.047162	0.044262	0.022503	0.004512	0.066244	0.017744	0.001703	0.001326	0.004673	0.001298	0.006958	0.02245	0.01978
	V-	0	0.004039	0.004092	0.001081	0.000132	0.000696	0.000357	0.002738	0.002694	0.002097	0.011747	0.0082	0.001229	0.000996	0.002871
VIKOR	Xj+	0.191	0.087	0.151	0.107	0.029	0.015	0.129	0.048	0	0	0	0	0.02	0.059	0.069
	Xj-	0	0	0	0	0	0	0	0	0.014	0.01	0.043	0.03	0	0	0

Table 4

Morphometric based assigned ranking and compound value in Kangsabati basin.

SB	Rbm	D	S	Fs	T	Lo	If	C	Ff	Bs	Cc	Rc	R	Rn	Rr	Compound parameter	Final priority
1	2	9	34	18	11	9	17	30	4	35	37	2	28	5	16	17.13	13
2	3	8	28	22	16	8	16	31	9	30	28	11	27	2	2	16.07	10
3	6	16	38	20	8	16	21	23	3	36	32	7	35	15	37	20.87	22
4	12	20	36	14	18	20	19	19	11	28	36	3	37	18	35	21.73	27
5	5	17	37	13	12	17	13	22	5	34	38	1	31	12	36	19.53	19
6	9	10	29	15	9	10	12	29	7	32	30	9	23	1	5	15.33	7
7	7	12	32	16	6	12	15	27	6	33	29	10	26	6	14	16.73	12
8	14	11	26	11	19	11	9	28	22	17	26	13	34	11	15	17.80	14
9	17	19	30	19	21	19	22	20	15	24	33	6	36	17	26	21.60	26
10	19	13	27	21	22	13	20	26	27	12	27	12	38	20	24	21.40	24
11	11	14	35	17	17	14	18	25	10	29	35	4	33	13	30	20.33	21
12	8	15	33	12	13	15	11	24	8	31	34	5	32	14	31	19.07	17
13	4	6	31	7	2	6	7	33	2	37	31	8	16	4	25	14.60	6
14	13	3	18	5	3	3	4	36	13	26	4	35	14	7	7	12.73	3
15	16	21	24	8	7	21	8	18	14	25	6	33	18	10	8	15.80	8
16	10	1	20	2	5	1	2	38	17	22	25	14	20	8	17	13.47	4
17	22	7	19	6	14	7	6	32	34	5	9	30	29	16	1	15.80	8
18	15	5	15	3	4	5	3	34	24	15	5	34	15	9	4	12.67	2
19	23	37	21	27	29	37	29	2	37	2	22	17	30	36	6	23.67	31
20	26	28	14	25	26	28	26	11	29	10	13	26	19	23	18	21.47	25
21	30	34	23	30	30	34	31	5	20	19	14	25	25	30	22	24.80	35
22	27	26	22	29	27	26	28	13	19	20	7	32	24	28	29	23.80	32
23	1	2	25	1	1	2	1	37	1	38	21	18	1	3	34	12.40	1
24	36	33	17	33	34	33	33	6	28	11	19	20	22	31	19	25.00	36
25	27	24	4	24	25	24	24	15	30	9	8	31	7	25	13	19.33	18
26	18	18	5	10	15	18	10	21	21	18	16	23	3	19	27	16.13	11
27	37	31	6	32	36	31	32	8	36	3	12	27	17	32	20	24.00	34
28	20	23	3	9	20	23	14	16	35	4	15	24	11	24	28	17.93	15
29	23	27	12	26	24	27	25	12	25	14	1	38	9	27	23	20.87	22
30	25	22	7	23	23	22	23	17	31	8	17	22	11	22	12	19.00	16
31	38	36	10	38	38	36	38	3	33	6	11	28	21	35	11	25.47	37
32	31	25	2	28	28	25	27	14	32	7	23	16	8	26	9	20.07	20
33	32	32	13	34	33	32	34	7	16	23	24	15	4	33	21	23.53	29
34	34	38	8	36	35	38	37	1	26	13	3	36	6	37	10	23.87	33
35	32	29	11	31	32	29	30	10	23	16	10	29	5	34	32	23.53	29
36	29	4	1	4	10	4	5	35	38	1	2	37	10	21	3	13.60	5
37	34	35	16	37	37	35	36	4	18	21	18	21	13	38	33	26.40	38
38	20	30	9	35	31	30	35	9	12	27	20	19	2	29	38	23.07	28

In VIKOR model, sub-basin priority has incorporated with ES, which showed that least compound values of SB 23 (0.144), 16 (0.155), 18 (0.156) are given highest rank as 1 to 3 for their high ES, while more compound values of SB 31 (0.999), 24 (0.811), 37 (0.808) are given lowest rank as 38 to 36 for their low ES (Fig. 5a). In TOPSIS model, there is positive relation between compound value and assigned rank value based on erosion sensitivity. In particular, maximum values of SB 23 (0.734), 16 (0.689), 18 (0.58) are assigned highest rank as 1 to 3 for most sensitivity to erosion, while minimum values of SB 21 (0.205), 24 (0.213), 37 (0.526) are assigned lowest rank as 38 to 36 for their least sensitivity to erosion (Fig. 5b). Based on positive relation between final scores and ES, in SAW model, maximum final scores of SB 23 (0.767), 16 (0.682), 18 (0.6127) are given highest rank as 1 to 3, while SB 31 (0.276), 37 (0.278), 21 (0.287) are given lowest rank as 38 to 36, respectively. It is a point that sub-basins with more complex score have high sensitivity to erosion, and sub-basins with least final score have low sensitivity to erosion (Fig. 5c). On the other hand, according to CF model, least compound value of SB 23 (12.40), 18 (12.67), 14 (12.73) are assigned rank as 1 to 3 for their high ES, in contrast, more compound values of SB 37 (26.40), 31 (25.47), 23 (25) are given rank as 38 to 36 for their least sensitivity to erosion, respectively (Table 4). Based on ES and loss of natural resources, all priority rank of MCDM models (except CF) are well categorized into four classes i.e. very high (0.75–1), high (0.5–0.75), moderate (0.25–0.5) and low (<0.25) in thirty-eight sub-basins (Ameri et al., 2018). In this study, four categories (very high, high, moderate and low) in VIKOR model, three categories (high, moderate and low) in TOPSIS model, and three categories (very high, high and moderate) in SAW model are identified (Fig. 5a, b, c). On the other hand, based on compound values from selected morphometric parameters, CF model has classified into four categories i.e. very high (>100), high (50–75), medium (25–50) and low (>25). In this study, only two CF categories of moderate and low are identified in the thirty-eight sub-basins (Fig. 5d). According to rank values, MCDM priority classes showed that moderate class (55%) is most dominant than low (19%), high (18%) and very high classes (8%) in VIKOR model; high class (71%) is most dominant than moderate (13%) and low classes (16%) in TOPSIS; and maximum moderate class (79%) has dominated in SAW model than high (18%) and very high classes (3%), respectively. In contrast, low class (92%) has reaches dominant position than medium class (8%) in CF model. Therefore, ES distribution is not same in four different MCDM methods in the entire Kangsabati basin. However, MCDM method denoted that moderate and low erosion classes are predominant categories throughout the basin.

Fig. 5

Classification of sub-watershed to erodibility using MCDM models: a VIKOR, b TOPSIS, c SAW and d CF.

MCDM models validation

Two different evaluation methods of percentage change and intensity change were used for the determination of MCDM models efficiency. Percentage change result showed best efficiency of CF method with high accuracy value (68.42) than VIKOR (57.89), TOPSIS (42.11), and SAW methods (31.28) (Table 5 ), whereas intensity change denoted that VIKOR method is more efficient with high accuracy (3.165) than TOPSIS (3.097), CF (3.087) and SAW methods (3.085), respectively (Table 6 ). According to percentage change and intensity change of MCDM validation test, it is the point that CF and VIKOR have the best efficiencies with high accuracy values than other methods.

Table 5

Percentages of change in models.

	CF	VIKOR	TOPSIS	SAW	SUM
CF	0	94.736	89.473	84.21	68.421
VIKOR	94.736	0	84.21	78.947	57.894
TOPSIS	89.473	84.21	0	68.421	42.105
SAW	84.21	78.947	68.421	0	31.578

Table 6

Intensity of change in models.

	CF	VIKOR	TOPSIS	SAW	SUM
CF	0	1.044	1.025	1.018	3.087
VIKOR	1.053	0	1.06	1.052	3.165
TOPSIS	1.033	1.061	0	1.003	3.097
SAW	1.026	1.056	1.003	0	3.085

SWAT model

Sensitivity, calibration and validation analysis on flow discharge and sediment load

In general, SWAT model performance is evaluated by three different analysis of sensitivity, calibration, and validation from assigned parameters (Welde, 2016). In this study, sensitive analysis was done in a Latin hypercube sampling at 12 intervals, where twenty parameters were considered for flow discharge and sediment load estimation. 1000 iterations were used for the determination of following parameters, where 500 iterations for flow parameters (10 × 50 per iteration), and 500 iterations for sediment parameters (10 × 50 per iteration) were required with the help of SWAT-CUP software. Table 7 denotes maximum, minimum and fitted values of ten sensitive parameters for flow discharge simulation. Table 8 denotes maximum, minimum and fitted values of ten parameters for sediment load simulation. Sensitivity analysis reveals that CN2.mgt (number of SCS runoff curve in moisture condition II) is more sensitive followed by GWQMN.gw (occurring of return flow for threshold depth of water) in flow discharge calibration and validation. While R_SPEXP.bsn (exponent parameter for calculating sediment recent rained in channel sediment) is more sensitive followed by R__SLSUBBSN.hru (average slope length) in sediment load calibration, but HRU_SLP.hru (average slope steepness) is more sensitive followed by R_SPEXP.bsn in sediment load validation. After the successful run of calibration procedure on sensible parameters, Table 7, Table 8 both are indicated that the SWAT model has sensibly explained of hydrological processes in Kangsabati basin. During calibration period (January 2010-Decembar 2013), observed versus model simulated coefficient of regression (R2) values of 0.86 and 0.87 for flow discharge and sediment load (Table 9 ; Fig. 6a, b, c, d), respectively; Nash-Sutcliffe efficiency criterion (NSE) values of 0.75 and 0.69 for flow discharge and sediment load, respectively; denotes well acceptable range (Moriasi et al., 2007). During validation period (January 2014-Decembar 2015), observed versus model simulated R2 values of 0.79 and 0.79 for flow discharge and sediment load, respectively; NSE values of 0.55 and 0.76 for flow discharge and sediment load, respectively; showed recommendable range (Moriasi et al., 2007). In addition, PBIAS values of −19.9 and 34.6 for flow discharge and sediment load calibrations, respectively; −29 and 10.2 for flow discharge and sediment load validations, respectively; demonstrated that average performance level of SWAT model has underestimated of flow discharge by 19.9% and 29% during calibration and validation period, respectively, however; overestimated of sediment load by 34.6% and 10.2% during calibration and validation period (Table 9), respectively. In summary results of calibration test, mean and standard deviation values of observed and simulated flow discharge are 25.45 and 35.63, and mean and standard deviation values of observed and simulated sediment load are 11,270.16 and 14,224.11, respectively. In summary results of validation test, mean and standard deviation values of observed and simulated flow discharge are 20.32 and 30.77, and mean and standard deviation values of observed and simulated sediment load are 14,655.22 and 18,284.36, respectively. Fig. 6 indicated that the observation value of monthly flow discharge becomes more or less similar with simulated monthly flow discharge, whereas observed values of sediment load (ton/month) has greater than the simulated value at Mohanpur station, during the calibration and validation phases. Therefore, based on all statistical indices of calibration and validation, simulated results are considered as goodness of fit with observed data.

Table 7

Sensitive parameters and its fittest values of flow discharge after calibration using SUFI 2.

Parameter name	Definition	Reference	Fitted value	Minimum value	Maximum value
R_CH_N2.rte	Manning's ‘n’ value for the main channel	Mittal et al. (2016); Himanshu et al. (2019); Markhi et al. (2019)	−0.05	−0.068	0.11
R_TLAP.sub	Temperature lapse rate	Asl-Rousta and Mousavi (2019)	−6.46	−7.22	−3.65
R_GWHT.gw	Initial groundwater height (m)	Markhi et al. (2019)	5.01	−1.89	5.72
R_PLAPS.sub	Precipitation lapse rate	Asl-Rousta and Mousavi (2019)	−182.56	−199.1	186.83
R_TIMP.bsn	Snow pack temperature lag factor	Asl-Rousta and Mousavi (2019)	0.78	0.59	0.98
R_SNOCOVMX.bsn	Maximum snow water content that corresponds to 100% snow cover	Asl-Rousta and Mousavi (2019)	140.9	107.65	230.16
R_CN2.mgt	Number of SCS runoff curve in moisture condition II	Romagnoli et al. (2017); Mittal et al. (2016); Himanshu et al. (2019)	0.103	0.08	0.13
V_ALPHA_BF.gw	Ground water recession curve for base flow alpha factor	Mittal et al. (2016); Sridhar et al. (2018); Himanshu et al. (2019)	1.26	1.08	1.41
V_GW_DE LAY.gw	Groundwater delay time from bottom root zone to shallow aquifer	Romagnoli et al. (2017); Mittal et al. (2016); Himanshu et al. (2019)	42.89	−56.13	121.61
V_GWQMN.gw	Occurring of return flow for threshold depth of water	Romagnoli et al. (2017); Mittal et al. (2016); Sridhar et al. (2018); Himanshu et al. (2019)	1.65	1.39	2.18

Table 8

Sensitive parameters and its fittest values of sediment discharge after calibration using SUFI 2.

Parameter name	Definition	Reference	Fitted value	Minimum value	Maximum value
R_SPCON.bsn	Calculation of linear factor for maximum amount of channel sediment re-entrained	Bokan (2015); Welde (2016); Markhi et al. (2019)	0.0033	−0.0037	0.0059
R__USLE_K (..).sol	USLE soil erodibility factor	Bokan (2015); Himanshu et al. (2019); Welde (2016)	0.32	0.049	0.35
R_SPEXP.bsn	Exponent parameter for calculating sediment recent rained in channel sediment	Bokan (2015); Welde (2016)	1.19	1.08	1.28
R_USLE_C {..}.plant.dat	USLE cropping factor for strip cropped fields	Bokan (2015); Welde (2016)	0.27	0.04	0.33
R__BIOMIX.mgt	Biological mixing efficiency	Bokan (2015); Himanshu et al. (2019)	0.85	0.60	1.31
R__RSDIN.hru	Initial residue cover [kg/ha]	Bokan (2015)	4518.76	381.1	10,161
R__CH_COV1.rte	Channel erodibility factor	Bokan (2015)	0.35	0.25	0.64
R__CH_COV2.rte	Channel cover factor	Bokan (2015); Welde (2016); Himanshu et al. (2019)	1.08	0.55	1.2
R__SLSUBBSN.hru	Average slope length	Bokan (2015); Himanshu et al. (2019)	−44.18	−48.11	54.19
R_HRU_SLP.hru	Average slope steepness	Bokan (2015); Welde (2016)	0.005	−0.59	0.31

Table 9

Model performance for simulation of flow discharge and sediment concentration at Mohanpur station (2010–2015).

Parameters	Flow discharge		Sediment load
	Calibration	Validation	Calibration	Validation
R2	0.86	0.79	0.87	0.79
NSE	0.75	0.55	0.69	0.76
PBIAS	−19.9	−29	34.6	10.2
Mean	25.45	20.32	11,270.16	14,655.22
Standard deviation	35.63	30.77	14,224.11	18,284.36
Time step	Monthly	Monthly	Monthly	Monthly
Iterations	500		500

Fig. 6

Monthly calibration and validation at Mohanpur station: a calibration of flow discharge (m3s−1), b validation of flow discharge (m3s−1), c calibration of sediment load (ton/month) and d validation of sediment load (ton/month).

Sediment yield analysis and its sub-basin prioritization

Suspended sediment load and its movement in stream has entirely depended on several hydraulic variables like stream discharge, watershed slope, including flow and sediment regime characteristics (Sridhar et al., 2018; Bhattacharya et al., 2020a). Huge sediment load helps to accumulation of sediment yield at sub-basin level where SWAT model is applied to estimate the amount of yield following some considerable hydrological parameters (Welde, 2016). In Kangsabati basin, monthly simulated sediment yield of SWAT model is classified into five categories i.e. very low (<0.56 m ton/ha), low (0.57–0.83 m ton/ha), middle (0.84–1.08 m ton/ha), high (1.09–1.38 m ton /ha) and very high (1.39–1.88 m ton/ha), respectively. Fig. 7 showed that most of the basin area come under low to medium sediment yield priority class, whereas high and very high sediment yield classes are mainly concentrated at outlets or confluence points. Generally, absences of flow discharge, low drainage density and lower catchment area helps to generate low sediment yield, while presence of confluence points, high drainage density and large catchment area helps to accumulate of high sediment yield (Bhattacharya et al., 2020a). In terms of ES at sub-basin level, monthly sediment yield priority is classified into four classes i.e. low (<0.69 m ton/ha), moderate (0.70–0.93 m ton/ha), high (0.94–1.17 m ton/ha) and very high (>1.17 m ton/ha). According to assigned rank of simulated sediment yield, very high erosion priority class mainly concentrated in SB 13, 16, 18, 19, 23, while high erosion priority class concentrated in SB 2, 6, 14, 17, 20, 22, 25, 26, 28 and 31, respectively (Fig. 8 ). In addition, SB 1, 3, 4, 5, 7, 8, 9, 10, 12, 15, 24, 29, 30, 32 and 34 having monthly sediment yield of 0.70–0.93 m ton/ha comes under moderate erosion priority class, and SB 1, 3, 4, 5,7, 8, 9, 10, 12, 15, 24, 29, 30, 32 and 34 having monthly sediment yield of 0–0.69 m ton/ha comes under low erosion priority class. Therefore, it is pointed that amount of sediment yield are fully dependent with ES level throughout the basin.

Fig. 7

Sediment Yield (SY) distribution zone using SWAT model.

Fig. 8

Classification of sub-watershed to erodibility using SWAT model.

Comparison between MCDM model and SWAT model

MCDM results showed that most of the sub-basins fall under moderate ES class, predicted by VIKOR (55%), TOPSIS (71%) and SAW methods (79%), in contrast, CF method predicted that nearly 92% area are considered as low ES class in Kangsabati basin. Validation result of MCDM methods demonstrated that prediction of VIKOR and CF methods are within acceptable range; however, there is no single model that can be acceptable for ES. In this context, based on assigning the priority rank from simulated sediment yield in thirty-eight sub-basins, SWAT model is used to select the effective MCDM methods on ES. In this study, SWAT model predicted that moderate to low ES classes are predominant (62%) than high and very high erosion susceptible classes (38%) throughout the basin where annual stream discharge and annual average sediment load are 19.23m3/s and 12,312 ton, during 2010–2015. After successful validation of MCDM methods, results indicated that VIKOR and CF methods are in acceptable range; in contrast, calibration and validation results of SWAT model are good fit with observed data. In order to ES priority rank, VIKOR model nearly matches with SWAT model, according to both model predictions, moderate to low erosion classes are considered as more dominant classes throughout the basin. These susceptible classes are major sediment source sites where huge sediment load supply helps to vast accumulation of sediment yield. Therefore, VIKOR method is an effective MCDM model, which helps to identify the erosion-prone sites, sediment load supply, and sediment yield accumulation throughout the basin.

Discussion

Morphometric parameters are not considered as only key factors for the determination of ES, but other dominant geo-environmental parameters such as LULC, soil erodibility, slope length-steepness, lithology, geomorphic set-up, etc. are also important (Chauhan et al.,2016; Bhattacharya et al., 2019b). In comparison between MCDM methods and SWAT model, effective MCDM methods are helpful to identify the significant role of morphometric parameters on ES; furthermore, linear and aerial aspects are positively correlated with ES (Nooka et al., 2005; Bhattacharya et al., 2019b), but shape aspects are inversely correlated with ES (Patel et al., 2012, Patel et al., 2013). In contrary, when preparing sediment yield distribution map from seven considerable factors: LULC, slope, soil, ground water depth, run off volume, precipitation lapse rate and channel hydro-morphogenetic properties (Markhi et al., 2019; Himanshu et al., 2019), SWAT model also helps to assess ES in each sub basin. Therefore, in this study, LULC and morphometric parameters are considered for ES assessment following the validation performance of MCDM and SWAT models, respectively. In terms of response to erosion priority, land covers like laterite with barren land, wasteland, and settlement are positively related to ES, while agricultural land, dense forest, pasture land are inversely related to ES (Altaf et al., 2014).

Based on MCDM validation performances; TOPSIS, SAW and CF does not give satisfactory results due to disadvantages of unavailable data sets for all decision-making problems and assumption based relative weight assignment in each variables (Khosravi et al., 2019). In contrast, VIKOR method has given better results than other MCDM methods, which are validated by SWAT model for its advantages: hierarchical formulating issue, pair wise comparison using expert quantitative and qualitative knowledge, and assessment of compatibility and incompatibility decision (Saaty, 1980; Arabameri et al., 2019).

SWAT and VIKOR helps to understand the role of LULC and morphometric properties on ES in Kangsabati basin. In spite of dense vegetation cover and pasture land, SB 13, 16, 18, 23 are more susceptible to erosion due to maximum values of linear morphometric parameters (mean bifurcation ratio, drainage density, drainage texture, stream frequency) and relief morphometric parameters (slope, ruggedness index, relief ratio, basin relief) as well as minimum values of basin form parameters (form factor, shape factor, compact coefficient and circular ratio) (Fig. 9 ). In contrast, despite the presence of barren land, double crop practice and dense settlement, SB 38, 36, 35 faces are of low ES category due to minimum values of linear and relief parameters as well as maximum values of basin parameters. Moreover, moderate erosion susceptible class in rest thirty one sub-basins has dominant susceptibility due to maximum coverage of barren land, pasture land, and generic agricultural practices, on the other hand, morphometric parameters allowed same susceptible class in those sub-basins. As a result, maximum sediment yield deposition is found in SB 13, 16, 18, 23 with an average value of 1.27 m ton/ha, and minimum sediment yield deposition is found in SB 38, 36, 35 with an average value of 0.5 m ton/ha, respectively. Furthermore, rest of the sub-basins reaches in moderate sediment yield class with an average value of 0.8073 m ton/ha. Thus, based on the above discussion, it can be said that among all the available MCDM methods for ES, VIKOR model is more pragmatic for ES in respect to simulated sediment yield rank as predicted by SWAT. In order to ES and sediment yield deposition at sub-basin level, in this study demonstrated that VIKOR is fittest MCDM model for suitable selection of best and worst values from morphometric parameters using normalize decision matrix of AHP, maximum group utility index, minimum regret index of opposite group value, as well as assessment of compatibility and incompatibility decision among effective morphometric parameters. Therefore, VIKOR model is useful MCDM method to prepare rational sub-basin prioritization from linear normalize ranking of all morphometric parameters. Moreover, VIKOR method helps to identify the five critical sub-basins that are most sensible to erosion due to presence of responsible morphometric parameters in line with the results of Bhattacharya et al. (2020a). On the other hand, SWAT parameters helps to determine the significant role of morphometric parameters on ES in response to sediment yield deposition throughout the basin. Therefore, it can be stated that morphometric parameters are considered as crucial contributing parameters for ES generation in this plateau fringe basin, followed by LULC patterns, climate, and soil characteristic. These agreements are validated with the findings of Biswas et al. (1999), Hembram and Saha (2018), Sridhar et al. (2018). Based on literature review in introduction, previous researchers have presented sub-basin prioritization in order to ES using MCDM methods; however, there is no hydrological model to measure the significant role of hydrological parameters like runoff volume, sediment load, sediment concentration, etc. on sub-basin prioritization. In this context, there is research gap in previous studies of sub-basin prioritization. Present study tried to address this research gap using a comparison between MCDM and SWAT models. In term of ES, present study humbly argues that effective morphometric and hydrological parameters are prerequisite for sub-basin prioritization in any region in the world as assigned by VIKOR and SWAT models. Moreover, both models reveal that all morphometric and hydrological parameters are not equally significant in every sub-basin as they have own characteristics.

Fig. 9

Priority base soil erosion susceptibility location in SB 13, 2, 23, 36 of Kangsabati basin.

The findings in this research might help to identify critical sub-basins where ES has been found most severe with the presence of responsible factors, thus VIKOR and SWAT could provide important tools for planner or policy maker to take rational strategies for soil and water conservation in watershed management.

Conclusion

The present study demonstrated that ES is a sensitive criterion to determine the sub-basin prioritization using a comparison among MCDM methods (VIKOR, TOPSIS, SAW, and CF). On the other hand, effective MCDM methods are selected by applying of SWAT model. Extraction of morphometric parameters from SRTM DEM is required to assign the priority rank on ES for MCDM models, on the other hand, the spatial distribution of LULC, slope, and soil, including twenty sensitive parameters, are required for SWAT model estimation that is predicted the simulated sediment yield as well as ES in thirty-eight sub-basins of Kangsabati basin. To evaluate the MCDM model validation performance on ES, percentage and intensity change are tested, and performance of SWAT model is determined using calibration and validation test on flow discharge and sediment load. According to priority class of ES predicted by MCDM, ES is divided into four different classes including very high, high, medium and low class. In MCDM methods, VIKOR method predicted that moderate and lower erosion classes are dominated than high and very high classes, whereas TOPSIS method predicted that high and very high classes are most dominated than medium and lower classes, respectively. In addition, SAW model predicted that moderate class is predominated class than high and very high classes, but lower class is most dominated in CF model. Although, lower class is completely absent in SAW method, while very high and high classes are almost absent in CF method. In contrary, SWAT model predicted that moderate and low erosion classes are much more dominant than high and very high classes using simulated sediment yield distribution map where superior sediment yield in each sub basins are corresponded with high erosion priority class. Finally, MCDM validation results demonstrated that VIKOR and CF model have more acceptable range, while calibration and validation results of SWAT model has acceptable range that means simulated results are well fitted with observed data. Among all the effective MCDM methods, prediction of VIKOR model matches with SWAT speculation. It is seen that moderate erosion class represents dominant position throughout the basin when we apply it on fifteen morphometric parameters, and SWAT model applying on twenty sensitive parameters of flow discharge and sediment load. This susceptible class leads to the supplied average value of sediment load as well as accumulation of sediment yield. Therefore, VIKOR model is an effective MCDM model, which helps to identify erosion-prone sites and recommends on mitigation measure to protect erosion especially in critical sub-basins.

CRediT authorship contribution statement

Raj Kumar Bhattacharya: Conceptualization, Methodology, Software, Data curation, Writing - original draft. Nilanjana Das Chatterjee: Supervision, Writing - review & editing. Kousik Das: Software, Validation, Visualization, Investigation.

Declaration of competing interest

We are no conflict of interest in this work.