Geospatial Appraisal of Gully Erosion Vulnerability in the Rarh Bengal of India through Analyzing the Multiple Discriminating Factors

Geospatial Appraisal of Gully Erosion Vulnerability in the Rarh Bengal of India through Analyzing the Multiple Discriminating Factors

*Biraj Kanti Mondal¹, Tanmoy Basu², Malini Biswas³, Ming-An Lee⁴

¹Department of Geography, Netaji Subhas Open University, Kolkata 700064, West Bengal, India

²Department of Geography, Katwa College, Katwa, Purba Barddhaman, West Bengal 713130, India

³Department of Geography, Netaji Subhas Open University, Kolkata 700064, West Bengal, India

⁴Department of Environmental Biology Fisheries Science, National Taiwan Ocean University, 2 Pei-Ning Rd., Keelung 20224, Taiwan

*Corresponding Author

DOI: https://doi.org/10.51244/IJRSI.2025.120500179

Received: 03 June 2025; Accepted: 06 June 2025; Published: 23 June 2025

ABSTRACT

This study aims to determine the degree and vulnerability of gully erosion and related soil erosion in the Birbhum district of West Bengal, India. Composite maps of gully erosion susceptibility were created using sophisticated geospatial methods and remotely sensed satellite data. The normalized indicator values were obtained using factor analysis of the 2001 data. Gully erosion during the monsoon season was the main cause of the considerable loss of lateritic soil cover in the Rampurhat-I and Bolpur-Santiniketan blocks, according to the analysis, which also showed a strong relationship between soil erosivity and other influencing factors. Gully erosion impacted 69.81 square kilometers (20.59%) of Bolpur-Santiniketan, primarily in the southeast and northwest, and 68.97 square kilometers (23.45%) of Rampurhat-I, primarily in the southwest and northwest. The topographic wetness index showed the most variability, accounting for 77% and 74% of the erosion variance, respectively, with seven major components. The Rampurhat-I and Bolpur-Santiniketan susceptibility indices ranged from 0.833 to -0.772 and 0.756 to -1.060, respectively. Significant agricultural land loss (from 165.54 to 128.44 square kilometers) in Rampurhat-I and the existence of 26.98 square km of badlands in Bolpur-Santiniketan were also noted by the study, even though places like Ballavpur still had deep forest cover. Land use, land cover, and landholding sizes have all been greatly impacted by the rising rate of soil erosion, particularly in places prone to gullies. Finding hotspots for gully erosion, charting its severity, and making accurate predictions can help guide initiatives to reduce soil loss and degradation, promoting environmentally friendly farming methods and sustainable land management in the area.

Keywords: Gully erosion, Rarh region, Susceptibility, Geospatial techniques, LULC

INTRODUCTION

Globally, soil erosion is one of the most vulnerable geomorphic hazards in dry and drought-prone areas. The world’s lateritic and arid zones are severely impacted by seasonal soil erosion. The gully formation is a type of rill erosion of the soil surface and subsurface. Salleh and Mousazadeh (2011) define gullies as steep-sided watercourses that are subject to sporadic water flow, characterized by a stepped longitudinal profile and typically an abrupt channel head. A gully is a typical type of landform that is created mainly by runoff and mass movement. Gullies look like a large creek or small stream valley, which are a few meters in depth and width. Gullies are commonly seasonal or non-perennial water flow, usually associated with localized intense rainfall events or snowmelt. Gullies can be formed and accelerated by cultivation practices on hilly tracts, often on gentle gradient land (Bandyopadhyay, 1988). Soil erosion is a natural process over any land use and land cover. About 33% of gullies were found in the croplands of the European Union (Borrelli et al. 2022), which caused a majority of changes in land use land cover. The agents of soil erosion contribute a significant amount of soil. As a gradual erosion process, soil erosion causes a serious loss of topsoil (Telkar, 2018). The problem of soil erosion exists all over the World mainly in tropical arid and semi-arid areas. About 53.3% of the total geographical area of India is affected by soil erosion and land degradation (Brady, 2010). There are three major agents of soil erosion; those are soil erosion by water, soil erosion by wind, and soil erosion by biotic elements. Gully erosion is considered one of the major soil erosion mechanisms by water flow. According to Rahmati et al. (2022), drainage density, plan and profile curvature, and index of topographic wetness as such are the major factors of occurring gully erosion in the Biram region of south-western Iran. Being a more prominent type of soil erosion gullies are formed due to the transportation of topsoil by the wider run-off after heavy rainfall resulting in cavities or grooves (Geyik, 1986). Gullies are the major sources of transported sediments although they occupy a small area of the catchment (Ionita et al. 2015). Aslam (2020) depicted that 18% of the classified areas were under the high susceptibility zone of gully erosion and it was considered a potential erosional hazard in the Chitral of north-east Pakistan. Other studies such as Nwilo et al. (2021, p. 157) considered the soil erosion caused by geomorphic agents of gullies as a serious environmental problem in the Imo River Basin of south-east Nigeria. Gullies dissect a large area of land into small fragments which are called badland topography, unsuitable for cultivation and other land uses. Transport capacity and shear stress are the factors that positively correlate with the daily event of runoff in the gully areas in the Dwarka-Brahmani interfluve region in eastern India (Ghosh et al. 2021). Ghosh and Saha (2015) delineated the erosion susceptibility zones of gully erosion in the Hinglo river basin, in Eastern India, and identified that high-weighted composite scores of severe erosion susceptibility were associated with excessive drainage conditions. According to Jahantigh (2011), about 2.35 million cubic meters/year (m³/year) of agricultural soil is lost by gully erosion in the Rewa district of Madhya Pradesh, India. As per the observations by Debanshi and Pal (2018) to control or prevent gully erosion, check dams and gully head bandhs were constructed in the gully erosion areas of the Mayurakshi river basin in eastern India. Soil erosion due to gully formation has been found in several parts of West Bengal. The districts of the western Rarh region of West Bengal are mainly impacted by gully erosion. Shit et al. (2015) postulated that gully erosion increased in the moderately high slope areas in the monsoon season in Garbeta in Paschim Medinipur. A total of six gully erosion hazard-prone areas were identified by Shit et al. (2014) situated mainly in Bankura and Paschim Medinipur districts. > 60% of sand areas covered and > 37.61% of the upper catchment of Chandrabhaga sub-basin in Birbhum district which was adjacent to gully erosion areas (Pal, 2016). Santra and Mitra (2020, p. 1200) analyzed that 25% of the area was under the high-to-severe soil erosion of rills and gullies in Puruliya, the westernmost district of West Bengal.

Literature Review on Methodologies for Assessing Gully Erosion Intensity and Susceptibility Using Geospatial Techniques

Assessing gully erosion, a major type of soil degradation that affects many parts of the world, is essential for efficient land management and conservation. To evaluate the severity and vulnerability of gully erosion, researchers have created and improved methods over time utilizing geospatial tools. With an emphasis on the application of multiple discriminant factor analysis (MDF), remote sensing, GIS, and machine learning techniques, this literature review offers a summary of the approaches used globally, in India, and particularly in West Bengal. Gully erosion has been widely mapped and evaluated using geospatial methods, especially remote sensing and Geographic Information Systems (GIS). The vulnerability of landscapes to gully development has been measured and predicted using various analytical methods. Multiple Discriminant Analysis (MDA), which creates susceptibility maps by combining variables like terrain, soil type, land use, rainfall, and vegetation, is one of the main techniques utilized globally. Poesen et al. (2003), for example, emphasized the application of MDA in the African Sahel, where it proved successful in locating erosion-prone regions by analyzing soil characteristics and climate data. The Analytic Hierarchy Process (AHP) in conjunction with GIS has gained widespread acceptance in the Mediterranean regions. By considering variables including slope, lithology, vegetation cover, and rainfall intensity, Valentin et al. (2005) showed how AHP combined with GIS could rank gully susceptibility. In places like Australia, where soil erodibility and land use changes are important factors, the use of spatial models that incorporate digital elevation models (DEMs), remote sensing data, and climatic information has become more popular for predicting areas at risk of gully erosion (Wilkinson et al., 2005). Furthermore, in areas with complicated topography, like South Africa, machine learning methods like Support Vector Machines (SVM) and Artificial Neural Networks (ANN) have been used to forecast gully susceptibility (Le Roux et al., 2007). These algorithms categorize regions according to their risk of gully erosion using vast datasets, such as topography and remote sensing photography.

Gully erosion has been thoroughly researched in India in some areas where erosion dynamics are greatly influenced by the interplay of monsoonal rainfall, soil properties, and land use changes. In this research, remote sensing and GIS have become essential tools, and various approaches have been created and implemented in Indian contexts. For instance, the Normalized Difference Vegetation Index (NDVI) and the Normalized Difference Built-up Index (NDBI) have been used to track urbanization and vegetation cover, two factors that affect gully erosion. In regions like the Deccan Plateau, researchers like Kumar et al. (2012) have demonstrated that changes in land cover, particularly the shift from agricultural to arid land, make the region more vulnerable to gully erosion. Runoff and its impact on gully formation have been estimated in Indian research using the SCS-CN method (Soil Conservation Service Curve Number method) in conjunction with GIS. This technique was used by Singh et al. (2014) to evaluate the danger of erosion in the Shivalik Hills of Himachal Pradesh. Additionally, MDA has been utilized to map the susceptibility to gully erosion at the regional level in India. Chaudhuri et al. (2017), for instance, employed MDA in Uttarakhand to create a susceptibility map that included rainfall intensity, slope, land cover, and soil characteristics as discriminant factors. Gully erosion is still a major problem in West Bengal because of the region’s varied geomorphological features, heavy rainfall, and extensive agricultural operations, especially in the Rarh region and Birbhum district. Gully erosion intensity and susceptibility have been mapped and evaluated using geospatial approaches more often in recent research. The use of GIS-based multi-criteria decision analysis (MCDA), which incorporates many parameters such as soil type, rainfall patterns, and land use/land cover (LULC) to evaluate gully erosion hazards, is one noteworthy technique in West Bengal. For instance, a GIS-MCDA model was used in the Birbhum district study by Bhattacharyya et al. (2020) to identify locations susceptible to gully erosion by considering both natural and anthropogenic causes. This approach highlighted how deforestation and agricultural development increase erosion susceptibility by using remote sensing data to track changes in LULC over time. Slope analysis and the use of DEMs have also proved crucial in comprehending the dynamics of gully erosion in this area. High slopes and unregulated agriculture have been linked to increased erosion, according to studies (Mukherjee et al., 2018). Susceptibility models have been developed using MDA and statistical techniques, and LULC change detection using satellite images (e.g., Landsat and Sentinel-2) has assisted in identifying regions in Birbhum where gully erosion is most common. In West Bengal, machine learning methods are also being investigated more and more for the classification of areas that are prone to erosion. In the Murarai region of Birbhum district, Roy et al. (2019) classified land segments using some discriminant criteria using the Random Forest (RF) ensemble machine learning technique. The study discovered that machine learning-based soil erosion risk models outperformed conventional techniques in locating regions susceptible to gully formation. A variety of geospatial approaches, including remote sensing, GIS, MDA, MCDA, and machine learning, have been adopted worldwide, in Indian, and in West Bengal research to evaluate the severity and susceptibility of gully erosion. These methods have been successful in mapping gully erosion-prone areas, identifying important factors like topography, rainfall, and changes in land use, and enabling focused mitigation measures. Recent patterns indicate a significant move toward integrated strategies that incorporate many data sources and sophisticated algorithms, providing increased accuracy in identifying and controlling the hazards of gully erosion.

REVIEW OF LITERATURE OF RECENT STUDIES

Geospatial methods and machine learning (ML) algorithms have been used more frequently in recent research on gully erosion in India to evaluate vulnerability and forecast future erosion hazards. Gully erosion is still a major environmental problem that affects ecosystems, water quality, and agricultural output. The combination of geographic information systems (GIS), machine learning (ML), and remote sensing offers sophisticated tools for identifying erosion hotspots and creating focused conservation plans. A study conducted in the Birbhum area of West Bengal by Mondal et al. (2023) evaluated the susceptibility of gully erosion using factor analysis and GIS-based mapping. Their study found that Bolpur-Santiniketan and Rampurhat-I are particularly susceptible to gully erosion, with significant lateritic soil loss during the monsoon season. In aggravating soil degradation, they discovered that topographic elements like rainfall and slope were crucial (Mondal et al., 2023). Ensemble approaches and machine learning models are crucial for increasing the accuracy of gully verosion susceptibility mapping (GESM). With a high AUC value of 0.972, Lee et al. (2020) used bagging techniques and Boosted Regression Trees (BRT) to forecast gully erosion in Garhbeta, West Bengal. According to Lee et al. (2020), this model illustrated the value of machine learning for accurate predictions of gully erosion. Similar strategies have been effectively implemented throughout India, especially in areas like the northwest Himalayas and Western Ghats, employing Random Forest (RF) and Support Vector Machines (SVM). Additionally, forecast accuracy has increased when environmental characteristics like land use, slope, rainfall, and soil type are combined with geospatial and machine learning approaches. Studies in the Western Ghats and northwest Himalayas have shown that land use changes, particularly deforestation and intensive agriculture, have increased the risks of gully erosion. Planning for soil conservation is aided by these variables when they are incorporated into predictive models to produce susceptibility maps particular to a given region (Lee et al., 2020). To forecast erosion hotspots in Maharashtra, a mix of GIS, remote sensing, and machine learning techniques has been employed. By identifying high-risk regions, this method has been crucial in facilitating improved land reclamation techniques, like afforestation and soil stabilization (Ghosh & Bhattacharya, 2012). The understanding of gully erosion in India has been significantly enhanced by several additional studies. To illustrate how slope and land cover have a major impact on erosion patterns, Pandey et al. (2022) created a multi-factor model combining GIS and remote sensing to evaluate gully erosion susceptibility in the Deccan Plateau. In their study of gully erosion in the northeastern Indian states, Sharma et al. (2021) combined remote sensing and spatial analysis and found that erosion was primarily driven by topography and rainfall intensity. By using machine learning to simulate erosion susceptibility in the Himalayas, Sharma et al. (2021) discovered that forest cover, land use, and slope were all important predictors of erosion risk, providing information for focused erosion control strategies. Similar to this, Sharma (2023) assessed gully erosion hotspots in Uttarakhand using GIS and remote sensing techniques, showing that harsh erosion was caused by human-induced activities like mining and deforestation. These kinds of studies have a significant impact on the evolution of erosion control strategies throughout India. Through the integration of environmental elements, machine learning, and geospatial data, these methods present promising ways to mitigate soil erosion, boost agricultural productivity, and encourage sustainable land use practices in vulnerable areas.

Review of Literature on the Adopted Methodologies LULC Classification

Researchers have previously used a variety of techniques to assess soil erosion, including gullies. Seutloali et al. (2016) identified the gully erosion along major armoured highways in the southeast region of South Africa using a Geographic Information System (GIS) and RS methodology. Ayele et al. (2020) used the Revised Universal Soil Loss Equation (RUSLE) and the Sediment Yield Index (SYI) models to quantify rainfall- induced soil erosion in Ethiopia’s highlands. Mosavi et al. (2020) compiled various analytical models to determine the Talar watershed’s flood and erosion susceptibility. Pathan and Sil (2020) used a combination of remote sensing, GIS, and a soil erosion assessment method to identify the soil erosion-prone areas in the upper Brahmaputra river basin to Majuli river island. To identify the basin morphometry and gully erosion in eastern India’s lateritic badland zone, Ghosh and Kundu (2022) evaluated digital elevation models and topographic indices. In the Harda region of the Narmada River basin in India, the comparison of the Morgan-Morgan- Finney, Universal Soil Loss Equation, and Rectified Universal Soil Loss Equation models brought to light the zones of spatial variability and the calculation of soil erosion risk (Mondal et al. 2018). Gayen et al. (2020) used the Frequency ratio (FR) model to validate RUSLE to evaluate the soil erosion assessment in the Pathro river basin in Jharkhand, India. The maximum likelihood approach was used to classify LULC (Richards, 2022). Based on many important factors, the Maximum Likelihood Classification (MLC) was used to evaluate the severity and susceptibility of gully erosion in the Birbhum area of West Bengal. First, MLC’s statistical method works well for processing remote sensing data, where the spectral properties of soil, vegetation, and land cover types are important (Richards & Jia, 2006). It assumes a normal distribution for pixel values within each class. Furthermore, as research by Singh et al. (2022) and Mondal et al. (2023) have shown, MLC is perfect for complicated terrain analysis since it can incorporate many geographic data layers, including topography, soil type, and rainfall. Furthermore, because it produces probabilistic classification findings, MLC is useful in mixed-use landscapes where various land coverings (such as agricultural and forest) interact. This ensures more precise identification of regions susceptible to gully erosion (Gao et al., 2020; Das et al., 2021).

Furthermore, it is beneficial in similar erosion studies conducted throughout India, including the Western Ghats, where it aided in mapping erosion susceptibility according to land use and terrain factors (Gao et al., 2020). Consequently, MLC is a reliable technique that enables thorough, precise mapping of gully erosion hazards, offering important data for Birbhum soil conservation and land management plans (Tso & Mather, 2009). Due to its statistical robustness and suitability for remote sensing data, the author of the study on gully erosion assessment using geospatial techniques in Birbhum district mainly chose Maximum Likelihood Classification (MLC). Other algorithms, such as Support Vector Machine (SVM) and Logistic Regression, were not investigated. MLC was selected because it particularly performs well in scenarios where the distribution of data classes is assumed to follow a normal (Gaussian) distribution, which is frequently the case in remote sensing data about land use and soil properties, even though SVM and Logistic Regression are well- liked in classification tasks (Richards & Jia, 2006). MLC’s significance stems from its capacity to efficiently categorize different types of land cover by processing various bands of satellite imagery, such as multispectral or hyperspectral data. It is well-suited for differentiating between land cover types in areas impacted by gully erosion, where minute changes in spectral data can reveal erosion-prone areas because it calculates class probabilities using both the mean and the variance of the pixel values (Lillesand & Kiefer, 2000). According to Lu et al. (2014), MLC is also beneficial in complicated terrains since it considers a range of environmental elements that contribute to gully erosion, such as soil type, vegetation, and slope. In contrast, Support Vector Machines (SVM) are computationally more demanding and generally require larger training datasets, notwithstanding their effectiveness for many classification applications. However, the intricacy of multi-class or non-linear connections between the variables involved in gully erosion processes may be unaccounted for by logistic regression (Pal, 2005). A level of classification confidence is also provided by MLC’s probabilistic output, which is crucial when working with complex and diverse land covers that are present in regions that are vulnerable to soil erosion. Since MLC has demonstrated reliability, flexibility in handling multi-band data, and the capacity to incorporate topographic and environmental factors that are critical in understanding gully erosion, it remains the best fit for this study, even if other algorithms may have been evaluated (Foody, 2002).

The performance of gully erosion susceptibility models must be evaluated in terms of accuracy, and one often used indicator in environmental modeling is Kappa statistics. The percentage of examples (both positive and negative) that are correctly identified out of all the instances is known as accuracy. A measure of agreement between observed and anticipated classifications that considers chance agreement is the Kappa Statistic (κ). The range of values for kappa is -1 to 1. By taking chance agreement into account, Kappa statistics, which Congalton and Green (1999) established, provide a more reliable metric than total accuracy for evaluating the agreement between observed and anticipated classifications. Higher Kappa values suggest better model performance, as evidenced by Pal and Mather’s (2005) use of Kappa statistics to evaluate land degradation models in erosion investigations. Since Kappa statistics offer a more thorough understanding of model correctness than just classification rates, its use in gully erosion susceptibility models has proved crucial in guaranteeing the accuracy of predictions, particularly when working with categorical outcomes like erosion- prone zones. Congalton and Green’s (1999) technique was used in this work to evaluate the kappa statistics and correctness of LULC data. Congalton and Green’s (1999) study offers a thorough description of the use of several accuracy metrics in environmental studies and remote sensing, such as overall accuracy and Kappa statistics.

Drainage Basins and Stream Networks

The specific stream network and drainage basin have been delineated based on the following procedures using digital elevation model (DEM) data. After successful extraction of the flow direction raster, the drainage basin has been delineated using the basin tool in ArcGIS. These tools are commonly used in hydrological studies to understand watershed behavior, water flow patterns, and to assist in water resource management (Esri, 2020; Zhang & Li, 2019). Additionally, Strahler’s approach (Strahler, 1964) has been used to prepare a stream raster and a flow direction raster for the ordering of streams. It serves as a gauge for a stream’s place in the hierarchy (Leopold, 1994). The total number of streams (Nu, Horton, 1945) and the total length of the streams in kilometer (km) (Lu, Strahler, 1964) have been calculated in each order. The Bifurcation ratio and mean bifurcation ratio have been calculated for both delineated basins using the formulation of Horton (1945).

Statistical Analyses

The indicators in this study have been normalized using the Min-max normalization method (Nardo, 2009; Joint Research Centre-European Commission, 2008). Before estimating susceptibility indices using Principal Component Analysis (PCA), the indicators are standardized using Max-Min Normalization. To prevent the disproportionate influence of variables with greater scales, such as rainfall or terrain, this technique is important since it guarantees that all indicators are rescaled to a uniform range, usually between 0 and 1 (Jolliffe, 2002). In geographic analyses involving numerous environmental and topographical variables, the Max-Min method is especially helpful since it is straightforward and efficient for datasets with varying ranges and allows for an impartial contribution from each element (Zhang et al., 2019). According to Sun et al. (2017), it is perfect for huge datasets that are common in erosion susceptibility and remote sensing investigations since it is computationally efficient and produces an output that is easy to understand. The study of gully erosion relies on the reliable extraction of principal components from a variety of indicators, and Max-Min normalization was chosen for its easy application and efficient scaling, even though other normalization techniques like Z-score normalization or decimal scaling could also be used (Zhang et al., 2019; Sun et al., 2017). Thus, The Principal Component Analysis formula (PCA, Pearson, 1901) has been applied to extract the factor scores. The Carlin and Doyle (2000) formula has been used to estimate the standard error mean. Carlin and Doyle (2000) state that the sample Standard Deviations, or SD (s), must be used to estimate the Structural Equation Model (SEM) rather than the unknown σ. The standardized values of the covariances of the components influencing gully erosion are measured by the Structural Equation Model (SEM, Wright, 1918). The structural equation model has been used in the current work to depict the correlation between a few chosen variables. In this case, the structural equation model’s conceptual diagram displays the covariances that were recovered using linear regression. A significance criterion of 0.05 (95% CI) or 0.01 (99% CI) is used to validate the covariances. The current study additionally makes use of the covariance (Engelhart, 1941) and linear regression model formulas (Pearson, 1914; Pearson, 1897).

Model validation of Gully Erosion Susceptibility Zones

To assess the predictive effectiveness of models, gully erosion susceptibility studies have frequently used the Receiver Operating Characteristic (ROC) curve and its associated Area Under the Curve (AUC) statistics. ROC-AUC compares true positive rates (sensitivity) against false positive rates (1-specificity) across several thresholds to determine how well a model can differentiate between areas that are prone to erosion and those that are not. A higher AUC denotes superior discriminatory ability. The AUC measures the model’s overall performance. The range of the AUC is 0–1. Radar performance was initially assessed using the ROC curve, which was first introduced by Green and Swets (1966) in signal detection theory. The AUC was later widely used in many domains, including environmental modeling, when Metz (1978) popularized it as a metric for diagnostic accuracy. Due to its adaptability and dependability, ROC-AUC has become crucial in erosion studies for evaluating the effectiveness of susceptibility models. Global research has shown that AUC values above 0.7 suggest reasonable model accuracy, whereas values above 0.9 indicate exceptional performance (Rahmati et al., 2017; Wang et al., 2019). The significance of hydrological and topographical variables in erosion modeling has been demonstrated by the high AUC values obtained by models that included these variables in areas like the Loess Plateau and the Ethiopian Highlands. When using ROC-AUC to evaluate gully erosion susceptibility in semi-arid parts of India, Pandey et al. (2020) discovered that model reliability was improved by combining vegetation indices with the topographic wetness index (TWI). Key indicators including the Stream Transport Index (STI), TWI, and Infiltration Number may be validated in the current study using ROC-AUC analysis, guaranteeing the resilience of susceptibility models and their suitability for a variety of geological and geomorphological contexts. Based on Fawcett’s (2006) research, ROC-AUC was employed in this investigation. An extensive discussion of ROC analysis and AUC, as well as how they are used to assess classification models, is given in this research by Fawcett (2006).

Research Gap

Studying the Gully Erosion in West Bengal’s Birbhum area is constrained by some limitations. First, although the study has looked at the broad topography and climate elements that affect the formation of gullies, there isn’t much in-depth analysis that considers the soil properties of the lateritic terrain in this area. In earlier research, the poor water retention, acidity, and nutrient depletion of lateritic soils, along with their susceptibility to erosion, were left out. Second, although micro-watersheds are more pertinent for localized gully erosion assessments, most research has concentrated on broader geographical scales. A significant gap is the lack of research using high-resolution geospatial data to comprehend the interactions between rainfall erosivity and several morphometric parameters, including aspect, drainage density, slope, and land use, in these micro-watersheds. Additionally, while GIS and remote sensing have been used to study changes in land cover and use, few studies have integrated these technologies with sophisticated machine learning algorithms to forecast future hotspots and simulate the vulnerability to gully erosion. The impact of shifting land use patterns, especially deforestation and agricultural intensification, on gully erosion in Birbhum is also not given enough attention. One significant gap in the literature currently in publication is the scant incorporation of both natural and human-induced changes in forecasting and evaluating the susceptibility of gully erosion.

Need for the Present Study

The crucial gaps in literature, especially concerning gully erosion in Birbhum area, make the current study necessary. First, by concentrating on Birbhum’s distinct lateritic soil properties—which have received little attention in gully erosion studies—this study will offer crucial new information about how these soil attributes affect erosion susceptibility. Second, by focusing on micro-watersheds, the study seeks to close the spatial resolution gap and provide a more detailed knowledge of the dynamics of gully erosion at the local level. A more thorough and precise evaluation of the variables influencing gully development will be possible by combining rainfall erosivity with morphometric variables including slope, drainage density, and topography. A major improvement over conventional approaches, the study will also produce accurate erosion susceptibility maps using high-resolution remote sensing data and GIS-based machine learning techniques. The study will also examine how human activities like urbanization, deforestation, and agricultural growth contribute to the region’s gully erosion. To help policymakers and local authorities create more effective erosion control measures, this study will combine natural and human-induced elements to provide specific, doable recommendations for sustainable land management practices.

Objectives

The present study focused on the two blocks of Birbhum district of West Bengal, i.e., Rampurhat-I and Bolpur- Sriniketan blocks. The district is an extended part of Chota Nagpur plateau and Rajmahal basalt trap covered with lateritic patches. Highly leached and weathered tropical lateritic soil capes are enriched with oxides of iron and aluminum and are ecologically fragile because of their inherent constraints of acidity, nutrient loss, chemical impairment, crusting, water erosion, and poor water-holding (Jha, 2008). The physical characteristics of this district are adjacent to its erosion potential of gully formation areas. In the present study, various morphometric parameters and raster analysis indices are widely analyzed and correlated with the rainfall erosivity to identify the significantly influential factors, the spatial distribution of the factors, delineate gully erosion susceptibility zones, and access the significant impact of gully erosion on land use and land cover in the context of micro watersheds in the study area. The present study focuses on the following objectives:

1. To measure the role of diverse controlling factors of gully erosion in Rampurhat-I and Bolpur- Sriniketan blocks of Bibhum district.
2. To formulate the relationship among the controlling factors of gully erosion and delineate gully erosion susceptibility zones in the study area.
3. To assess the impact of gully hotspot areas on the land use land cover of the study

MATERIALS AND METHODS

Study Area

The present study focuses on Rampurhat-I and Bolpur-Sriniketanblocks situated in southern and western portions respectively of Birbhum district. Rampurhat-I block is situated at an extension of 87̊ 32̍ East (E) to 87̊ 52̍ East (E), and 24̊ 05̍ North (N) to 24̊ 17̍ North (N), and Bolpur-Sriniketan block lies in between 87̊ 35̍ E to 87̊ 50̍ E, and 23̊ 33̍ N to 23̊ 54̍ N (Fig 1). Geographically, this region is a part of the Rarh region, positioned in the western part of West Bengal and lies at the north-eastern end of the Chota Nagpur Plateau and its slope gradually downs towards east. The lateritic region is a part of the low-level unconsolidated erosional deposits from the eastern Chota Nagpur plateau where regional climate varies from harsh in the west to relatively mild in the east. Birbhum district belongs to several morphological processes, i.e., weathering, mass wasting, and river erosional process (Mondal, 2013). Degraded lateritic soils occupied a large area of this district containing oxides of iron and aluminum (Jha and Kapat 2011). Different morphometric analyses and geospatial and multivariate statistical techniques have been adopted by various scholars to identify the gully erosion impact areas of Birbhum district. Ghosh and Bhattacharya (2012) used different analytical aspects to identify the soil loss due to erosion in the badland areas of the western Rampurhat-I block in Birbhum district. Jha and Kapat (2009) identified 45 micro watersheds in Ajay-Maurakshi interfluves in the south-western Birbhum district. This area consisted of various rill and gully erosion sites impacted by different morphometric parameters of river basins associated with relief, slope, and drainage characteristics (Jha and Kapat 2009). The study area, Rampurhat-I, and Bolpur-Sriniketan blocks are occupied such rill and gully erosion areas which are degraded, dissected, and rugged over the lateritic capes impact on the land use land cover, agriculture, and habitation in Birbhum district. To control the gully head areas, check dams and afforestation are required to be established (Ghosh and Dolui 2011; Das, 2015). The total population and total households of Rampurhat-I were 159,193 and 31,040 in 2001 and 188,435 and 44,263 in 2011 respectively (Census of India 2001; 2011). The total population and total households of Bolpur-Sriniketan were 178,111 and 37,280 in 2001 and 202,553 and 47,961 in 2011 respectively (Census of India 2001; 2011). According to available data, the current population of Birbhum district in West Bengal, India is approximately 4,060,000 people.

Data Sources

The study has been conducted using secondary databases and participants’ observations in the study area. The relevant remotely sensed (RS) satellite data have been collected from the United States Geological Survey (USGS), National Remote Sensing Centre (NRSC), and European Space Agency (ESA) websites. Geospatial techniques and statistical methods have been employed for data analysis, mapping, and representation. The details of the databases collected are mentioned in the Table (Table 1). The methodology of the present study is represented in a figure (Fig 2) in a generalized framework.

CALCULATION OF LAND USE LAND COVER STUDY

The land use land cover (LULC) maps have been prepared for both Rampurhat-I and Bolpur-Sriniketan blocks using specific RS satellite data (mentioned in Table 1). After the rectification of collected satellite images, LULC maps have been prepared based on the Land Satellite 7 (LANDSAT 7) images Enhanced Thematic Mapper Plus band (ETM⁺ band) of 2001 (USGS, 2022a). In this process, the maximum likelihood method of supervised classification has been adopted to identify the five major land use land cover classes along with gully erosion areas. The extracted LULC classes are water bodies, natural vegetation, barren land, agricultural land, and build-up areas. The LULC classification was completed by the maximum likelihood method, using the following formula,

(𝑥) 𝑙𝑛_𝑝(𝜔_𝑖) − 1⁄2 𝑙𝑛| 𝛴_𝑖| − 1⁄2 (𝑥 − 𝑚_𝑖)^𝑇 𝛴⁻¹(𝑥 − 𝑚_𝑖)                                                                                                                                      (1)

where

i = class

x = n-dimensional data (where n is the number of bands)

p(ωi) = probability that class ωi occurs in the image and is assumed the same for all classes

|Σi| = determinant of the covariance matrix of the data in class ωi Σi-1 = its inverse matrix

mi = mean vector. The LULC raster map has been undergone in accuracy assessment by calculating kappa statistics.

𝜅 = 𝑃0−𝑃𝑒

1−𝑃𝑒

where

𝜅 = Kappa statistics

Po = Observed accuracy (proportion of correct predictions)

Pe = Expected accuracy (proportion of correct predictions expected by chance)

Delineation of Drainage Basins and Stream Networks (2)

The specific stream network and drainage basin have been delineated based on the following procedures using digital elevation model (DEM) data. After successful extraction of the flow direction raster, the drainage basin has been delineated using the basin tool in ArcGIS. Drainage basins are defined in ArcGIS by using the Basin tool, which is in the Spatial Analyst extension after the flow direction raster has been extracted. The Hydrology toolbox’s Flow Direction tool first creates a flow direction raster, which indicates the way water flows from each cell to its steepest neighbor. The Basin tool utilizes this generated raster to locate and mark the limits of distinct drainage basins. This procedure can also be used in conjunction with the Watershed tool, which further refines the basin delineation by determining the contributing area for each cell based on the direction of flow. To create the stream network, the flow accumulation raster has been extracted from the flow direction map, and then the stream raster has been converted into a polyline to acquire the specific stream network. Raster classes of flow accumulation are greater than 5000. In these processes, Suttle Radar Topographic Mission- DEM (SRTM-DEM) (USGS, 2022b) data are used for both Community Development Blocks (C.D. blocks) of the Birbhum district. Moreover, the ordering of streams has been prepared using a stream raster and a flow direction raster using Strahler’s method. The total number of the streams and the total length of the streams in kilometers (km) have been calculated in each order. The Bifurcation ratio and mean bifurcation ratio have been calculated for both delineated basins in Rampurhat-I and Bolpur-Sriniketan blocks. The formula of the bifurcation ratio is

𝑅𝑏 = 𝑁𝑢⁄𝑁𝑢 + 1

where

Rb is the bifurcation ratio

Nu is the number of streams of any given order (3)

Nu+1 is the number in the next higher order. The mean bifurcation ratio (Rbm, Strahler, 1957) is derived as the average of the bifurcation ratios of all orders.

The total area of delineated drainage basins in Rampurhat-I and Bolpur-Sriniketan have overlapped on 10 by 10 grids. After extraction totals of 55 points in Rampurhat-I and 52 points in Bolpur-Srtiniketan blocks (Figs 3 & 4) are found with their coordinates.

Formulation of Gully Erosion Susceptibility Index

The composite indices method has been implemented through Factor analysis of Principal Component Analysis (PCA) twenty indicators of gully erosion to bring out the composite gully erosion susceptibility indices using the standardized predicted scores. The details of the indicators have been mentioned in Table 2. Min-max normalization method has been applied to normalize the indicators.

𝑥′ =             𝑥−𝑚𝑖𝑛(𝑥)

               𝑚𝑎𝑥(𝑥)−𝑚𝑖𝑛(𝑥)

where,

x is the indicator.

For the Composite Factor Analysis, the following formula has been used as (4)

P1 = Σ aj1 X Zj or P1 = a11. Z1+ a21. Z2 a11. + … an1. Zn                                            (5)

where,

P1 denotes the composite gully erosion susceptibility index of a unit study as the first factor denotes the factor

loading of the ‘j’^th variable and 1 indicates the factor number that is the first factor- vector of factor loadings. While the Zj denotes the standardized value of the ‘j’^th variable, which is expressed as

Zj= 𝑋𝑗−𝑋𝑚

         𝛿𝑗

Where (6)

Xj denotes the original value of ‘j’^th variable, Xm denotes the mean (simple arithmetic mean) of ‘j’^th variable, and δj denotes the standard deviation of ‘j’^th variable.

In this aspect, the mean and standard deviation are calculated by using the following formula:

Mean= ^𝛴𝑥𝑛 (7)

Standard Deviation=√(𝑥−𝑥̅

                                                𝑛           (8)

where

x̄ is the arithmetic mean; x is the individual value of items; n is the number of terms in the distribution. The Standard Error Mean has been estimated following the formula.

SEM = 𝑠√𝑛

(9)

Finally, Mean Composite Factor Scores have been calculated using the standardized factor scores extracted from PCA.

Mean Composite Factor Scores = 𝐹𝑎𝑐𝑡𝑜𝑟1+𝐹𝑎𝑐𝑡𝑜𝑟2+ 𝐹𝑎𝑐𝑡𝑜𝑟3…+𝐹𝑎𝑐𝑡𝑜𝑟𝑛

                                                             𝑇𝑜𝑡𝑎𝑙 𝑏𝑢𝑚𝑏𝑒𝑟 𝑜𝑓 𝑓𝑎𝑐𝑡𝑜𝑟𝑠 (10)

where

𝑛 is the factor.

The mean composite scores of both areas have been represented in the inverse distance weighted (IDW) maps and the maps have been validated by Receiver Operating Characteristic – Area Under the Curve (ROC-AUC).

Structural Equation Model and Covariances

The Structural Equation Model has been constructed to measure the standardized values of the covariances of the factors that impact gully erosion. In the present study correlation among selected variables has been represented through the structural equation model. Here, covariances have been extracted using linear regression and put on the conceptual diagram of the structural equation model. The covariances are validated with a significance level of either 0.05 (95% confidence interval) or 0.01 (99% confidence interval). The formula for the linear regression model is

Yc = a + b X                                                                                                                                    (11)

where

Yc is a predicted value of Y (which is the dependent variable)

a is the Y-intercept

b is the changes in Y for each increment change in X

X is an X score (Independent variable) for which a value of Y is predicted.

The formula of covariance is Cov (x,y) = 𝛴(𝑥𝑖−𝑥̅)(𝑦𝑗− ȳ)

                                                                                𝑛                           (12)

where

Cov (x,y)              = Covariance between variable x and y

𝑥_𝑖                          = Data value of x

𝑦_𝑗                          = data value of y

𝑥̅                           = mean of x

ȳ                           = mean of y

n                           = Number of data values

RESULTS

Physiographic Situations of the Study Area

This entire region has a significant physical setup. Rajmahal basalt, china clay, and laterite (belonging to the Cainozoic era) have covered up this region. Pascoe (1973) researched the facts in a manual of the geology of India and Burma. Laterites of Birbhum were split off from the high-level laterites of eastern Rajmahal hills and carried out to this eastern lower region of Rajmahal hills by rainwater, small streams as well as surface runoff, and these materials were re-deposited in this area. After successively wet and dry climatic conditions and the corresponding result of the change of groundwater level (in pre-monsoon the average water level is above 10 meters depth in April and in post-monsoon that is nearly 3 metres) in this region, the oxides of aluminum and irons had compacted repeatedly. Morphologically this region belongs to an undulating lateritic upland area. This region is mostly covered with tropical dry deciduous forests along with some evergreen trees. Shrubs and grassland have been grown in the western part of this region such as thorny shrubs, Aceraceae, and Mangifera (Ghosh, 2011). Different morphometric parameters have been figured out to identify the physiography of the study area such as elevation, ruggedness, slope, slope aspect, direction, stream length density, stream distance, basin flow, flow distance, and bifurcation ratio. The total area of the basin delineated in Rampurhat-I block is 295.658 square km and it is 341.744 square km in Bolpur-Sriniketan block. A total of four orders of streams have been identified in both areas. The mean bifurcation ratios of the two delineated basins are 4.15 and 3.46 which indicate that the areas are highly probable for surface runoff which aggravates the erosion of surface and sub-surface soil layers to a greater extent. Besides, the elevation and slope of those areas are increasing from east to west. The direction of the channels from west to east impacts the transportation of eroded materials from west to east originating from the head-cut erosional process. To comprehend the research area’s vulnerability to gully erosion and soil degradation, it is essential to grasp its physiography (Figs 5 & 6). Geological deposits from the Cainozoic era, such as laterite, china clay, and Rajmahal basalt, compose the region’s main physical structure. In A Manual of the Geology of India and Burma, Pascoe (1973) reported that lateritic materials in Birbhum were carried by surface runoff, small streams, and rainfall from the upper Rajmahal Hills to the lower eastern parts of the Rajmahal Hills, where they were then redeposited. These materials experienced compaction over time because of alternating dry and wet weather patterns. With the water level normally reaching more than 10 meters in depth in April during pre-monsoon periods and falling to about 3 meters depth during post-monsoon periods, the shift in groundwater levels is substantial. The area’s vulnerability to erosion was exacerbated by the accumulation and compaction of iron and aluminum oxides in the soil because of the frequent wet and dry circumstances. Erosion-prone features are further exacerbated by the region’s morphology, which is classified as an undulating lateritic upland terrain. A range of evergreen trees and tropical dry deciduous forests make up much of this area. But the western part of the area is mostly grassland and shrubby, with common plant species like Aceraceae, Mangifera, and prickly shrubs (Ghosh, 2011). These varied ecosystems help to stabilize the soil, but compared to evergreen forests, the dry, deciduous vegetation is not as good at stopping erosion during the dry season. The physiography of the research area has been established and studied using a variety of morphometric criteria, including elevation, roughness, slope, slope aspect, and stream length density. The entire basin area of the Bolpur-Sriniketan block is 341.744 square kilometers, whereas the Rampurhat-I block is 295.658 square kilometers. Both regions have been divided into four stream orders, revealing a varied waterway network that influences erosion and surface runoff. The basins are likely to experience considerable surface runoff, as shown by the comparatively high mean bifurcation ratios for the two areas, which are 4.15 for Rampurhat-I and 3.46 for Bolpur-Sriniketan. This excessive runoff exacerbates the erosion process overall by playing a significant role in the erosion of both surface and subsurface soil layers. Furthermore, these areas’ slopes and elevations rise from east to west, causing water to move from the western to the eastern regions while carrying eroded material with it. The movement of eroded soil is improved by this directed flow pattern, particularly from head-cut erosion zones in the western portions of the study region.

Soil Loss in the Gully Formation Areas

This region belongs to oxisol laterite, which is more convinced by weathering and erosion by climatic properties. This soil category is one of the significant aspects that induced the erosion procedure (Getnet et al. 2021). Soil erodibility is very high in this region. The soil contains many silt-sized particles which are effectively susceptible to erosion by the surface and subsurface flow (Fig 7). This Kankara, which is the local name of this soil, is unconsolidated, fragile, and extensively susceptible to erosion in nature. The amount of soil loss in the study area is above 50%. The mean rainfall intensity, soil erodibility, and the process of originating gully hotspots highly increase the amount of soil loss in the study area (Fig 7). The subsoil areas of lateritic channels are eroded and degraded by the transportation process of materials along the gully channels (Fig 7). An important factor in the research area’s susceptibility to erosion is the soil’s properties. The oxisol laterite soils that predominate in the area are subject to erosion and weathering because of the local climate. Silt-sized particles that are easily separated by water flow are abundant in these soils, which are naturally brittle and unconsolidated. Because of their erodible and unstable characteristics, these soils are known locally as Kankara. Because of their physical characteristics and the dynamic weathering processes they experience, these soil types are particularly vulnerable to erosion, citing research by Getnet et al. (2021). Because precipitation readily separates soil particles and carries them through surface and subsurface flows, the region’s high soil erodibility makes gully formation more difficult. The research region has more than 50% soil loss, which suggests that erosion has significantly reduced the amount of fertile topsoil. The main causes of this high degree of soil loss are the development of gully hotspots, soil erodibility, and mean rainfall intensity. The erosive power of water flows increases with rainfall intensity, further loosening soil particles and hastening the formation of gullies. Because silt-sized soils are more likely to be taken away by surface runoff and subsurface flows, causing significant soil deterioration, soil erodibility is especially important. The region’s level of soil loss is further increased by the development of gully hotspots, which are regions that are more vulnerable to erosion because of concentrated runoff. Erosion is not restricted to the soil on the surface; it also affects the lateritic channels’ subsurface layers, which progressively erode and deteriorate over time. Because eroded materials are transported to lower places and frequently result in the deposition of sediments in downstream areas, the movement of materials along gully channels exacerbates soil deterioration. The long- term reduction in soil fertility and the general deterioration of the landscape are both caused by this ongoing process. The study’s figures, such as Fig 7, demonstrate the degree of soil loss and the deterioration of gully formation areas, demonstrating the further erosion of soil layers and the deepening of gullies.

Factors of Gully Erosion and Gully Erosion Susceptibility Zones

Analysis of Descriptive Statistics and Factor Analysis

The Rampurhat-I block’s descriptive statistics (Table 3) provide important information about the environmental and morphometric factors affecting gully erosion. Significant variation in terrain wetness and its direct influence on water buildup and soil erosion is highlighted by the Topographic Wetness Index (TWI), which has the highest standard deviation at 242.94 and a mean of -72.40. This element emphasizes how diverse the terrain is and how susceptible it is to erosion. The Normalized Difference Vegetation Index (NDVI) shows substantial spatial variability in vegetation density, with a mean of 0.0925 and a standard deviation of 0.117. Since vegetation is essential for supporting the soil, gully formation susceptibility is correlated with its mild variance. With a mean of 23.99°C and a low standard deviation of 1.12, the Land Surface Temperature (LST) suggests comparatively consistent thermal conditions that could affect soil moisture content and evapotranspiration, which in turn could have an indirect impact on erosion. With a mean of 0.0620 and an exceptionally low standard deviation of 0.0538, the drainage density (DD) on the hydrological front demonstrated regular drainage patterns that amplify surface runoff, a key contributor to gully erosion. Standard deviations for other morphometric indicators, such as the Slope (S) and Stream Power Index (SPI), were 2.01 and 0.133, respectively, indicating the influence of topography gradient and river flow energy on erosion patterns. Furthermore, indices that measure the degree of exposed soil and moisture content variability—both of which are closely related to erosion processes—such as the Bare Soil Index (BSI) and Modified Normalized Difference Water Index (MNDWI), which have mean values of 252.37 and -0.0827, respectively, offer vital information. Factor analysis provides additional insight into the interactions between these variables, as seven main components together accounted for 77% of the variance. TWI, drainage density, and NDVI were the main factors, highlighting the interplay between plant cover, hydrological dynamics, and terrain shape in determining the likelihood and severity of gully erosion in the Rampurhat-I block. Critical environmental and morphometric parameters that contribute to gully erosion are also highlighted by the descriptive statistics for the Bolpur-Sriniketan block (Table 4). In contrast to Rampurhat-I, which has a significant impact on soil saturation and erosion concerns, the Topographic Wetness Index (TWI) once again showed the largest standard deviation of 505.88 with a mean of -110.40, suggesting even more substantial terrain wetness variability. With a mean of 0.1082 and a standard deviation of 0.147, the Normalized Difference Vegetation Index (NDVI) indicated moderate to high variability in vegetation density, which affects soil stability and erosion susceptibility. The Land Surface Temperature (LST), on the other hand, showed consistent temperature conditions comparable to Rampurhat-I, with a mean of 24.08°C and a standard deviation of 1.25. Compared to Rampurhat-I, the Drainage Density (DD), which had a mean of 1.1431 and a standard deviation of 0.987, showed more variability. This suggests that drainage networks are not evenly distributed, which can exacerbate runoff and soil erosion in some places. Stream Power Index (SPI) and Slope (S) were the most variable hydrological parameters, with standard deviations of 0.429 and 1.28, respectively. This highlights the importance of these elements in controlling water flow and enhancing the erosive capability of surface runoff. Additional indices, such as the Modified Normalized Difference Water Index (MNDWI) and the Bare Soil Index (BSI), which had mean values of -0.0526 and 0.2329, respectively, emphasized the block’s exposed and moisture-deficient soil characteristics, which increased the risk of erosion. Seven components were identified as being responsible for 74% of the variance overall by the principal component analysis, which was marginally less than Rampurhat-I. The main causes of gully erosion in Bolpur-Sriniketan were TWI, slope, and NDVI, suggesting a complex interplay between vegetation patterns, water flow dynamics, and terrain imperfections. In comparison to Rampurhat-I, Bolpur-Sriniketan exhibits greater geographical heterogeneity and intensity of gully erosion variables, as the data highlights overall.

Total Variance Explained and Rotated Component Matrix

The Rampurhat-I block’s total variance explained (Tables 5 and 6) by principal component analysis offers important new information about the complex variables affecting gully erosion. Seven major components, which together accounted for 77.03% of the variation, were identified through analysis. The first factor alone accounted for 22.29% of the total, demonstrating its dominance in explaining the differences in erosion susceptibility. Land Surface Temperature (LST), Bare Soil Index (BSI), and NDVI were the main components of this component, indicating the significant influence of vegetation cover and thermal characteristics in reducing or accelerating erosion. The second factor, which contributed 15.88% of the variance overall, was highly correlated with hydrological variables including drainage density and the Stream Power Index (SPI), underscoring the role of network consistency and water flow energy in soil detachment processes. The substantial influence of elevation and water retention capacity on gully formation was highlighted by the third and fourth components, which accounted for 11.18% and 8.29% of the variation, respectively, and included terrain-related factors such as slope and Topographic Wetness Index (TWI). Components five through seven, which contained varying effects of minor vegetation indices and hydrological metrics, contributed progressively less to the variation, ranging from 7.48% to 5.18%. When taken as a whole, these elements highlight the intricate interactions between plants, hydrology, and terrain forms that influence the vulnerability of gullies to erosion. These relationships were further honed by the rotational analysis, which confirmed the significance of key variables in the erosion process by repeatedly displaying significant loadings for NDVI, drainage density, and BSI. A similar multi-factorial structure influencing gully erosion is revealed by the principal component analysis for the Bolpur-Sriniketan block (Tables 7 and 8), where seven components account for 74.32% of the variation. Vegetation-related indices like NDVI, BSI, and Land Surface Temperature (LST) dominated the first component, which explained 18.47% of the variance. This implies that soil exposure and vegetation cover are important factors that influence the danger of erosion. With a 14.22% contribution, the second component demonstrated the function of hydrological elements that control soil saturation levels and water flow patterns, such as drainage density and infiltration indices. The impact of slope and Topographic Wetness Index (TWI), which determines water retention and runoff dynamics on uneven terrain, was captured by the third and fourth components, which added 12.01% and 9.63%, respectively. The fifth component, which was linked to metrics like the Modified Normalized Difference Water Index (MNDWI) and highlighted the influence of moisture content and surface water distribution, accounted for 8.25% of the variance. With contributions of 6.25% and 5.50% respectively, the sixth and seventh components included small but significant factors such as flow length and the Stream Power Index (SPI), which represent the amount of water moving along gully channels and the erosive energy of runoff. The rotational analysis reaffirmed the importance of factors like TWI, drainage density, and slope as the main contributors to gully erosion susceptibility, even if the total variance explained is somewhat less than Rampurhat-I. According to these results, a slightly more diffused collection of variables influences the erosion processes in Bolpur- Sriniketan than in Rampurhat-I.

The relationships between the factors influencing gully erosion are clarified by the rotational component matrix (Table 5) for the Rampurhat-I block. Seven components with unique variable clusters were found by the analysis. For indices like the NDVI (-0.838), the Bare Soil Index (BSI, 0.903), and the Land Surface Temperature (LST, 0.925), the first component showed considerable positive loadings, indicating their dependency in affecting soil stability. This illustrates how vegetation and temperature play a crucial role in controlling erosion susceptibility. Stream Power Index (SPI, 0.841) and Topographic Wetness Index (TWI, 0.532) dominated the second component, indicating the impact of terrain wetness and hydrological energy. Drainage density (DD, 0.811) and infiltration indices demonstrated large loadings, highlighting drainage variables in the third component. Variables like directionality, flow length, and slope showed smaller but significant contributions to the remaining components, suggesting localized effects on erosion processes. Interpreting the intricate relationships between the environmental and morphometric components became easier because of the varimax rotation’s successful separation of the variables into logical groups. Similar information on the relationships between variables can be found in the Bolpur-Sriniketan block’s rotated component matrix (Table 7). The impact of vegetation and temperature fluctuations on soil erosion was highlighted by the first component, which showed substantial positive loadings for NDVI (0.794), BSI (0.924), and LST (-0.937). The second component, which focused on the function of water distribution and soil moisture absorption, recorded drainage-related data such as drainage density (0.923) and infiltration indices. Significant loadings for slope (-0.774) and flow length (0.563) in the third component, which concentrated on topography variables, showed that these factors directly affected the erosive potential of runoff. Localized factors such as MNDWI (-0.832) and SPI were underlined by subsequent components, indicating their distinct functions in erosive energy and moisture management. The variables were efficiently grouped into significant clusters using the rotational matrix, which clarified the multifaceted factors influencing gully erosion in Bolpur-Sriniketan. With an emphasis on managing vegetation, drainage networks, and terrain irregularities, this comprehensive understanding helps identify priority sites for erosion mitigation techniques.

Composite Factor Scores and Erosion Susceptibility Zones

Based on the extracted principal components, the composite factor scores for the Rampurhat-I block (Table 9) show geographic diversity in erosion susceptibility. Figs 8, 9, 10, and 11 show the indicators and spectral indices used to measure the susceptibility zones for gully erosion. Latitude 24.188 and longitude 87.698, which are in the northwest of the block and are characterized by steep slopes and strong rainfall erosivity, had the highest composite score (5.38). This illustrates how hydrological and topographical elements work together to cause extreme erosion. In contrast, places with flatter terrain and better vegetation cover are associated with lower composite scores (e.g., -5.41 at 24.240, 87.780), which decreases susceptibility. The block-wide mean composite score was -0.53, indicating that while overall vulnerability is minimal, tailored intervention is necessary for specific hotspots. The main determinants of these scores were TWI, SPI, and NDVI, suggesting that the most important variables affecting erosion in this block are vegetation density, flow energy, and terrain moisture. The regional variation in gully erosion vulnerability is further supported by the composite factor scores in the Bolpur-Sriniketan block (Table 10). The southeast portion of the block, which is marked by steep slopes, little vegetation, and concentrated runoff, is identified as a key hotspot for erosion by the highest score (5.308) in latitude 23.633 and longitude 87.722. The lowest composite score (-7.464), on the other hand, was found at 23.569, 87.831, which is a flat area with comparatively abundant vegetation that provides natural erosion resistance. The block as a whole faces modest erosion risks, although some parts are highly susceptible, according to the mean composite score of -1.066, which is marginally lower than Rampurhat-I. The interrelationship of vegetation coverage, elevation, and water retention in forming erosion patterns was highlighted by the significant effects of NDVI, slope, and TWI on these scores. To stabilize the most susceptible places and preserve overall resistance against gully erosion, these findings highlight the necessity of specialized management techniques. The Rampurhat-I block’s detailed composite factor scores are shown in Table 9, which sheds light on the regional variability of gully erosion susceptibility. Seven principal components are used in the analysis to get composite scores that represent the interplay of several morphometric and environmental elements. The northwestern part of the block, at latitude 24.188, and longitude 87.698, had the highest composite score (5.38) of all the extracted locations. Because of the area’s steep slopes, high rainfall erosivity, and lack of vegetation, there is a lot of runoff and soil erosion. This hotspot demonstrates how erosion processes are made worse by hydrological factors and compounded terrain imperfections. On the other hand, latitude 24.240, and longitude 87.780, which is in a flat region with higher vegetation cover and less intense runoff, had the lowest composite score (-5.41). With a mean composite score of -0.53 for the entire block, the area is somewhat susceptible to gully erosion overall. The spatial variability is noteworthy, though. Zones with positive composite scores (above the mean) make up about 40% of the block (calculated from clustering of high scores), designating them as extremely vulnerable areas in need of urgent erosion control measures. However, over 60% of the blocks are in zones with negative scores, which are below the mean and indicate regions with comparatively reduced risks of erosion. The Topographic Wetness Index (TWI), Stream Power Index (SPI), and Normalized Difference Vegetation Index (NDVI) are important variables that affect these values. TWI highlights the uneven topography and water retention patterns due to its significant variability (SD = 242.94). In regions with greater slopes and less vegetation, the energy of runoff flow is captured by SPI and is directly correlated with the severity of erosion. The moderately variable NDVI shows how vegetation stabilizes the soil and lowers the risk of erosion. The interplay of these elements is highlighted by the composite ratings in Table 9, which focus attention on the block’s most vulnerable areas for focused actions such as afforestation and slope stabilization. The Bolpur-Sriniketan block’s composite factor scores are shown in Table 10, which highlights the regional variability in the susceptibility of gully erosion. Latitude 23.633, and longitude 87.722, which is in the southeast part of the block, have the highest composite scores (5.308). This hotspot is crucial for erosion control because of its steep slopes, little vegetation, and high concentration of runoff. On the other hand, the location with the lowest composite score (-7.464) was latitude 23.569, and longitude 87.831, which indicates a flatter topography with more dense vegetation, hence decreasing the sensitivity to soil erosion. Compared to Rampurhat-I, the block’s mean composite score was -1.066, indicating that Bolpur-Sriniketan has a comparatively higher percentage of low-susceptibility zones. However, around 35% of the block is made up of zones with positive composite scores (above the mean), which are mostly found in the southeast and northwest. These places are in line with those that were found to be extremely vulnerable in earlier evaluations because of their high hydrological energy, uneven terrain, and inadequate vegetative cover. With negative ratings, the remaining 65% of the block is classified as a reduced susceptibility zone, emphasizing areas with flatter slopes and higher vegetation densities. The Topographic Wetness Index (TWI), slope, and NDVI are major factors that affect these results. The NDVI emphasizes the vital function that vegetation plays in preventing erosion because of its moderate variance (mean = 0.1082; SD = 0.147). Because vegetation absorbs surface water and stabilizes the soil, areas with higher NDVI are often less susceptible. Erosion processes are significantly impacted by TWI, which shows significant fluctuation in water retention and flow concentration with a standard deviation of 505.89. Another important factor is slope since steeper terrain has higher runoff velocities, exacerbating gully erosion. The composite factor scores in Table 10 indicate the regional extent of erosion vulnerability, which also clearly shows the priority regions for interventions like check dam construction, contour farming, and reforestation. The trends in Figs 12 and 13 demonstrate how these scores support the need for localized management strategies in high-susceptibility zones while maintaining the stability of lower-risk areas. According to Figs 12 and 13, the geographical distribution of gully erosion susceptibility in the Rampurhat-I and Bolpur- Sriniketan blocks categorize regions into five classes: very high, high, moderate, low, and very low. The majority of the extremely sensitive areas are in the southwest and northwest, making up about 23.45% of Rampurhat-I’s total size (68.97 square km). These areas are characterized by high runoff intensity, scant vegetation, and steep slopes, all of which make soil erosion worse. The high SPI and TWI scores of Rampurhat-I’s northwest correspond to areas of extremely high susceptibility. However, over 20.59% of the Bolpur-Sriniketan block’s land (69.81 square kilometers) is categorized as extremely vulnerable. The hotspots have been identified in the block’s southeast and northwest, where erosion risks are increased by concentrated runoff and less vegetation. Particularly susceptible to erosion are regions with steep slopes, limited vegetation cover, and high rainfall erosivity, as indicated by the regional variability in susceptibility in both blocks. This spatial study highlights the necessity of location-specific interventions to reduce gully erosion, such as slope stabilization in areas with steep gradients and afforestation in areas with scant vegetation.

With an average ROC-AUC value of 0.97 for Rampurhat-I, the model performs exceptionally well, demonstrating its capacity to reliably differentiate between the two classes in 97% of cases. This high score indicates that the model exhibits nearly complete separability between positive and negative classes, indicating that the predictor variable, Mean Composite Value, is very successful for classification. Practically speaking, this type of model is dependable for forecasting since it demonstrates outstanding sensitivity (high true positive rate) and specificity (low false positive rate) across a range of thresholds. This is a robust result because the model’s discriminative power increases with the AUC’s proximity to 1. A model’s capacity to differentiate between two classes at different threshold levels is graphically represented by the Bolpur- Sriniketan Receiver Operating Characteristic (ROC) curve. The model’s ROC-AUC (Area Under the Curve) value in this analysis is 0.975, which shows that it can effectively distinguish between the positive and negative classes. Positive occurrences are routinely given higher scores by the model than negative ones, according to an AUC value near 1. For binary classification tasks, this result demonstrates the model’s great reliability by achieving a high true positive rate while retaining a low false positive rate across thresholds. There is little overlap in the prediction scores for the two classes, indicating high predictive performance, as indicated by the small departure from a perfect score of 1.0.

Correlation Analysis of Gully Erosion Indicators

The interdependence of important parameters influencing gully erosion is examined in the correlation analysis for the Rampurhat-I block (Table 11). Rainfall erosivity (R) and drainage density (DI) showed a strong negative association (-0.363, p < 0.01), indicating that locations with larger drainage densities are less vulnerable to rapid runoff as erosivity rises, potentially because of improved water dispersion. The Modified Normalized Difference Water Index (MNDWI, -0.821, p < 0.01) and the Bare Soil Index (BSI, -0.620, p < 0.01) showed a substantial negative connection with NDVI, indicating that denser vegetation reduces bare soil exposure and maintains soil moisture, hence decreasing erosion vulnerability. Similar to this, the relationship between slope (S) and stream power index (SPI, 0.310, p < 0.05) showed that steeper terrains typically have higher erosive energy from water flow, which leads to the construction of more noticeable gullies. Stream Power Index and Topographic Wetness Index (TWI) had a negative correlation (-0.305, p < 0.05), suggesting that erosion dynamics may be affected by wetter areas’ lower flow velocity. The interplay of these elements highlights the vital roles that flora, hydrology, and terrain play in this block’s erosional processes. Significant correlations between the erosion indicators are also highlighted by the correlation analysis for the Bolpur- Sriniketan block (Table 12). Higher rainfall intensities lead to concentrated runoff rather than widely dispersed flow via existing drainage networks, as evidenced by the strong negative connection between rainfall erosivity (R) and drainage density (DI, -0.781, p < 0.01). The stabilizing effect of vegetation on soil and water retention was demonstrated by the significant negative correlation of the Modified Normalized Difference Water Index (MNDWI, -0.807, p < 0.01) and the significant positive correlation of the NDVI with the Bare Soil Index (BSI, 0.835, p < 0.01). Furthermore, there was a positive association between slope (S) and the Modified Normalized Difference Vegetation Index (0.471, p < 0.01), indicating that places with moderate slopes are likely to have more vegetative cover, which helps reduce erosion. The strong correlation between drainage parameters like drainage density (DD) and infiltration (IF) (0.922, p < 0.01) highlights how they work together to control surface water flow and lower the danger of erosion. Additionally, there was a negative correlation between rainfall erosivity and the Topographic Wetness Index (TWI) (-0.579, p < 0.01), suggesting that wetter regions tend to mitigate the immediate erosive effects of heavy rainfall. According to this analysis, Bolpur- Sriniketan’s hydrological, topographical, and vegetation-related interactions are intricate and multifaceted, necessitating the use of integrated management techniques to reduce erosion. The interrelationships between 20 important parameters influencing gully erosion in Rampurhat-I block are examined using the correlation analysis in Table 11. These relationships demonstrate how morphometric, hydrological, and environmental factors interact to affect soil erosion processes. Rainfall Erosivity (R) and Drainage Density (DI) showed a significant negative correlation (-0.363, p < 0.01), indicating that areas with a dense drainage network are less likely to experience concentrated runoff as rainfall intensity increases because of improved water dispersion throughout the landscape. On the other hand, a positive association between Slope (S) and Stream Power Index (SPI, 0.310, p < 0.05) highlights how steeper topography boosts water flows’ erosive power. This is crucial because steep slopes make up about 25% to 30% of Rampurhat-I, which raises SPI values and causes substantial soil erosion. Strong negative relationships were found between the vegetation indices, including the Normalized Difference Vegetation Index (NDVI), and markers of soil exposure and runoff, including the Modified Normalized Difference Water Index (MNDWI, -0.821, p < 0.01) and the Bare Soil Index (BSI, – 0.620, p < 0.01). These results demonstrate that regions with greater vegetation densities have better moisture retention and less soil exposure, which reduces their vulnerability to erosion. The protective function of vegetation, which stabilizes soil and lessens gully formation, is seen in the link between NDVI and BSI. These correlations highlight the fact that between 50% and 60% of regions with low NDVI values are more vulnerable to erosion. The Stream Power Index (SPI, -0.305, p < 0.05) showed a negative connection with the Topographic Wetness Index (TWI), indicating that wetter regions tend to slow down water flow, lowering its erosive potential. The variation in erosion intensity throughout Rampurhat-I can be explained by this, as flatter areas with higher TWI values are less vulnerable. DD demonstrated high positive associations with both flow length (Fl, 0.364, p < 0.01) and infiltration variables (IF, 0.926, p < 0.01), suggesting that regions with well- distributed drainage patterns effectively regulate runoff, reducing the susceptibility to erosion. Table 11 shows that the interplay of rainfall intensity, slope, plant cover, and drainage efficiency have a major impact on gully erosion in Rampurhat-I. Roughly 35–40% of the block is extremely prone to erosion, with steep slopes and limited vegetation creating localized hotspots. To reduce soil loss, this emphasizes the necessity of focused interventions such as afforestation and drainage network upgrades. The correlation matrix for the Bolpur- Sriniketan block is shown in Table 12, which clarifies the intricate connections between the 20 gully erosion indicators. Rainfall Erosivity (R) and Drainage Density (DI) have a substantial negative connection (-0.781, p < 0.01), indicating that regions with higher rainfall intensities are more likely to have concentrated runoff because of limited drainage dispersion. This demonstrates how 30% to 40% of the block is susceptible to significant gully erosion, especially in areas with inadequate drainage systems. Vegetation indices were important; the NDVI had a negative connection with the Modified Normalized Difference Water Index (MNDWI, -0.807, p < 0.01) and a substantial positive correlation with the Bare Soil Index (BSI, 0.835, p < 0.01). This suggests that densely vegetated areas have more stable soil and higher moisture retention, both of which help prevent erosion. However, the absence of vegetation in low-NDVI areas—which make up roughly 45% to 50% of the block—exacerbates soil exposure and runoff, making the area more vulnerable to erosion. Wetter regions are less impacted by rainfall-induced erosion because water retention lessens the effect of runoff, according to the Topographic Wetness Index (TWI), which showed a substantial negative connection with Rainfall Erosivity (R, -0.579, p < 0.01). This result is consistent with findings showing lesser erosion hazards are experienced by 25%–30% of the block’s terrain, which is defined by high TWI values. Furthermore, a strong association was found between slope (S) and hydrological measures such as the Stream strength Index (SPI, 0.471, p < 0.01), highlighting the erosive strength of water flows being exacerbated by steep gradients. The crucial role drainage networks play in controlling runoff and minimizing erosion is demonstrated by the substantial interdependencies (0.922, p < 0.01) and positive correlations with flow length (Fl, 0.349, p < 0.05) and slope of drainage parameters, such as Drainage Density (DD) and Infiltration (IF). Nonetheless, 20% to 25% of the block is still made up of sections with inadequate drainage distribution, which makes them extremely vulnerable. In conclusion, Table 12 shows that the interaction of rainfall intensity, vegetation cover, slope, and drainage efficiency mostly affect gully erosion in Bolpur-Sriniketan. The remaining portions of the block show resilience because of efficient drainage patterns and superior vegetation covering, although 35%–40% of the block is extremely vulnerable to erosion because of inadequate drainage and sparse vegetation. As part of erosion mitigation methods, these findings highlight the need to increase vegetation density, strengthen drainage networks, and lower slope instability. Rainfall and Drainage: There is a strong negative correlation between rainfall erosivity (R) and drainage density (DI) in both blocks; however, the relationship is slightly stronger in Bolpur-Sriniketan (-0.781) than in Rampurhat-I (-0.363), suggesting that Bolpur-Sriniketan is more susceptible to runoff. plant and Erosion: Bolpur-Sriniketan had stronger NDVI associations with BSI and MNDWI (0.835, -0.807) than Rampurhat-I (-0.620, -0.821), indicating more spatial variability in plant cover and its protective roles. Hydrology and Slope: Both blocks exhibit a considerable link between SPI and slope, highlighting their combined influence on erosion intensity. This link demonstrates how steep terrain in 35%–40% of Bolpur-Sriniketan and Rampurhat-I is vulnerable. These studies highlight the complexity of gully erosion and the requirement for integrated management strategies catered to the unique features of every block.

Covariance among Gully Erosion Indicators (SEM)

The covariance matrix for the main indicators of gully erosion in the Rampurhat-I block is shown in Table 13, providing information on the intricate interactions between several hydrological, topographical, and environmental factors. Rainfall Erosivity (R) and Drainage Density (DI) have a substantial negative correlation of -0.669, indicating that areas with high rainfall erosivity which are defined by heavy rainfall events frequently have low drainage densities. Gully formation is made worse by concentrated surface runoff brought on by ineffective drainage. Erosion risk is raised in these places because of inadequate drainage networks, which make it difficult for water that falls during intense rain to disperse. Conversely, regions with a strong drainage system aid in directing surplus water and lessening the effects of rainy erosivity. The correlation of 0.310 between the Stream Power Index (SPI) and the Topographic Wetness Index (TWI) is another important discovery. Areas with higher moisture retention (higher TWI values) typically have higher water flow energy (higher SPI), which makes them more vulnerable to gully formation during periods of heavy rainfall. TWI is a measure of water accumulation and drainage. The protective function of vegetation against erosion is highlighted by the positive correlation (0.219) between the Normalized Difference Vegetation Index (NDVI) and the Bare Soil Index (BSI). In Rampurhat-I’s gully-prone sections, high plant coverage stabilizes the soil and helps stop soil erosion by limiting the amount of bare soil exposed to water flow. Additionally, the cooling impact of vegetation on the land surface is shown by the negative covariance of -0.259 between the NDVI and LST. Densely vegetated areas usually have cooler surface temperatures, which helps retain moisture and lowers evaporation. By preventing the soil from drying out, this moisture retention enhances soil structure and lessens erosion vulnerability. greater drainage density areas also have greater rates of water infiltration, which means that more water is absorbed by the soil rather than flowing over the surface as runoff, according to the covariance of -0.481 between Drainage Density (DI) and Infiltration Factor (IF). This procedure lessens erosion by stopping the water from moving quickly, which would otherwise damage the soil. These covariance correlations show that Rampurhat-I contains areas that are particularly vulnerable to erosion due to a combination of circumstances such as heavy rainfall, steep slopes, and a lack of vegetation. Targeted erosion control measures are necessary in these susceptible locations, as approximately 40–50% of the block is in areas with increased erosion risks because of these interconnected causes. The covariance study for the Bolpur-Sriniketan block is shown in Table 14, which also highlights the intricate relationships between the main parameters affecting gully erosion in this area. The significant negative correlation of -0.781 between Drainage Density (DI) and Rainfall Erosivity (R) is a noteworthy finding that emphasizes a crucial link. In regions with heavy rainfall, the concentration of runoff is made worse by inadequate drainage systems, which speed up erosion. Because water cannot effectively spread out and penetrate the soil in these regions, there is a high surface runoff that separates soil particles and worsens gully formations, making these areas vulnerable to severe gully erosion. This correlation emphasizes how crucial it is to improve drainage infrastructure in areas with heavy rainfall to manage runoff and lessen the likelihood of erosion. Slope (S) and Stream Power Index (SPI) have another important covariance of 0.254, which shows that steeper terrain often has more erosive water flow because of the higher runoff velocity. The force of running water is increased on steep slopes, increasing the likelihood of gully development and soil erosion. This relationship emphasizes the necessity of slope stabilization measures in sensitive regions of the block, as locations with higher slope values are more prone to these erosive forces. Gully erosion is more likely to occur in areas with low plant cover (and hence high BSI values), according to the covariance between the Normalized Difference plant Index (NDVI) and the Bare Soil Index (BSI) of 0.232. Without the stabilizing influence of plants, the bare soil in these places is more vulnerable to rainfall, making it easier for runoff to separate the soil. Another important discovery is that the Topographic Wetness Index (TWI) and Stream Power Index (SPI) have an inverse connection (- 0.343). Greater water retention capacity is indicated by higher TWI values, which would typically slow down surface runoff’s pace and erosive power. However, in regions with high TWI, the buildup of water may cause the soil to become more saturated, which increases the soil’s susceptibility to erosion when significant amounts of water are released. Ultimately, the covariance of -0.433 between the Infiltration Factor (IF) and Drainage Density (DI) indicates that areas with higher drainage densities are better equipped to control water flow, enabling more water infiltration and lessening the severity of surface runoff. Therefore, these places are less susceptible to erosion than those with inadequate drainage systems. About 35% to 45% of Bolpur-Sriniketan is extremely susceptible to gully erosion, according to the covariance correlations in Table 14, especially in regions with steep slopes, little vegetation, and inadequate drainage. Reducing erosion hazards in the block requires addressing these problems through better slope management, drainage infrastructure, and vegetation cover.

A solid mathematical foundation for comprehending how morphometric, hydrological, and environmental elements contribute to gully erosion in the area is provided by the Structural Equation Modeling (SEM) analysis for the Rampurhat-I block, as shown in Fig 14. To determine the factors that contribute to and hinder erosion, SEM is very helpful in assessing the direct and indirect correlations between important variables. The Topographic Wetness Index’s (TWI) high standard deviation (SD = 242.94) highlights the landscape’s variability in moisture retention and flow accumulation, making it one of the SEM model’s most significant findings. The variability of TWI emphasizes how important it is in regulating runoff dynamics and soil moisture, both of which have a direct impact on the formation of gullies. Greater water accumulation in areas with higher TWI values is likely to cause soil saturation and increase erosion vulnerability. Another important component of the SEM model is the Stream Power Index (SPI), which measures the energy of runoff. SPI has a strong correlation with both rainfall erosivity and moisture retention (TWI), indicating that regions with high rainfall erosivity and TWI have more concentrated and powerful water flows, which speeds up the formation of gullies. Furthermore, the SEM model revealed a substantial correlation between the Normalized Difference Vegetation Index (NDVI) and the Bare Soil Index (BSI), indicating the protective function of vegetation in halting soil erosion. More vegetation cover is generally found in areas with higher NDVI values, which stabilizes the soil and lessens its vulnerability to erosion by shielding it from the direct effects of rainfall. According to the model, by promoting water infiltration into the soil and lowering the volume and velocity of surface runoff, Drainage Density (DD) and Infiltration Factor (IF), which had smaller standard deviations (SD= 0.0534 and 0.0874, respectively), consistently help to mitigate runoff. The SEM also highlights the indirect impact of terrain characteristics like elevation and slope on erosion, mainly through their influence on hydrological elements like drainage efficiency and runoff. To put it another way, slope, and elevation affect how water moves across the landscape and how vulnerable some locations are to erosion, even though they do not directly cause erosion. The complicated nature of gully erosion in Rampurhat-I is explained by the combination of these components, where successful erosion control requires consideration of both hydrological circumstances (such as water retention and flow) and environmental variables (such as vegetation cover). As shown in Table 14, the SEM analysis for the Bolpur-Sriniketan block builds on concepts like those found in Rampurhat-I, but it also identifies some significant distinctions in the way morphometric and environmental factors interact in this block. Like Rampurhat-I, the Topographic Wetness Index (TWI) had the largest standard deviation (SD = 505.89) and was found to be a dominant factor in the SEM model. This implies that TWI has a similar function in Bolpur-Sriniketan by highlighting regions that are extremely vulnerable to erosion due to their high-water retention and runoff potential. Nonetheless, Bolpur-Sriniketan’s noticeably greater TWI variability indicates that the block’s moisture retention varies more drastically, which may lead to more specialized erosion hotspots based on landscape characteristics and rainfall intensity. In Bolpur-Sriniketan, erosion dynamics are also significantly impacted by the Slope and Stream Power Index (SPI), which have standard deviations of 1.28 and 0.43, respectively. The significance of topography in influencing the velocity and erosive potential of water flows is highlighted by the covariance between these factors. Steeper slopes lead to higher SPI values, which raise water’s erosive energy and exacerbate gully formation and soil detachment. Strong negative associations between the Normalized Difference Vegetation Index (NDVI) and the Bare Soil Index (BSI) were discovered by the SEM model in terms of vegetation. This suggests that places with dense vegetation, or high NDVI values, assist limit bare soil exposure and hence avoid erosion. This connection emphasizes how crucial vegetation is for stabilizing soil and lessening its vulnerability to erosion brought on by water. The model also emphasizes the relationship between Infiltration Factor (IF) and Drainage Density (DI), which were found to be highly associated. This suggests that regions with advanced drainage systems are better equipped to manage surface runoff, which lowers the risk of erosion. Gully development is largely influenced by concentrated flow, which is reduced when water is dispersed throughout the landscape with a high drainage density. The SEM findings for Bolpur-Sriniketan show that a combination of vegetation cover, moisture retention, and topography features significantly affects the block’s vulnerability to gully erosion. This block’s localized erosion highlights the necessity of specialized management approaches that take into consideration particular hydrological circumstances and vegetation distribution.

The accuracy of the sample means in determining the population mean for different indicators is shown by the standard error of the mean (SEM) for the Rampurhat-I and Bolpur-Sriniketan blocks. Rainfall (R) has a higher SEM of 1.02322, indicating significant variability in rainfall patterns throughout the block, while SEM values in Rampurhat-I vary greatly among indicators. The region’s large variety of wetness conditions, which may influence the susceptibility of gullies to erosion, is also indicated by the Topographic Wetness Index (TWI), which shows an even greater SEM of 32.75845. In comparison to Rampurhat-I, Bolpur-Sriniketan displays comparatively lower SEM values, with Slope (S) at 0.17708 and Rainfall (R) at 1.15830, indicating less variation in slope and rainfall. Although Bolpur-Sriniketan shows more homogeneity in these areas, Rampurhat-I’s higher SEM values generally indicate greater uncertainty in the mean values for several variables, reflecting the more varied nature of the landscape and climate.

Gully Hotspots and Their Impact on Land Use Land Cover

This region is highly influenced by weathering, rill, and gully formation as well as soil loss and soil degradation. Gully erosion hotspot areas are identified from the two blocks shown in Fig 16. Gully erosion has severely escalated the degradation of agricultural land through an increase in surface runoff and a decrease in the groundwater level. The gully-affected area is about 31.01 square km in Rampurhat-I block. Agricultural land was 165.54 square km in 2011 (BAES, 2011), which has been reduced to 128.44 square km at present. Contrarily Rampurhat-I block belongs to a very low productivity zone, whereas Bolpur-Sriniketan block is under the high productivity zone of this district (Saha and Rudra 2019). Observable gully-affected badland areas (locally noted as Khoai) are 26.98 square km, originating in Bolpur-Sriniketan block. As per the impact of gully erosion, many agricultural workers have changed their occupation from agriculture to stone crackers, construction workers, and small retailers. The middle, southern, and south-western portions of Rampurhat-I block are less populated as they belong to the high gully susceptibility zone, compared to the rest of the areas. Similarly, the southeastern portion of Bolpur-Sriniketan block comes under the highly susceptible zone, and therefore, it is less populated. Gully erosion has gradually decreased the forest coverage area of both C. D. blocks. Besides, Chhora, Nrayanpur, Tejhati, and Hazarpur of Rampurhat–I block are highly deforested; the central and southern parts (Kurumgram, Kluha, Ambhaetc) of this block have been covered with part of the vegetation patches. On the other hand, Ballavpur, Bonerpukur Danga, and Golbarietc are densely covered by forest areas. The eastern portion of this block is also deforested due to the impact of gully erosion and head- cut channel formation (Fig 17). The present study shows that the existing area under gully erosion is greater in Rampurhat-I than Bolpur-Sriniketan block, and it also impacts on the average landholding sizes in the two

C.D. blocks. The Rampurhat-I consisted of 1.05 hectares of average landholding sizes, and Bolpur-Sriniketan consisted of 1.13 hectares of average landholding sizes in 2001 (Agricultural Census 2001 of India), which was decreased to 0.87 hectares in Rampurhat-I and 1.07 hectares in Bolpur-Sriniketan block. As per participants’ observation, the drained-out soil and sediments transported through the degraded gully channels and deposited in the extended erosional land and thus impacting the neighboring areas during and after severe rainfall and deteriorating the fertility of the soil and agricultural lands. Fig 18 shows the situation of the cross profile of the gully channels along with the gully erosion land cover areas in Rampurhat-I and Bolpur- Sriniketan, respectively. Identification of Hotspots: Rampurhat-I and Bolpur-Sriniketan both have sizable hotspots for gully erosion, with the highest scores of 5.38 and 5.308, respectively. These sites are primarily found in areas with significant runoff, little vegetation, and steep slopes. Composite Mean Scores: Bolpur- Sriniketan’s mean composite score is -1.066, while Rampurhat-I’s is -0.53, suggesting that Rampurhat-I is comparatively more susceptible. Areas of High Susceptibility: About 40% of Rampurhat-I and 35% of Bolpur- Sriniketan are high-susceptibility zones, which call for focused interventions. Altogether, these findings provide a strong framework for organizing erosion control tactics, highlighting the necessity of flexible approaches catered to the unique patterns of susceptibility in every block. Gully erosion-related changes to land use and land cover (LULC) in both blocks show considerable losses in forest cover and substantial deterioration in agricultural fields. The gully-affected regions in Rampurhat-I cover 31.01 square km, which means that agricultural land has decreased by 22.4% from 165.54 square km in 2011 to 128.44 square km now. The rising gully erosion, which reduces soil fertility and renders vast tracts unusable for cultivation, is the direct cause of this decline in agricultural land. Similar to this, gully-affected areas in Bolpur-Sriniketan span 26.98 square kilometers, resulting in noticeable deforestation, especially in the block’s central and southern regions. Hotspots like Ballavpur in Bolpur-Sriniketan and Chhora and Nrayanpur in Rampurhat-I have experienced severe vegetation loss as a result of erosion. The wider socioeconomic effects of erosion are shown in the change in land usage from agriculture to alternate forms of income like stone-cracking and construction. Landholding sizes decreased from 1.05 hectares in 2001 to 0.87 hectares in Rampurhat-I and from 1.13 hectares to 1.07 hectares in Bolpur-Sriniketan, respectively. To prevent additional LULC alterations and enhance agricultural output, these modifications underscore the pressing need for integrated land management methods that prioritize erosion prevention, sustainable land use practices, and the repair of degraded areas. Different levels of model accuracy are revealed by the Kappa statistic values for LULC accuracy evaluation and gully erosion susceptibility in Rampurhat-I and Bolpur-Sriniketan blocks. Strong agreement between the observed and predicted classifications is indicated by a higher Kappa value (e.g., 0.75 for Rampurhat-I), indicating that the gully erosion susceptibility model is highly accurate in this block and that there is a clear correlation between erosion-prone areas and changes in land use, such as deforestation and agricultural land loss. A lower Kappa value (such as 0.60) for Bolpur-Sriniketan, on the other hand, indicates moderate agreement and suggests that the model’s predictions for gully erosion and the LULC changes that occur from it are less accurate. This lower Kappa indicates greater uncertainty in the Bolpur-Sriniketan erosion assessment, requiring additional model refining to increase accuracy. The need for region-specific approaches to regulating gully erosion and reducing its effects on land use is highlighted by these disparate Kappa values.

Major Findings

In summary, Gully erosion occurs mainly when the surface runoff concentrates strongly within a particular channel which detaches the soil particle and makes its channel. This type of erosion is very common in Rampurhat-I and Bolpur-Sriniketan blocks. Gully-susceptible zones are identified in these two blocks based on erosion vulnerability. The southwestern and north-western parts of the Rampurhat–I block are highly susceptible to gully erosion. This consists of a 68.97 square km area of this block. That means 23.45% of this block is badly affected by gully erosion. In the case of Bolpur-Sriniketan block, this value is 20.59% and the area coverage is 69.81 square km area. It is a clear indication that gully erosion in both blocks is going to be a matter of concern because this area’s coverage is increasing day by day. In the study, a total of twenty factors have been selected to formulate a composite gully erosion susceptibility zone in the case of both Rampurhat- I and Bolpur-Sriniketan. A total of seven principal components are extracted which cumulatively justified the gully erosion with 77% in Rampurhat-I and 74% in Bolpur-Sriniketan blocks. The highest deviation (SD = 242.94) of the factor topographic wetness index and lowest deviation (SD = 0.0534) of the factor drainage density are being followed; whereas in Bolpur-Sriniketan the highest deviation (SD = 505.89) is of the factor topographic wetness index and lowest deviation (SD = 0.054) is of the factor flow length. The highest standard error of the mean has been followed in the case of topographic wetness index in both of the areas, and the lowest standard error of the mean has been followed in the case of drainage density in Rampurhat-I and flow length in Bolpur-Sriniketan. The highest standard deviation values of the topographic wetness index depict that the impact of topographic elements on the hydrological elements of gully-originated river basin areas of the two blocks is variable, and, in this condition, the differential erosional process has occurred in those areas. Here, drainage density and flow length are mostly consistent in affecting the gully formation and sediment transportation in the study areas. They highlighted the bivariate trends of the relationship of rainfall erosivity with its predictors in Rampurhat-I and Bolpur Sriniketan, respectively. The most influential factors are slope, Normalized Difference Vegetation Index (NDVI), Modified Normalized Difference Water Index (MNDWI), Normalized Difference Buit-up Index (NDBI), Bare Soil Index (BSI), and Land Surface Temperature (LST) analysis in the scenarios of both Rampurhat-I and Bolpur-Sriniketan. The spatial variation of the selected factors which impact gully erosion formation and processes in the study area in Rampurhat-I, and Bolpur- Sriniketan. The composite gully erosion susceptibility indices show that very high (0.833-1.230), high (0.431- 0.832), moderate (0.0298-0.431), low (-0.370 to -0.0297), and very low (-0.772 to -0371) gully erosion susceptibility areas are found in the north-west; north-west; north-west, west, middle and part of north, south, and east; mainly north, east and part of the south and north-east respectively in Rampurhat-I block. While, the areas of very high (0.393-0.756), high (0.0283-0.392), moderate (-0.335 to -0.0282), low (-0.669 to -0.336), and very low (-1.060 to -0.700) gully erosion susceptibility are found in the part of the south-east; south-east and north-west; part of the south-east, middle and north-east, parts of north-west and rarely in the north-east and middle-east portion of Bolpur-Sriniketan block respectively. The figures show the predicted zones of gully erosion susceptibility, respectively, in the two community development blocks. The structural equation model culminates in the correlation among the influential factors of the gully erosion process, respectively, in the two blocks. The broader categories are terrain analysis factors, drainage basin morphometry, and raster analysis, which assemble the composite gully erosion susceptibility indices. The value of covariances among the factors and statistical significances have been mentioned in a tabular format of Rampurhat-I and Bolpur- Sriniketan blocks, respectively. The overlapping map of gully erosion susceptibility isolines on rainfall erosivity index maps shows that higher rainfall erosivity creates high soil erosion in the gully susceptibility areas.

DISCUSSION

Indicators of Gully Erosion Susceptibility

The Rampurhat-I and Bolpur-Sriniketan blocks differ greatly in their susceptibility to gully erosion because of variations in terrain, hydrology, and land use. The importance of the Topographic Wetness Index (TWI) as a crucial governing factor can be explained by examining the correlation of indicators to determine the dominating components in each block. Studies by Rahmati et al. (2017) and Wang et al. (2019) found that high STI values are strong drivers of erosion because they can quantify sediment transport capacity. In Rampurhat-I, the STI is the best predictor of gully erosion susceptibility (r = 0.61). Demir and Kisi (2016) have highlighted the importance of Drainage Density (r = 0.51) and Average Slope (r = 0.51) as additional crucial indicators in their study of runoff-driven erosion. By efficiently channeling water, a denser drainage system produces concentrated flow routes that increase erosive force. As Pal et al. (2020) also point out, steeper slopes increase the danger of erosion by decreasing infiltration and speeding up soil separation. The Normalized Difference Built-up Index (NDBI) shows the effects of urbanization, with less vegetation increasing runoff (r = 0.47). According to Zhang et al. (2021), impermeable surfaces in metropolitan areas make them more susceptible to erosion. This finding is consistent with their findings. The rough topography, which directly controls runoff and sediment transport, explains why TWI has a weak association (r = 0.057) in this case. In the Ethiopian highlands, Gessesse et al. (2016) reported similar findings, showing that in steep terrain, slope overshadowed TWI. The dominant factors are different in Bolpur-Sriniketan. The most important indicator is the Infiltration Number (r = 0.48), which confirms the findings of Tebebu et al. (2010), who showed that decreased infiltration increases the risk of erosion by accelerating runoff. Chen et al. (2017) found that vegetation stabilizes soil and lessens the impact of raindrops. This protective role is supported by the NDVI (r = 0.45). Water channels are important, as indicated by hydrological indicators such as Flow Direction (r = 0.43) and Stream Frequency (r = 0.40). Rahimi et al. (2020) discovered that stream frequency played a significant role in sediment mobilization. TWI has growing significance in areas with less steep topography, as evidenced by its moderate correlation (r = 0.16). According to studies like Sørensen et al. (2006) and Beven and Kirkby (1979), TWI directly measures soil saturation and water buildup, both of which are essential for gully start. Because the mild slopes in Bolpur-Sriniketan let water concentrate rather than immediately drain off, they strengthen the role of TWI.

The Topographic Wetness Index (TWI) calculates the likelihood of water accumulation by combining drainage area and slope. Because saturated soils are more likely to detach, TWI is a universal indicator of erosion vulnerability. TWI is essential in explaining erosion because it connects topographic and hydrological processes, according to studies by Rahmati et al. (2017) and Ali et al. (2022). The role of TWI is less evident in steep regions like Rampurhat-I because water storage is overshadowed by quick runoff. However, TWI becomes more prevalent in moderate terrains like Bolpur-Sriniketan because of higher soil saturation and slower drainage. This dual behavior supports the wider conclusions of Beven and Kirkby (1979) and emphasizes TWI’s versatility as a crucial framework for simulating gully erosion susceptibility across a range of environments.

Rugged Terrain and Active Fluvial Erosion

The majority of the granite-gneiss and metamorphic rock formations that make up Rampurhat-I’s geological setting are susceptible to chemical weathering in tropical climates. Because of their weathering, these minerals create loose sediments that are very transportable. Steep slopes and small valleys in the rough terrain increase runoff velocity, which makes the environment perfect for severe gully erosion. The Stream Transport Index (STI), which has the highest influence (r = 0.61), emphasizes how sediment transport capability contributes to erosion. This behavior is also seen in the Ethiopian Highlands, where steep slopes and concentrated flow greatly increase sediment mobility (Gessesse et al., 2016). A well-developed stream network’s erosive power is highlighted by drainage density (r = 0.51), where concentrated flows more efficiently erode soils. In India’s Chambal Ravines, comparable patterns have been documented (Mishra et al., 2015). The Himalayan foothills have a similar pattern to the average slope (r = 0.51), which promotes rapid runoff generation, decreases infiltration, and increases soil detachment (Sharma et al., 2018). Steep slopes and roughness, which more directly control erosion processes, are the main causes of the Topographic Wetness Index’s (TWI) weak association (r = 0.057). Similar results in the mountainous regions of northern India and Nepal highlight how slope-related factors frequently take precedence over TWI in these environments (Pandey et al., 2020). With rolling hills and vast lateritic soils, Bolpur-Sriniketan offers a striking geomorphological and geological backdrop. Laterites, which are high in iron and aluminum, are more prone to erosion because they lose structural strength when they become soaked. The significance of soil saturation indicators like TWI is increased by the more substantial water buildup made possible by the moderate slopes and lower drainage density. The major component is the Infiltration Number (r = 0.48), suggesting that low-infiltration locations are especially vulnerable. Lower permeability increases runoff and erosion in semi-arid areas of Maharashtra and the Western Ghats, which is consistent with this finding (Tebebu et al., 2010; Gadgil et al., 2016). The stabilizing and erosion-reducing properties of plant cover are highlighted by the Normalized Difference plant Index (NDVI) (r = 0.45). In China’s Loess Plateau, where forest cover considerably reduces the risk of erosion, similar patterns have been documented (Zhao et al., 2019). An important factor in this block is the Topographic Wetness Index (TWI), which has a moderate association (r = 0.16). Permeable soils and moderate slopes let water build up, which affects soil saturation levels and erosion processes. Similar dynamics are shown in studies conducted in the lateritic terrains of southern Africa and West Bengal, where TWI becomes crucial in flat-to-moderate terrain (Lal et al., 2021).

Geological and Geomorphological Context: Global and Indian Comparisons

TWI has been identified as a crucial component of gully erosion research worldwide, especially in North America (Wilson & Gallant, 2000) and Europe (Sørensen et al., 2006). However, as shown in the Ethiopian Highlands, slope and drainage network characteristics predominate in steep terrains such as Rampurhat-I (Gessesse et al., 2016). On the other hand, semi-arid areas of India are like mild terrains like Bolpur-Sriniketan, where vegetation and soil permeability affect erosion (Pandey et al., 2020). These results are corroborated by Indian research conducted in the Western Ghats (Gadgil et al., 2016) and the Chambal Ravines (Mishra et al., 2015). Regions with undulating topography place more emphasis on vegetation and soil saturation than steep areas do on slope and drainage density. As a crucial measure of gully erosion susceptibility, TWI’s dual role of responding to terrain-specific dynamics highlights its global relevance. An examination of the Rampurhat- I and Bolpur-Sriniketan blocks’ susceptibility to gully erosion reveals how important topographical, hydrological, vegetative, and land use/land cover (LULC) elements are in determining the erosion dynamics of these areas. Though Bolpur-Sriniketan has moderate susceptibility, which is influenced by factors like soil saturation and vegetation cover, Rampurhat-I exhibits higher overall susceptibility, especially in areas with steep slopes, sparse vegetation, and high rainfall intensity. Using sophisticated statistical and modeling methods, this study offers insightful information that can direct land management and erosion control plans specific to each region’s circumstances.

Topographic and Hydrological Factors

The key elements causing erosion in Rampurhat-I were found to be slope, drainage density (DD), and stream power index (SPI). Slope (S) and SPI showed strong positive associations (0.310), suggesting that steeper terrain increases surface runoff energy and causes more noticeable erosion. This result aligns with international research, such as Gessesse et al. (2016) in the Ethiopian Highlands and Mishra et al. (2015) in the Chambal Ravines, which found that high slopes were a major factor in the creation of gullies. Similar to this, Rampurhat-I’s Drainage Density (DD) demonstrated a substantial correlation with both Flow Length (Fl) and Infiltration Factor (IF), indicating that the effectiveness of a dense drainage network in distributing water can either improve or decrease runoff dynamics. These results are consistent with those of Pandey et al. (2020), who demonstrated that drainage networks had a major impact on soil erosion and runoff dispersion in comparable semi-arid areas. Topographic Wetness Index (TWI), despite its modest correlation in Rampurhat- I (r = 0.057), is a significant element in Bolpur-Sriniketan where the terrain is less steep. In Bolpur-Sriniketan, where mild slopes and comparatively flat terrain improve moisture retention, TWI showed a stronger link with erosion processes (r = 0.16). Several studies have emphasized the significance of TWI in places with moderate slopes, such as Beven and Kirkby (1979), who identified TWI as a critical element in situations where water buildup causes soil saturation, increasing erosion susceptibility. In a similar vein, Sørensen et al. (2006) observed that TWI has a major impact on soil detachment processes in regions with moderate topography by measuring water retention.

Vegetation and Soil Protection

A significant outcome of this study has been the importance of vegetation in reducing gully erosion. The Normalized Difference Vegetation Index (NDVI) in Rampurhat-I showed significant negative correlations with both the Modified Normalized Difference Water Index (MNDWI) (-0.821) and the Bare Soil Index (- 0.620). This suggests that denser vegetation lessens soil exposure and improves moisture retention, both of which lessen erosion susceptibility. This aligns with Zhao et al. (2019) in China’s Loess Plateau, where they discovered that vegetation significantly reduced soil erosion by improving water retention and lessening the effects of rainfall. Afforestation and vegetation restoration are necessary in places with low vegetation cover since Chen et al. (2017) pointed out that NDVI and soil erosion are inversely correlated. The NDVI-BSI relationship in Bolpur-Sriniketan (0.835) emphasizes the importance of vegetation in stabilizing the soil. The stabilizing effect of vegetation in forested areas, where high plant density considerably mitigated erosion, was also shown by the study conducted in the Western Ghats by Gadgil et al. (2016). However, in some parts of Bolpur-Sriniketan, the low NDVI values increase the risk of erosion since bare soil is more vulnerable to surface runoff and rainfall. Accordingly, Bolpur-Sriniketan’s sparsely vegetated areas need urgent measures to boost plant cover to reduce erosion.

Soil and Moisture Dynamics

The Bare Soil Index (BSI) and Infiltration Factor (IF) are crucial measures of soil exposure and permeability. The negative correlation between Infiltration (If) and Drainage Density (DD) in Rampurhat-I (-0.481) indicates that better water infiltration occurs in locations with higher drainage densities, reducing runoff and soil erosion. In support of this, Tebebu et al. (2010) discovered that in semi-arid areas, high drainage density enhances water retention and infiltration while lowering runoff and the erosive potential of rainfall. The opposite is true in Bolpur-Sriniketan, where runoff is concentrated in low drainage density areas, which speeds up erosion. To lessen the negative effects of heavy rainfall, drainage systems in high-risk locations must be improved. In Rampurhat-I, the Topographic Wetness Index (TWI) and SPI had a negative correlation (-0.305), suggesting that wetter regions, which retain more water, typically have lower flow velocities, which lowers the water’s erosive capacity. However, if the soil is too saturated during periods of heavy rainfall, regions with high TWI may still have an increase in erosion. This result is consistent with research by Pandey et al. (2020), who found that if the moisture content in a certain topography is more than what the soil can absorb, water retention may make the soil more vulnerable.

Spatial Distribution of Erosion Susceptibility

According to the spatial distribution of gully erosion susceptibility in both blocks, Bolpur-Sriniketan has 20.59% of its area classed as extremely sensitive, while Rampurhat-I has roughly 23.45%. Steep slopes, little vegetation, and intense runoff are characteristics of these regions that make soil erosion worse. These results are consistent with the findings of Gessesse et al. (2016) and Mishra et al. (2015), who discovered that areas with comparable topography and hydrology were extremely vulnerable to gully erosion. The identification of hotspots in both blocks calls for specific erosion management strategies, including plant restoration, slope stabilization, and drainage network enhancement. The changes in land use and land cover (LULC) that have been noted in both blocks are in line with overall patterns in India and international research, especially the decline in agricultural land brought on by gully erosion. Bolpur-Sriniketan and Rampurhat-I both saw a 22.4% decrease in agricultural land, which was indicative of a move away from agriculture and toward non- agricultural occupations like construction and stone-cracking. The findings from the Chambal Ravines (Mishra et al., 2015) and Western Ghats (Gadgil et al., 2016), where traditional farming methods have been abandoned by rural people due to erosion-induced soil degradation, are in line with these changes. This emphasizes how erosion has a socioeconomic influence, causing changes in local livelihoods in addition to deteriorating land quality. This study highlights the complex relationship between gully erosion susceptibility in Bolpur-Sriniketan and Rampurhat-I. The pattern of Bolpur-Sriniketan is a little more complicated, influenced by moderate slopes, increased soil saturation, and variable drainage networks, whereas Rampurhat- I is more vulnerable to erosion because of its steep slopes and sparse vegetation. The findings support the necessity of integrated land management plans that take into consideration the dynamics of the local vegetation, hydrology, and topography. Restoring degraded lands and reducing soil loss in these areas depend heavily on erosion control techniques such as afforestation, enhanced drainage systems, slope stability, and soil conservation. Constructing adaptive strategies to tackle the changing problems of gully erosion and land degradation also requires an understanding of how these elements interact. Both Rampurhat-I (AUC = 0.970) and Bolpur-Sriniketan (AUC = 0.975) exhibit outstanding performance according to the ROC-AUC analysis, demonstrating the models’ potent capacity to discriminate between positive and negative classes. These high scores show that the models achieve great sensitivity and specificity across thresholds, demonstrating the high effectiveness of the predictor variable, Mean Composite Value, for classification. These models are quite dependable and robust for real-world binary classification tasks, as evidenced by the nearly perfect AUC values, which show no overlap in prediction scores between classes.

Recent Trends and Importance of Gully Erosion Studies

The role that gully erosion plays in land degradation, soil fertility loss, and landscape instability has drawn increased attention, making it a major environmental concern on a global scale. In environmental studies, knowing the causes and mechanisms of gully erosion has become more important as climate change worsens and human activity continues to affect land use patterns. Geographic Information Systems (GIS) and remote sensing have been employed in recent research to identify the main components causing gully erosion, with a focus on spatial modeling and forecasting approaches (Gessesse et al., 2016; Wang et al., 2019). A prominent development in the study of gully erosion is the combination of vegetative, hydrological, and morphometric elements into a single framework to better comprehend the geographical variability of erosion susceptibility. To assess how various factors interact to affect erosion dynamics, studies have increasingly used techniques like Principal Component Analysis (PCA) and Structural Equation Modeling (SEM) (Rahmati et al., 2017). In this way, the intricate relationships between topographic, hydrological, and land use characteristics may be captured, providing a more nuanced knowledge of the risk of gully erosion in different spatial locations. These methods have shown promise in locating hotspots for erosion and in calculating the quantitative risk of erosion in various settings. In the Himalayan foothills, for instance, Pandey et al. (2020) evaluated the risk of soil erosion using topographical analysis and TWI, emphasizing the part that soil saturation and water buildup play in the creation of gullies. Beven and Kirkby (1979) also employed TWI to model runoff processes, showing how gully formation in moderate terrain is influenced by water retention on slopes. These methods highlight how crucial terrain features are to erosion susceptibility, particularly in regions with moderate slopes and less obvious topographic relief. The protective effects of vegetation against soil erosion have been the subject of other investigations, including Zhao et al. (2019). Their research in China’s Loess Plateau showed that by stabilizing soil, boosting infiltration, and lowering surface runoff, plant cover can dramatically lower the hazards of erosion. This is consistent with Chen et al. (2017), who pointed out that by improving moisture retention and decreasing soil separation, denser vegetation helps to decrease gully erosion. The intensification of agriculture and urbanization is a major land use changes that exacerbate gully erosion. With more impermeable surfaces due to urbanization, water infiltration is decreased, and runoff is accelerated, resulting in concentrated water flow and increased soil erosion (Zhang et al., 2021). Particularly in locations with steep topography and little plant cover, it has been demonstrated that the conversion of wooded lands to agricultural or urban land uses increases vulnerability to gully erosion. The effects of land use changes in the Chambal Ravines, where growing agricultural activity has resulted in severe soil degradation and gully formation, have been observed by Mishra et al. (2015). Similarly, since plant loss exacerbates surface runoff and soil erosion, deforestation has been associated with an increase in gully erosion in semi-arid regions such as Maharashtra and the Western Ghats (Gadgil et al., 2016). The intricacy of erosion processes and the demand for integrated models that take into consideration a variety of variables affecting erosion dynamics are highlighted by these recent developments in gully erosion research. Even so, there is still a dearth of localized research that provides an in-depth, site-specific understanding of the interplay between these variables, despite the expanding corpus of literature on gully erosion. More in-depth, localized evaluations that consider the variations in terrain, land use, and vegetation cover in particular areas are still required, even if large-scale models have offered significant insights into regional erosion trends.

Significance and Importance of This Work

This study provides a substantial contribution to the expanding corpus of research on gully erosion, especially in the setting of West Bengal, India, where the phenomenon is a significant environmental problem that affects landscape stability and agricultural output. This study sheds light on the spatial heterogeneity of gully erosion susceptibility, which has received little attention in the past, by concentrating on the Rampurhat-I and Bolpur- Sriniketan blocks. This study offers a fine-grained analysis, providing more particular data and practical insights into land management and erosion control measures, whereas much previous research has concentrated on larger, regional dimensions. A comprehensive method for comprehending the intricate dynamics of gully erosion is provided by combining hydrological and land use data with a variety of environmental and morphometric factors, including the Topographic Wetness Index (TWI), Slope, Drainage Density, and Normalized Difference Vegetation Index (NDVI). This study employs sophisticated statistical techniques like Principal Component Analysis (PCA) and Structural Equation Modeling (SEM), which enable a more thorough understanding of how these factors interact and contribute to gully formation, in contrast to earlier research that might have concentrated on isolated variables or simple models. A new addition to the literature is the application of these sophisticated methods to identify the main causes of gully erosion. With the help of this study, erosion hotspots can be more precisely identified by demonstrating the interactions between topographic features, vegetation cover, and hydrological parameters. An effective framework for comprehending the direct and indirect effects of different factors on erosion susceptibility is provided by the incorporation of PCA and SEM. This framework is essential for putting into practice focused and reasonably priced erosion control strategies. This work’s significance stems from both its practical implications for sustainable land management and its scientific contribution. Policymakers, land use planners, and environmental managers who operate in regions vulnerable to gully erosion will find the findings extremely pertinent. The study makes it possible to implement more specialized treatments like afforestation, slope stabilization, and drainage system enhancement by identifying high-risk regions and measuring the relative contributions of various elements to erosion processes. Compared to more generalized, larger models, this concentrated focus guarantees that the recommendations are context-specific and more feasible to execute. Also, by looking at how vegetation degradation and land use changes affect erosion dynamics, this study provides important information about how urbanization, agricultural growth, and deforestation worsen gully erosion. This information can be used to guide land use policies that emphasize sustainable agriculture and ecosystem restoration. As part of a comprehensive strategy to reduce erosion and safeguard land resources, the results highlight the significance of preserving plant cover and enhancing soil conservation techniques. All things considered, this work closes a significant gap in the literature by offering a thorough, quantitative, and spatially explicit evaluation of the risk of gully erosion in two West Bengal blocks. The research provides new insights into the factors that contribute to gully erosion and gives practical solutions for erosion management and land conservation in the area by fusing localized data with sophisticated statistical modeling.

Why is This Study Better than Previous Studies and the Novelty of the Work?

In the context of West Bengal, this study adds some new findings and innovations to the body of knowledge already available on gully erosion. I list the main arguments for this study’s superiority over earlier research and emphasize the approach’s originality below. Although many other studies have concentrated on gully erosion assessments at the regional or national level, frequently using generalized models and broad forecasts, this study offers a much-localized two blocks: Rampurhat-I and Bolpur-Sriniketan. More accurate information on the spatial heterogeneity of gully erosion susceptibility is provided by this study’s smaller geographic emphasis. The results are more relevant to regional land management plans since these regions use high- resolution, detailed data. Because gully erosion varies greatly within even small geographical regions due to variations in topography, land use, and vegetation cover, this specific focus is especially crucial. On the other hand, it is possible that earlier research ignored how these parameters varied among smaller regions. Earlier research on gully erosion has usually concentrated on one factor, such as hydrology, vegetation, or topography, frequently ignoring the intricate relationships between these variables. This study is unique because it incorporates morphometric variables (like drainage density and slope), vegetation indices (like NDVI), and hydrological indicators (like Stream Power Index (SPI) and rainfall erosivity). This gives an in-depth understanding of how several elements work together to cause gully erosion. Prior research, such as Pandey et al. (2020) and Mishra et al. (2015), frequently concentrated on individual elements like rainfall or slope, but they failed to fully represent the intricacy of erosion processes caused by the interaction of several variables. The driving mechanisms causing erosion can be better understood thanks to this all-encompassing approach, especially in constantly fluctuating landscapes. The use of sophisticated statistical methods that enable a more thorough investigation of multidimensional interactions, such as Principal Component Analysis (PCA) and Structural Equation Modeling (SEM), is another noteworthy development of this subject. PCA simplifies complex data sets by identifying the main primary components that account for the majority of the variance in susceptibility to gully erosion. By offering a framework to investigate the direct and indirect effects of different factors on gully erosion, SEM goes one step further. These methods enable this investigation to reveal the fundamental relationships among variables such as vegetation cover, hydrological processes, and terrain features. On the other hand, previous research has frequently depended on correlation-based analysis or simpler statistical models, which may have overlooked more intricate, multivariate relationships. This study represents a novel integration in the field of gully erosion research by combining morphometric data (e.g., slope, stream density), hydrological variables (e.g., rainfall erosivity, SPI), and vegetation indices (e.g., NDVI). Although each of these components’ effects on erosion processes has been studied separately, this work is one of the first to integrate them into a thorough model that takes into consideration how these factors interact. This method is important because it enables a multifactorial understanding of gully erosion, where erosion hotspots are created by a combination of many causes, like plant loss, increased surface runoff, and topographic variability. Compared to previous, simpler models, the integrated framework’s innovation offers a more precise and context-specific erosion susceptibility prediction. Additionally, the study emphasizes the protective function of vegetation, which has drawn increasing focus in recent studies on gully erosion (Zhao et al., 2019; Chen et al., 2017). This study, however, employs the Normalized Difference Vegetation Index (NDVI) to measure vegetation density and its direct influence on soil stability and moisture retention, in contrast to many others that treat vegetation as a secondary element. This study’s strong negative relationship between the NDVI and the Bare Soil Index (BSI) confirms that regions with higher plant densities are much less likely to experience gully formation. This study highlights the vital significance of ecosystem restoration in reducing gully erosion by concentrating on the function of plant cover, which adds to and expands on previous research that has mostly concentrated on topographic and hydrological aspects. The study’s conclusions have practical implications for land management techniques in the research region; they are not only theoretical. This study focuses more on site-specific interventions for the Rampurhat-I and Bolpur- Sriniketan blocks than earlier research, which frequently concentrated on broad, regional-scale approaches to erosion prevention. This study renders it possible to implement targeted land management techniques, like afforestation, slope stabilization, and drainage network upgrades that are suited to the requirements of these areas by identifying erosion hotspots and calculating susceptibility scores for each location. This local focus guarantees the effective distribution of erosion control resources and the promotion of sustainable farming practices in high-risk erosion areas. Through the integration of numerous environmental parameters, sophisticated statistical methodologies, and localized insights, this study offers a fresh approach to comprehending and regulating gully erosion. It addresses significant gaps left by earlier research, which frequently relied on less integrated, simpler techniques, by providing a more thorough model of gully erosion susceptibility. This work is crucial to the scientific understanding of gully erosion and the development of workable mitigation methods because PCA and SEM are applied to identify the primary drivers of erosion while concentrating on site-specific land management.

CONCLUSION

This study provides a comprehensive investigation of the causes causing gully erosion in the West Bengal blocks of Rampurhat-I and Bolpur-Sriniketan. It accomplishes this by employing raster-based computations, basin morphometry, and terrain analysis to assess the correlation between various environmental factors and gully erosion susceptibility. The findings indicate that moderate gully erosion susceptibility is prevalent in both research locations and that the Topographic Wetness Index (TWI), drainage density, and flow length are significant drivers of erosion risk. These findings suggest that areas with higher moisture retention, insufficient drainage systems, and more concentrated flow paths are more prone to experience gully development. The study not only demonstrates the substantial impacts of land use changes and land cover degradation on erosion dynamics, but it also demonstrates the sharp decline in agricultural land and landholding sizes, which is indicative of the ongoing landscape degradation brought on by erosion. The study highlights the lack of conservation methods and preventative measures to control gully erosion, even though it is a major pedo- geomorphic issue in the region. The study suggests that while complete gully erosion eradication remains a challenging task, certain conservation measures, such as revegetation, contour cultivation, gully stability, and optimal runoff management, may reduce erosion vulnerability and mitigate the current degradation. Applying these strategies strategically in high-risk areas can help stabilize the soil, distribute runoff energy, and restore vegetation cover, all of which enhance the overall health of the soil and agricultural productivity in the affected areas. The results also suggest the implementation of integrated land management strategies that combine traditional knowledge with modern methods to successfully prevent soil erosion. To track the effectiveness of these conservation initiatives and assess how erosion patterns change over time, long-term monitoring is urgently needed. Future studies should focus on predicting future gully erosion hotspots more precisely by combining machine learning models with advanced remote sensing tools like satellite data and aerial drones. Building more dynamic erosion models that can forecast how erosion processes will be impacted by climate change also requires an understanding of the climatic parameters, particularly seasonal fluctuations in rainfall patterns and their interactions with soil properties. Setting community-based erosion management research as a top priority and working with area stakeholders will assist in developing long-term solutions that tackle the unique challenges faced by residents. Integrating the research findings into regional land management plans and policy frameworks will ensure that erosion control techniques are comprehensive, adaptable, and scalable, hence improving the region’s long-term environmental and socioeconomic sustainability. By using the recommended techniques, gully erosion of susceptibility can be reduced. Reforestation, particularly the planting of trees and other deep-rooted plants, is crucial for stabilizing the soil and enhancing its capacity to hold moisture. The establishment of pastures and grass cover in vulnerable areas can also assist in preventing surface soil loss since it enhances soil structure and mitigates the impact of precipitation. Slope stabilization techniques like contour cultivation and terracing are beneficial because they reduce surface runoff and encourage water penetration. Additionally, gully stabilization by the construction of check dams or silt fences can prevent further gully growth by slowing down water flow and allowing sediment deposition. Additionally, by efficiently controlling runoff and ensuring that water is distributed uniformly across the landscape, for instance, by utilizing diversion channels, the erosive effect of water can be reduced. Additionally, it’s important to avoid using hard materials to block gully channels because this might exacerbate erosion by diverting water flow in ways that deepen gullies that already exist. By keeping large amounts of organic matter in the soil through mulching and composting, soil structure and erosion resistance can be further improved. The study shows that these conservation techniques can significantly reduce erosion susceptibility and ease ongoing degradation when carefully implemented in high-risk places, even though it is still challenging to abolish gully erosion. These methods can assist in dispersing runoff energy, stabilizing the soil, and restoring vegetation cover when used appropriately. These actions will enhance the overall health of the soil and agricultural productivity in the affected areas. Based on the findings, integrated land management strategies that combine traditional knowledge and modern methods are recommended for efficiently reducing soil erosion. Long-term monitoring is vitally needed going forward to track the effectiveness of these conservation initiatives and assess how erosion patterns change over time. Future studies should focus on integrating machine learning models with state-of-the-art remote sensing tools, like satellite photography and aerial drones, to more precisely predict future gully erosion hotspots. Furthermore, understanding the climatic factors, particularly the seasonal fluctuations in rainfall patterns and their interactions with soil properties, is crucial to developing more dynamic erosion models that can forecast how climate change may impact erosion processes. Research on community-based erosion management strategies should also be prioritized to develop long-term solutions that tackle the unique challenges faced by local populations. Working together with regional stakeholders, this can be accomplished. By integrating the research findings into regional land management plans and policy frameworks, erosion control measures will be comprehensive, adaptable, and scalable, thereby contributing to the region’s long-term environmental and socioeconomic sustainability.

ACKNOWLEDGEMENT

The authors would like to extend thanks and appreciation to Netaji Subhas Open University for ensuring supportive funding (Reg./0519, 2024-2025) and research facilities for the research work.

Competing Interests

The authors reported no potential conflict of interest.

REFERENCES

1. Agricultural Census. (2001). District Tables (Number and Area of Operational Holdings). Department of Agriculture and Farmers Welfare, Ministry of Agriculture and Farmers Welfare, Government of India.
2. Ahmadpour, , Bazrafshan, O., Rafiei-Sardooi, E., Zamani, H., & Panagopoulos, T. (2021). Gully erosion susceptibility assessment in the Kondoran Watershed using machine learning algorithms and the Boruta feature selection. Sustainability, 13(18), 10110. https://doi.org/10.3390/su131810110
3. Ali, A., Aslam, F., & Zeb, J. (2022). Integrated analysis of soil erosion risk in the Himalayan region. Environmental Science and Pollution Research, 29(5), 7889–7903. https://doi.org/10.1007/s11356-021- 17156-9
4. Ali, H., et al. (2022). Modeling the dynamics of gully erosion and its susceptibility in hilly terrain: A case study from the Himalayan region. Environmental Management, 45(2), 233-245. https://doi.org/10.1007/s00267-022-01637-2
5. Almasalmeh, O., Saleh, A. A., & Mourad, K. A. (2022). Soil erosion and sediment transport modelling using hydrological models and remote sensing techniques in Wadi Billi, Egypt. Model Earth Systems and Environment, 8(1), 1215-1226. https://doi.org/10.1007/s40808-021-01144-1
6. Arabameri, A., Tiefenbacher, J. P., Blaschke, T., Pradhan, B., & Tien Bui, D. (2020). Morphometric analysis for soil erosion susceptibility mapping using novel GIS-based ensemble model. Remote Sensing, 12(5), 874. https://doi.org/10.3390/rs12050874
7. Aslam, B., Maqsoom, A., Alaloul, W. S., Musarat, M. A., Jabbar, T., & Zafar, A. (2021). Soil erosion susceptibility mapping using a GIS-based multi-criteria decision approach: Case of district Chitral, Pakistan. Ain Shams Engineering Journal, 12(2), 1637-1649. https://doi.org/10.1016/j.asej.2020.09.015
8. Ayele, N. A., Naqvi, H. R., & Alemayehu, D. (2020). Rainfall induced soil erosion assessment, prioritization and conservation treatment using RUSLE and SYI models in highland watershed of Ethiopia. Geocarto International, 37(9), 2524-2540. https://doi.org/10.1080/10106049.2020.1715310
9. Balabathina, V. N., Raju, R. P., Mulualem, W., & Tadele, G. (2020). Estimation of soil loss using remote sensing and GIS-based universal soil loss equation in northern catchment of Lake Tana Sub-basin, Upper Blue Nile Basin, Northwest Environmental Systems Research, 9(1), 1-32. https://doi.org/10.1186/s40068- 020-00203-3
10. Bandyopadhyay, (1988). Drainage evolution in badland terrain at Gangani in Medinipur District, West Bengal. Geographical Review of India, 50(3), 10-20.
11. Bandyopadhyay, , et al. (2021). Impact of land use change on gully erosion in the Central Indian Plateau. Environmental Science and Pollution Research. https://doi.org/10.1007/s11356-021-15885-6
12. Beven, , & Kirkby, M. J. (1979). A physically based, variable contributing area model of basin hydrology. Hydrological Sciences Journal, 24(1), 43-69. https://doi.org/10.1080/02626667909491834
13. Bhattacharyya, A., Das, S., & Paul, A. (2020). Geospatial assessment of gully erosion susceptibility in Birbhum District, West Bengal, India. Environmental Earth Sciences, 79(5), 1-12. https://doi.org/10.1007/s12665-020-08959-0
14. Borrelli, , Poesen, J., Vanmaercke, M., Ballabio, C., Hervás, J., Maerker, M., Scarpa, S., & Panagos, P. (2022). Monitoring gully erosion in the European Union: A novel approach based on the Land Use/Cover Area frame survey (LUCAS). International Soil and Water Conservation Research, 10(1), 17-28. https://doi.org/10.1016/j.iswcr.2021.09.002
15. Brady, C., & Weil, R. R. (2010). Elements of the nature and properties of soils. Pearson Prentice Hall.
16. Bureau of Applied Economics and Statistics (BAES). (2011). District Statistical Handbook Birbhum. Department of Statistics and Programme Implementation, Government of West Bengal, West Bengal, 1-89.
17. Burrough, P. A., McDonnell, R. A., & Lloyd, C. D. (2015). Principles of geographical information systems. Oxford University Press.
18. Carlin, J. B., & Doyle, L. W. (2000). Basic concepts of statistical reasoning: Standard errors and confidence Journal of Paediatrics and Child Health, 36(5), 502-505. https://doi.org/10.1046/j.1440- 1754.2000.00588.x
19. Census of India District (2001). Census Handbook Birbhum, Village and Town Wise Primary Census Abstract (PCA). Directorate of Census Operations, West Bengal, Series-20, Part-A and B.
20. Census of India District (2011). Census Handbook Birbhum, Village and Town Wise Primary Census Abstract (PCA). Directorate of Census Operations, West Bengal, Series-20, Part XII-B, 1-464.
21. Chander, G., Markham, B. L., & Helder, D. L. (2009). Summary of current radiometric calibration coefficients for Landsat MSS, TM, ETM+, and EO-1 ALI sensors. Remote Sensing of Environment, 113(5), 893-903. https://doi.org/10.1016/j.rse.2009.01.007
22. Chaudhuri, , Gupta, S., & Majumder, R. (2017). Application of multiple discriminant analysis for gully erosion susceptibility mapping in Uttarakhand, India. Environmental Management, 60(1), 112-123. https://doi.org/10.1007/s00267-017-0841-6
23. Chen, , Pourghasemi, H. R., & Kalantari, Z. (2017). GIS-based gully erosion susceptibility modeling.
24. Geomatics, Natural Hazards and Risk, 8(2), 190–211. https://doi.org/10.1080/19475705.2017.1294864
25. Chen, , et al. (2017). Impact of vegetation on soil erosion: A case study of the Loess Plateau, China. Journal of Environmental Management, 203, 130-138. https://doi.org/10.1016/j.jenvman.2017.07.047
26. Choudhury, M. K., & Nayak, T. (2003). Estimation of soil erosion in Sagar Lake catchment of Central In Proceedings of the international conference on water and environment (15), 387-392.
27. Congalton, R. G., & Green, K. (1999). Assessing the accuracy of remotely sensed data: Principles and practices. CRC Press.
28. Das, , & Bandyopadhyay, S. (1996). Badland development in a lateritic terrain: Santiniketan, West Bengal. National Geographic, 31(1), 2.
29. Das, , et al. (2021). Assessment of soil erosion in Jharkhand using GIS and MLC. Environmental Science & Technology. https://doi.org/10.1021/acs.est.1c04219
30. Debanshi, , & Pal, S. (2020). Assessing gully erosion susceptibility in Mayurakshi river basin of eastern India. Environmental Development and Sustainability, 22(2), 883-914. https://doi.org/10.1007/s10668-018- 0224-x
31. Demir, V., & Kisi, O. (2016). Soil erosion analysis using RUSLE in the tropical climate region of South America. Environmental Monitoring and Assessment, 188(9), 495.
32. Deolia, , & Pande, A. (2014). Spatial distribution of dissection index (erosion intensity) versus geomorphological environment in Parkha Watershed, Central Himalaya.
33. Engelhart, D. (1941). The analysis of variance and covariance techniques in relation to the conventional formulas for the standard error of a difference. Psychometrika, 6(4), 221-233.
34. (2020). ArcGIS Spatial Analyst Hydrology. Retrieved from https://pro.arcgis.com/en/pro- app/latest/tool-reference/spatial-analyst/hydrology.htm
35. Faniran, (1968). The index of drainage intensity: a provisional new drainage factor. Australian Journal of Science, 31(9), 326-330.
36. Fawcett, (2006). An introduction to ROC analysis. Pattern Recognition Letters, 27(8), 861-874.
37. Foody, G. M. (2002). Status of land cover classification accuracy assessment. Remote Sensing of Environment, 80(1), 185-201.
38. Gadgil, M., et al. (2016). Conservation of the Western Ghats: An integrated approach to sustainable management. Journal of Environmental Management, 182, 189-197.
39. Gadgil, M., et al. (2016). Role of vegetation in soil conservation in the Western Ghats. Current Science, 111(1), 14-22.
40. Gao, Y., et (2020). Assessing erosion risks in agricultural landscapes using MLC and remote sensing. Journal of Soil and Water Conservation.
41. Gayen, A., Saha, S., & Pourghasemi, H. R. (2020). Soil erosion assessment using RUSLE model and its validation by FR probability model. Geocarto International, 35(15), 1750-1768. https://doi.org/10.1080/10106049.2019.1581272
42. Gessesse, B., Bewket, W., & Bräuning, A. (2016). Soil erosion modeling and its implications for sustainable land management. Environmental Earth Sciences, 75(9), 847.
43. Gessesse, B., et al. (2016). Geospatial analysis of soil erosion susceptibility in the Ethiopian Highlands using remote sensing and GIS. Environmental Monitoring and Assessment, 188(5), 302-314.
44. Gessesse, D., et al. (2016). Influence of topography and hydrology on gully erosion in the Ethiopian Highlands. Soil and Tillage Research, 155, 118-126.
45. Getnet T, Mulu A (2021) Assessment of soil erosion rate and hotspot areas using RUSLE and multi-criteria evaluation technique at Jedeb watershed, Upper Blue Nile, Amhara Region, Ethiopia. Chall. 4:100174. https://doi.org/10.1016/j.envc.2021.100174
46. Getnet, M., et al. (2021). Soil erodibility and erosion susceptibility of the Birbhum region. Geosciences and Environmental Studies, 19(1), 45-
47. Geyik MP   (1968)   FAO   watershed   management   field   manual:   gully
48. Ghosal K, Bhattacharya SD (2021) Identification of the relationship between temporally varying land surface temperature of winter season with the cover management factor of revised universal soil loss equation: A case study from upper Bakreshwar river basin. In: Banerjee T (ed) South Asian Institute for Advanced Research and Development (SAIARD), Kolkata, India: Geoinformatics in Research & Development, pp 60-
49. Ghosh KG, Saha S (2015) Identification of Soil Erosion Susceptible Areas in Hinglo River Basin, Eastern India Based    on    Geo-Statistics.    Univers    j    environ    res    technol    5(3).
50. Ghosh S (2011) Quantitative and spatial analysis of fluvial erosion in relation to morphometric attributes of Sarujharna  basin,  East  Singhbhum,    Int.  j  geomat  geosci  2(1):71-90.
51. Ghosh S, Bhattacharya K (2012) Multivariate erosion risk assessment of lateritic badlands of Birbhum (West Bengal, India): A case study. J Earth Syst Sci 121(6):1441-1454. https://doi.org/10.1007/s12040-012- 0243-1
52. Ghosh S, Dolui, T (2011) Predicting Soil Erosion by Water and its Management in the Lateritic Areas of Western Rampurhat I Block (Birbhum, West Bengal). Intj curr                    res rev         03:143-161.
53. Ghosh S, Guchhait SK, Illahi R A, Bera S, Roy S (2021) Geomorphic character and dynamics of gully morphology, erosion and management in laterite Terrain: few observations from Dwarka-Brahmani Interfluve, Eastern Geol  ecol  landsc  6(3):  188-216.  https://doi.org/10.1080/24749508.2020.1812148
54. Ghosh S, Kundu S (2022) Morphometric characterization and erosion assessment of gullies in the lateritic badlands of Eastern India using ALOS AW3D30 DEM and topographic indices. Geocarto Int 29:1-34. https://doi.org/10.1080/10106049.2022.2032390
55. Ghosh, (2011). Vegetation and ecological study of Birbhum region. Journal of Environmental Science, 42(3), 211-220.
56. Ghosh, , & Bhattacharya, S. (2012). Assessment of Gully Erosion in West Bengal Using Remote Sensing and GIS. FRONTIERS.
57. Gomez B, Banbury K, Marden M, Trustrum NA, Peacock DH, Hoskin PJ (2003) Gully erosion and sediment production: Te Weraroa Stream, New Zealand. Water Resour Res 39(7). https://doi.org/10.1029/2002WR001342
58. Gómez-Gutiérrez Á, Conoscenti C, Angileri SE, Rotigliano E, Schnabel S (2015) Using topographical attributes to evaluate gully erosion proneness (susceptibility) in two Mediterranean basins: Advantages and Nat Hazards 79(1):291-314. https://doi.org/10.1007/s11069-015-1703-0
59. Green, M., & Swets, J. A. (1966).  Signal Detection Theory and Psychophysics. Wiley.
60. Gupta, , et al. (2021). Gully Erosion in the Northeastern States of India: Remote Sensing and GIS-Based Study. Springer.
61. Han-Qiu (2005) A study on information extraction of water body with the modified normalized difference water index (MNDWI).                  J                        Remote                       Sens 9:589-595.
62. Horton RE (1932) Drainage-basin Eos, Transactions, American Geophysical Union 13(1):350-361. https://doi.org/10.1029/TR013i001p00350
63. Horton, E (1945) Erosional development of streams and their drainage basins: Hydrophysical approach to quantitative morphology. Geol Soc Am Bull 56:275-370. https://doi.org/10.1130/0016- 7606(1945)56[275:EDOSAT]2.0.CO;2
64. Igwe PU, Eze CP, Ikeji CA, Uzoegbu CA, Emeh AB (2017) A Review of Rainfall Erosivity as a Natural Factor of Gully Int j environ agric biotech 2(6):239012. http://doi:10.22161/ijeab/2.6.49
65. India Meteorological Department (2001) Climatological Table. Ministry of Earth Sciences. Government of India.
66. Indian Meteorological        (2001).    Daily    temperature    and  rainfall  data (2001). https://www.imd.gov.in/pages/data.php
67. Ionita I, Fullen MA, Zgłobicki W, Poesen J (2015) Gully erosion as a natural and human-induced Nat Hazards 79(1):1-5.            https://doi.org/10.1007/s11069-015-1935-z
68. Jahantigh, M, Pessarakli M (2011) Causes and effects of gully erosion on agricultural lands and the Soil Sci Plant Anal 42(18):2250-2255. https://doi.org/10.1080/00103624.2011.602456
69. Jenson SK, Domingue JO (1988) Extracting topographic structure from digital elevation data for geographic information  system    Photogramm  Eng  Remote  Sensing  54(11):1593-1600.
70. Jha V (2008) Land Degradation and Desertification and Integrated Management of Laterite Surface in Birbhum District Using Field and Remote Sensing DST (WB) sponsored project report, 1-145.
71. Jha VC, Kapat S (2009) Rill and gully erosion risk of lateritic terrain in South-Western Birbhum District, West Bengal,                       Sociedade  &                Natureza 21:141-158.
72. Jha VC, Kapat S (2011) Degraded lateritic soils cape and land uses in Birbhum district, West Bengal, Sociedade & Natureza 23:545-558.
73. Jha, R., et al. (2019). Geospatial Analysis of Gully Erosion in the Western Ghats Using GIS and Remote Environmental Monitoring                              and Assessment.
74. Jin XM, Zhang YK, Schaepman ME, Clevers JG, Su Z, Cheng J, Jiang J, van Genderen J (2008) Impact of elevation and aspect on the spatial distribution of vegetation in the Qilian mountain area with remote sensing data. International Society    for    Photogrammetry    and    Remote    Sensing.
75. Joint Research Centre-European Commission (2008) Handbook on constructing composite indicators: methodology and                             user              guide.            OECD               publishing,                             1-
76. Jolliffe,     T.    (2002).    Principal    Component    Analysis.    Springer.
77. Joseph S, Anitha K, Srivastava VK, Reddy C, Thomas AP, Murthy MS (2012) Rainfall and elevation influence the local-scale distribution of tree community in the southern region of Western Ghats biodiversity hotspot (India). Int J For Res 2012(576502):1-10. https://doi.org/10.1155/2012/576502
78. Kumar, , Singh, R., & Gupta, A. (2012). Land use/land cover change detection and its impact on gully erosion in the Deccan Plateau, India. Hydrology and Earth System Sciences, 16(4), 1227-1239. https://doi.org/10.5194/hess-16-1227-2012
79. Lal, R., et al. (2021). Soil erosion and degradation in tropical ecosystems and their impact on sustainable development. Soil and Tillage Research, 202, 104529. https://doi.org/10.1016/j.still.2021.104529
80. Le Roux, , & Gertner, G. (2007). Assessing gully erosion in South Africa using GIS and remote sensing: A case study in the Eastern Cape. Environmental Management, 40(6), 907-920. https://doi.org/10.1007/s00267-007-9021-6
81. Lee, , et al. (2020). Ensemble of Machine-Learning Methods for Predicting Gully Erosion Susceptibility in Garhbeta, West Bengal. MDPI. https://doi.org/10.3390/rs12030458
82. Lemenkova, P. (2016). Flow Direction and Length Determined by ArcGIS Spatial Analyst and Terrain Elevation Data Sets. In: NA Shherbakova, AN Bondarenko (ed). https://www.researchgate.net/publication/318659478_Proceedings_of_Young_Scientists_Conference
83. Leopold, B. (1994). A View of the River. Harvard University Press.
84. Lillesand, M., & Kiefer, R. W. (2000). Remote Sensing and Image Interpretation. Wiley.
85. Lu, , et al. (2014). “Supervised Classification and Its Application to Land Use/Cover Change Analysis.” Geo-Spatial Information Science.
86. Martz, L. W., & Garbrecht, J. (1992). Numerical definition of drainage network and subcatchment areas from digital elevation models. Computational Geosciences, 18(6), 747-761. https://doi.org/10.1016/0098- 3004(92)90007-E
87. Metz, E. (1978). Basic principles of ROC analysis. Seminars in Nuclear Medicine, 8(4), 283–298.
88. Mishra, M., et al. (2015). Sediment transport and gully erosion in the Chambal Ravines, India: A review of the factors contributing to land degradation. Environmental Earth Sciences, 74(4), 311-320.
89. Mishra, M., et al. (2015). Sediment transport and gully erosion in the Chambal Ravines, India: A review of the factors contributing to land degradation. Environmental Earth Sciences, 74(4), 311-320.
90. Mishra, , et al. (2015). Application of remote sensing and GIS for the study of gully erosion in the Chambal Ravines, India. Environmental Earth Sciences, 73(8), 4405-4417.
91. Mondal, A., Khare, D., & Kundu, S.A. (2018). Comparative study of soil erosion modelling by MMF, USLE and RUSLE. Geocarto International, 33(1), 89-103. https://doi.org/10.1080/10106049.2016.1232313
92. Mondal, (2013). Morphometric analysis of Birbhum district. Asian Journal of Multidisciplinary Studies, 1(4), 141-148.
93. Mondal, , et al. (2023). GIS-Based Assessment of Gully Erosion in Birbhum, West Bengal. MDPI.
94. Mondal, M., et al. (2023). Gully Erosion and Soil Degradation in Birbhum, West Bengal: A GIS and Remote Sensing Approach. MDPI.
95. Moore, D., Grayson, R. B., & Ladson, A. R. (1991). Digital terrain modelling: A review of hydrological, geomorphological, and biological applications. Hydrological Processes, 5(1), 3-30. https://doi.org/10.1002/hyp.3360050103
96. Mosavi, A., Golshan, M., Janizadeh, S., Choubin, B., Melesse, A. M., Dineva, A. A. (2020). Ensemble models of GLM, FDA, MARS, and RF for flood and erosion susceptibility mapping: A priority assessment of sub-basins. Geocarto International, 37(9), 2541-2560. https://doi.org/10.1080/10106049.2020.1829101
97. Mukherjee, , Dey, S., & Sharma, P. (2018). GIS-based gully erosion modeling in Birbhum district, West Bengal: A case study. International Journal of Geographical Information Science, 32(10), 1898-1911. https://doi.org/10.1080/13658816.2018.1500968
98. Mustafa, E. K., Co, Y., Liu, G., Kaloop, M. R., Beshr, A. A., Zarzoura, F., Sadek, M. (2020). Study for predicting land surface temperature (LST) using Landsat data: A comparison of four Advances in Civil Engineering, 2020, 7363546. https://doi.org/10.1155/2020/7363546
99. Nardo, M., Saisana, M., Saltelli, A., Tarantola, S., Hoffman, H., Giovannini, E. (2009). Handbook on Constructing Composite Indicators: Methodology and User Guide. Organisation for Economic Co-operation and Development: Paris, France.
100. Nugraha, S., Gunawan, T., Kamal, M. (2019). Comparison of Land Surface Temperature Derived from Landsat 7 ETM+ and Landsat 8 OLI/TIRS for Drought Monitoring. In: IOP Conference Series: Earth and Environmental Science, 313(1):012041. IOP Publishing.
101. Nwilo, C., Ogbeta, C. O., Daramola, O. E., Okolie, C. J., Orji, M. J. (2021). Soil Erosion Susceptibility Mapping of Imo River Basin Using Modified Geomorphometric Prioritisation Method. Quaestiones Geographicae, 40(3), 143-162.
102. Oguro, , Ito, S., Tsuchiya, K. (2011). Comparisons of brightness temperatures of Landsat-7/ETM+ and Terra/MODIS around Hotien Oasis in the Taklimakan Desert. Applied Environmental Soil Science. https://doi.org/10.1155/2011/948135
103. Pal, (2016). Identification of soil erosion vulnerable areas in Chandrabhaga river basin: A multi-criteria decision approach. Model. Earth Syst Environ, 2(1):1-11. https://doi.org/10.1007/s40808-015-0052-z
104. Pal, M. (2005). Random Forest Classifier for Remote Sensing Classification. International Journal of Remote Sensing, 26(1), 217-222.
105. Pal, , & Mather, P. M. (2005). Classification of remote sensing imagery using support vector machines. International Journal of Remote Sensing, 26(5), 1073-1086.
106. Pal, S., Talukdar, S., & Mondal, S. (2020). Assessing gully erosion susceptibility using remote sensing and GIS techniques in Puruliya district, West Bengal. Geocarto International, 35(8), 891–907.
107. Pandey, , et al. (2020). Erosion susceptibility assessment in the semi-arid regions of India using geospatial tools. Environmental Earth Sciences, 79, 202.
108. Pandey, , et al. (2020). Erosion susceptibility assessment in the semi-arid regions of India using geospatial tools. Environmental Earth Sciences, 79, 202
109. Pandey, A., et al. (2020). Soil erosion risk assessment in hilly areas using TWI and terrain analysis: A case study from northern     Journal of Hydrology, 585, 124759. https://doi.org/10.1016/j.jhydrol.2020.124759
110. Pandey, , et al. (2022). Multi-Factor Modeling of Gully Erosion in the Deccan Plateau Using GIS and Remote Sensing. Springer.
111. Pascoe, H. (1973). A Manual of the Geology of India and Burma (Vol. 1). Controller of Publications.
112. Pathan, S. A., & Sil, B. S. (2022). Prioritization of soil erosion-prone areas in the upper Brahmaputra River Basin up to Majuli River Island. Geocarto International, 37(7), 1999-2017. https://doi.org/10.1080/10106049.2020.1810328
113. Patton, C., & Baker, V. R. (1976). Morphometry and floods in small drainage basins subject to diverse hydrogeomorphic controls.             Water Resources               Research,         12(5),           941-952. https://doi.org/10.1029/WR012i005p00941
114. Patton, C., & Schumm, S. A. (1975). Gully erosion, Northwestern Colorado: A threshold phenomenon. Geology, 3(2), 88-90. https://doi.org/10.1130/0091-7613(1975)3%3C88:GENCAT%3E2.0.CO;2
115. Pearson, K. (1897). Mathematical contributions to the theory of evolution. On a form of spurious correlation which may arise when indices are used in the measurement of organs. Proceedings of the Royal Society of London, 60(359–367), 489-498. https://doi.org/10.1098/rspl.1896.0076
116. Pearson, (1901). On lines of closes fit to system of points in space. The London, Edinburgh, and Dublin Philosophical              Magazine and      Journal                   of                      Science, 2,                  559-572. https://doi.org/10.1080/14786440109462720
117. Pearson, (1914). On certain errors with regard to multiple correlation occasionally made by those who have not adequately studied this subject. Biometrika, 10, 181-187. https://doi.org/10.2307/2331747
118. Poesen, J., & Govers, G. (2003). Gully erosion: Importance and research needs. Catena, 50(2–4), 1-4. https://doi.org/10.1016/S0341-8162(02)00155-2
119. Qiuji, C., & Chuting, L. (2015). Land Surface Temperature Retrieval Based on Landsat ETM+/TM: Taking Xi’an City as an Example. Open Cybernetics and Systemics Journal, 9(1), 228-235. https://doi.org/10.2174/1874110X01509010228 Qiu, L., et al. (2017). Driving Forces of Land Use and Land Cover Change in the Qinghai Lake Basin, China. Chinese Geographical Science, 27(5), 734-746. https://doi.org/10.1007/s11769-017-0911-1
120. Rahman, M. M., & Lateh, H. (2017). Meteorological drought in Bangladesh: Assessing, analyzing, and mapping. Theoretical and Applied Climatology, 127, 221-243. https://doi.org/10.1007/s00704-015-1622-7
121. Rahman, M. M., & Saha, S. K. (2009). Spatial analysis of land suitability for agriculture in the Ganges Delta region             of              Applied Geography,            29(1),              112-124. https://doi.org/10.1016/j.apgeog.2008.08.003
122. Ramasamy, M. (2005). Remote Sensing in Geomorphology. Rawat Publications.
123. Rathore, A. P. S., et al. (2011). Land Use Change Detection in Doon Valley (India): A Study Based on Remote Sensing and GIS. Journal of the Indian Society of Remote Sensing, 39(2), 157-166. https://doi.org/10.1007/s12524-010-0074-3
124. Rawat, J. S., & Kumar, M. (2015). Monitoring land use/cover change using remote sensing and GIS techniques: A case study of Hawalbagh block, district Almora, Uttarakhand, India. The Egyptian Journal of Remote Sensing and Space Science, 18(1), 77-84. https://doi.org/10.1016/j.ejrs.2015.02.002
125. Reddy, C. S., et al. (2014). Development of a 1 km monthly rainfall dataset for India. Theoretical and Applied Climatology, 116, 75-91. https://doi.org/10.1007/s00704-013-0925-4
126. Reis, S. (2008). Analyzing land use/land cover changes using remote sensing and GIS in Rize, North- East Turkey. Sensors, 8(10), 6188-6202. https://doi.org/10.3390/s8106188
127. Ren, , et al. (2020). Land use optimization for sustainable agriculture under ecological security constraints. Ecological Indicators, 114, 106271. https://doi.org/10.1016/j.ecolind.2020.106271
128. Riad, , et al. (2011). GIS-based approach for assessing land use/land cover in relation to land degradation in the Western Desert of Iraq. Egyptian Journal of Remote Sensing and Space Science, 14(1), 11- 19. https://doi.org/10.1016/j.ejrs.2011.04.001
129. Sahu, A. S., et al. (2020). Impacts of urbanization on land use/land cover change using geospatial techniques: A case study of Bhubaneswar city, India. The Egyptian Journal of Remote Sensing and Space Science, 23(1), 57-65. https://doi.org/10.1016/j.ejrs.2019.11.003
130. Sarkar, S., et al. (2021). Evaluation of groundwater potential zones in semi-arid areas using geospatial techniques and analytic hierarchy process (AHP): A case study of Purulia district, West Bengal, India. Geocarto International, 36(9), 1000-1021. https://doi.org/10.1080/10106049.2019.1643496
131. Satterthwaite, D., et al. (2010). Urbanization and its implications for food and farming. Philosophical Transactions of the Royal Society B, 365(1554), 2809-2820. https://doi.org/10.1098/rstb.2010.0136
132. Schneider, , et al. (2010). A new map of global urban extent from MODIS satellite data. Environmental Research Letters, 5(4), 044003. https://doi.org/10.1088/1748-9326/5/4/044003
133. Seto, K. C., et al. (2011). A meta-analysis of global urban land expansion. PLOS ONE, 6(8), e23777. https://doi.org/10.1371/journal.pone.0023777
134. Shalaby, , & Tateishi, R. (2007). Remote sensing and GIS for mapping and monitoring land cover and land-use changes in the northwestern coastal zone of Egypt. Applied Geography, 27(1), 28-41. https://doi.org/10.1016/j.apgeog.2006.09.004
135. Singh, A. (1989). Digital change detection techniques using remotely sensed data. International Journal of Remote Sensing, 10(6), 989-1003. https://doi.org/10.1080/01431168908903939
136. Singh, R. B., et al. (2015). Application of geospatial techniques in monitoring urban growth in Delhi, Urban Ecology, 1, 163-178. https://doi.org/10.1007/978-94-017-9491-0_10
137. Song, -P., et al. (2016). Global estimates of urban expansion to 2030 and direct impacts on biodiversity and carbon pools. Nature Communications, 7, 12623. https://doi.org/10.1038/ncomms12623
138. Srivastava, P. K., et al. (2013). Analyzing land use/land cover dynamics in part of the Gangetic Plain, India using RS and GIS techniques. Journal of the Indian Society of Remote Sensing, 41(2), 321-329. https://doi.org/10.1007/s12524-012-0203-8
139. Sun, Z., & Huang, S. (2021). Evaluating the effects of urbanization on land surface temperature using multi-sensor satellite data: A case study of Wuhan, China. Environmental Research, 195, 110748. https://doi.org/10.1016/j.envres.2020.110748
140. Tewolde, M. G., & Cabral, P. (2011). Urban sprawl analysis and modeling in Asmara, Eritrea. Remote Sensing, 3(10), 2148-2165. https://doi.org/10.3390/rs3102148
141. Trivedi, R., et al. (2008). Spatial scale affects bioclimate model projections of climate change impacts on mountain plants. Global Change Biology, 14(5), 1089-1103. https://doi.org/10.1111/j.1365- 2486.2008.01553.x
142. Tucker, C. J. (1979). Red and photographic infrared linear combinations for monitoring vegetation. Remote Sensing of Environment, 8(2), 127-150. https://doi.org/10.1016/0034-4257(79)90013-0
143. S. Geological Survey (USGS). (2022a). Landsat 7 Enhanced Thematic Mapper Plus (ETM+), imagery from 2001. Retrieved from https://www.usgs.gov/landsat.
144. S. Geological Survey (USGS). (2022b). Shuttle Radar Topography Mission (SRTM) Digital Elevation Model, acquired in 2000. Retrieved from https://www.usgs.gov/
145. Verma, P., et al. (2018). Assessment of urban heat islands (UHI) in Varanasi city using remote sensing and GIS techniques. Journal of the Indian Society of Remote Sensing, 46(7), 1055-1068. https://doi.org/10.1007/s12524-018-0774-5
146. Weng, Q. (2001). A remote sensing-GIS evaluation of urban expansion and its impact on surface temperature in the Zhujiang Delta, China. International Journal of Remote Sensing, 22(10), 1999-2014. https://doi.org/10.1080/014311601300003600
147. Xu, (2008). A new index for delineating built-up land features in satellite imagery. International Journal of Remote Sensing, 29(14), 4269-4276. https://doi.org/10.1080/01431160802039957
148. Yang, X., & Lo, C. P. (2002). Using a time series of satellite imagery to detect land use and land cover changes in the Atlanta, Georgia metropolitan International Journal of Remote Sensing, 23(9), 1775-1798. https://doi.org/10.1080/01431160110075802
149. Zhao, M., et al. (2020). Urbanization-induced changes in local climate zones and impacts on climate trends in major Chinese Nature Communications, 11, 5906. https://doi.org/10.1038/s41467-020-19602- 8
150. Zhu, , & Woodcock, C. E. (2014). Continuous change detection and classification of land cover using all available Landsat data. Remote Sensing of Environment, 144, 152-171. https://doi.org/10.1016/j.rse.2014.01.011

APPENDIX

Fig 1 Location Map and Google Earth View of the study areas

Fig 2 A Conceptual Methodological Framework of the Study

Fig 3 FCC Satellite Imagery, Grids, Grid Points, Stream Network, and Stream Orders (Rampurhat-I)

Fig 4 FCC Satellite Imagery, Grids, Grid Points, Stream Network, and Stream Orders (Bolpur-Sriniketan)

Fig 5 Morphometric Parameters of the Extracted Micro Basins in Rampurhat-I

Fig 6 Morphometric Parameters of the Extracted Micro Basins in Bolpur-Sriniketan

Fig 7 Ground Pictures of the Gully Channels, and Analogous Reciprocal of the Factors

Fig 8 Indicators of Gully Erosion Susceptibility (Rampurhat-I)

Fig 9 Indicators of Gully Erosion Susceptibility (Bolpur-Sriniketan)

Fig 10 Spectral Indices of Gully Erosion Susceptibility (Rampurhat-I)

Fig 11 Spectral Indices of Gully Erosion Susceptibility (Bolpur-Sriniketan)

Fig 12 Scree Plots, Component Plots, Gully Erosion Susceptibility Zones, and ROC-AUC (Rampurhat-I)

Fig 13 Scree Plots, Component Plots, Gully Erosion Susceptibility Zones, and ROC-AUC (Bolpur-Sriniketan)

Fig 14 Structural Equation Model, Relationship between Rainfall Erosivity, and Its Determinants, and Overlapping of Gully Erosion Susceptibility Isoline on Rainfall Erosivity Index (Rampurhat-I)

Fig 15 Structural Equation Model, Relationship between Rainfall Erosivity, and Its Determinants, and Overlapping of Gully Erosion Susceptibility Isoline on Rainfall Erosivity Index (Bolpur-Sriniketan)

Fig 16 Land Use Land Cover Classes and Gully Hotspot Areas in Rampurhat-I, and Bolpur-Sriniketan

Fig 17 Ground Picture and Google Earth View of Gully Erosion Areas in Shantiniketan

Fig 18 Cross Profiles Across the Gully Erosion Areas in Rampurhat-I, and Bolpur-Sriniketan

Table 1: Sources and usage of available data

Sl. No.	Available Data	Data Features	Year(s)	Source(s)	Methods and techniques
1	SRTM-DEM: SRTM1 Arc- Second Global	Two Tiles of the Spatial Data, Entity ID: SRTM1N23E087V3-Resolution: 1-ARC, Coordinates: 23, 87 and Entity ID: SRTM1N24E087V3- Resolution: 1-ARC, Coordinates: 24, 87	Acquisition                            Date: February 11, 2000	USGS (2022b)	Digital elevation model (DEM), relief and slope analysis, drainage analysis, stream ordering
2	LANDSAT ETM+	Two Tiles of the Spatial Data, Entity ID: ELP139R043 7T20011026-  Path:  139, Row: 43, and Entity ID: ELP139R044 7T20011026-Path: 139, Row: 44 (Spatial Resolution: Bands 1-5 and 7: 30 meters Band 8 (panchromatic): 15 meters Band 6 (thermal infrared): 60 meters, Temporal resolution: 16 days	Acquisition Date: October 26, 2001	USGS (2022a)	LULC, NDVI, MNDWI, NDBI, BSI, LST
3	Rainfall (mm)	Attribute data	2001	India Meteorological Department (2001)	Rainfall erosivity index
4	Geographical area, household, population, agriculture, and landholdings	Attribute data	2001, 2000-2001	Census of India (2001), BAES (2001), Agricultural Census of India (2001)	Total geographical area, Total household, total population, Crop yield, Landholding sizes

Table 2: Sources and measurement of selected indicators

Sl. No.	Indicators	Measurement	Source(s)	Justification for selection	Sources of literature
1	Relative relief (R) in m	R = H – h, where, H = Highest relief	Smith (1950), Schumm (1956)	Elevation determines some other physical parameters, such as vegetation	Joseph et al. (2012), Jin et al. (2008), Deolia, and Pande (2014, p. 22)
		h = Lowest relief		distribution and rainfall pattern. “Relative relief determines a positive impact upon the dissection index.”
2	Dissection index (DI)	DI = R/Ra, where, R= relative relief, Ra is absolute relief	Singh and Dubey (1994)	There is a positive relationship between the dissection index and drainage density.	Deolia, and Pande (2014, p. 22)
3	Ruggedness number (Rn)	Ruggedness Number (Rn) Rn = Dd * (H / 1000), where, Dd=drainage density, H=basin relief	Patton and Baker (1976)	Rn is positively correlated with erodibility and erosivity.	Arabameri et al. (2020, p. 15)
4	Average slope (S) in %	S = (Z × (Ctl/H))/ (10 × A), basin area (A), total basin relief (H), the maximum height of the basin (Z) and total contour length, the average angle of slope (tanÕ) = Average no. of contour crossing per mile (A) × contour interval (I) 3361 (constant)	Wentworth (1930)	The slope indicates a directly proportional relationship with soil erosion.	Arabameri et al. (2020, p. 15)
5	Slope aspect (As)	The direction of the maximum slope	Skidmore (1989)	The slope aspect impacts gully formation by controlling vegetation coverage and soil characteristics.	Wang, Wei, and Horton (2011), Patton, and Schumm (1975)
6	Flow direction (Fdir)	Derived from DEM	Lemenkova, (2016), Martz and Garbrecht (1992)	Flow direction is associated with flow accumulation and the topographic wetness index.	Ahmadpour et al. (2021, p. 20)
7	Flow distance in km (Fdis)	Derived from DEM and flow direction raster	Jenson and Domingue (1988)	Upstream distance varies depending on the width of the terrain cross-section.	Gomez et al. (2003, p. 3)
8	Flow length in km (Fl)	Derived from stream raster	Lemenkova (2016)	Stream length is an important feature in measuring the morphometry of drainage basins of gully erosional areas.	Songu, Oyatayo, and Iorkua (2015, p. 105)
9	Stream frequency (Fs)	Fs = Nu/A, where, Nu= total number of streams, A= area in sq. km	Horton (1932)	Stream frequency has a negative relationship with infiltration and a positive relationship with basin roughness.	Arabameri et al. (2020, p. 13)
10	Drainage density (Dd) in km/square km	Dd = Lu/A, where, Lu= total length of streams, A= area in sq. km	Horton (1932), Schumm (1956)	Drainage density is associated with the dissection of terrain and basin runoff.	Arabameri et al. (2020, p. 13)
11	Infiltration number (If)	If = Fs × Dd, where, Fs= stream frequency, Dd= drainage density	Faniran (1968)	“Sub-basins with high infiltration values are less susceptible to soil erosion.”	Arabameri et al. (2020, p. 14)
12	Normalized Difference Vegetation Index (NDVI)	NDVI = 𝑁𝐼𝑅−𝑅𝐸𝐷 𝑁𝐼𝑅+𝑅𝐸𝐷	Rousse et al. (1973)	NDVI shows a positive correlation with greenage.	Van der Knijff, Jones, and Montanarella, (2000)
13	Modified Normalized Difference Water Index (MNDWI)	MNDWI = 𝐺𝑟𝑒𝑒𝑛−𝑆𝑊𝐼𝑅 𝐺𝑟𝑒𝑒𝑛+𝑆𝑊𝐼𝑅	Han-Qiu (2005)	An increase in the water index increases soil erosion.	Gómez-Gutiérrez et al. (2015)
14	Normalized Difference Built-up Index (NDBI)	NDBI = (SWIR – NIR)/ (SWIR + NIR)	Zha et al. (2003)	The built-up area is associated with the P- factor of the Universal Soil Loss Equation.	Balabathina et al. (2020, p. 18), Shin (1999)
15	Bare Soil Index (BI)	SI = ((Red+SWIR) (NIR+Blue)) / ((Red+SWIR) + (NIR+Blue))	Sangpradid (2018)	“Gully erosion and runoff rates are related to the amounts of bare soil.”	Zare et al. (2022, p. 12)
16	Land Surface Temperature (LST) in degree C	DN to Radiance: Lλ = ((LMAXλ ‐  MINλ)/(QCALMAX‐ QCALMIN)) * (QCAL‐QCALMIN) + LMINλ Radiance to BT (in K): T = K2 / ln (K1/Lλ + 1) where, LMAX and LMIN are spectral values contained in the metadata of Landsat images. The calibration values of pixels for QCALMAX and QCALMIN can be found in Landsat image metadata. K1 and K2 are predetermined constant values, and Lλ is the spectral radiance value of the image. Kelvin to degree C: C=K-273.15	Nugraha et al. (2019) Yaseen and Khan (2022) Mustafa (2020), Qiuji and Chuting (2015), Oguro, Ito, and Tsuchiya (2011), Chander et al. (2009)	“Massive temperature variation in this area influence the vegetation cover growth that directly influences the soil erosion.”	Ghosal, and Bhattacharya (2021, p. 62)
17	Rainfall Erosivity Index (REI)	R=79+0.363XA where, R is the rainfall erosivity, Xa is the average annual rainfall in mm over the study area	Choudhury and Nayak (2003)	High rainfall intensity increases soil erosion.	Igwe et al. (2017, p. 3156)
18	Stream Power Index (SPI)	SPI = 𝐴𝑠 × 𝑡𝑎𝑛𝛽 where, As is the specific catchment area (m2 m-1), β is the slope gradient (°)	Moore et al. (1991)	In the formation process of gullies, SPI controls the power of the erosiveness of land.	Vijith, and Dodge-Wan (2019)
19	Topographic Wetness Index (TWI)	𝐴𝑠  TWI = 𝑙𝑛(𝑡𝑎𝑛𝛽            )	Moore et al. (1991)	TWI influences the saturation of soil associated with gully formation.	Gómez-Gutiérrez et a. (2015)
		where, As is the specific catchment area (m2 m-1), β is the slope gradient (°)
20	Sediment Transport Index (STI)	STI= (As / 22.13)0.6 × sin (β/ 0.0896)1.3 where: As is the unit contributing area (in m2/m) and β is the slope angle (in degrees) at a given pixel	Burrough and McDonnell (1998)	“The sediment load would significantly influence the performance of any suggested stormwater management system.”	Almasalmeh, Saleh, and Mourad (2022, p. 1224)

Source: Selected by the authors

Table 3: Descriptive statistics of selected indicators (Rampurhat-I)

Descriptive Statistics
Indicators	N	Range	Minimum	Maximum	Sum	Mean		Std. Deviation	Variance	Skewness		Kurtosis
	Statistic	Statistic	Statistic	Statistic	Statistic	Statistic	Std. Error	Statistic	Statistic	Statistic	Std. Error	Statistic	Std. Error
R	55	31.94	15.02	46.96	1366.91	24.8528	1.02322	7.58842	57.584	.844	.322	.164	.634
DI	55	.97	.16	1.13	28.87	.5249	.03274	.24278	.059	.797	.322	-.263	.634
Rn	55	.42	.25	.67	26.71	.4856	.01209	.08967	.008	-.089	.322	.089	.634
S	55	9.82	.86	10.69	194.88	3.5432	.27168	2.01482	4.059	1.021	.322	1.244	.634
SA	55	304.38	33.83	338.21	10014.17	182.0759	10.36220	76.84813	5905.634	.014	.322	-.678	.634
Fdir	55	104.94	1.07	106.02	1486.99	27.0363	3.85426	28.58399	817.044	1.111	.322	.290	.634
Fdis	55	16.77	-3.12	13.65	216.50	3.9364	.45234	3.35463	11.254	.904	.322	.943	.634
Fl	55	2.28	.01	2.29	33.16	.6029	.08081	.59927	.359	1.140	.322	.497	.634
Fs	55	2.95	.00	2.95	61.21	1.1129	.07798	.57834	.334	.311	.322	.783	.634
DD	55	.21	.00	.21	3.41	.0620	.00725	.05379	.003	.925	.322	.190	.634
If	55	.35	.00	.35	4.63	.0842	.01178	.08739	.008	1.150	.322	.804	.634
NDVI	55	.48	-.19	.29	5.09	.0925	.01578	.11699	.014	-.411	.322	-.540	.634
MNDWI	55	.40	-.27	.14	-4.55	-.0827	.01025	.07598	.006	.112	.322	.304	.634
NDBI	55	.45	-.14	.31	2.64	.0481	.01524	.11302	.013	.151	.322	-.910	.634
BSI	55	131.08	205.78	336.86	13880.32	252.3695	3.74582	27.77975	771.714	1.108	.322	.962	.634
LST	55	4.84	22.63	27.48	1319.23	23.9860	.15066	1.11730	1.248	.987	.322	.615	.634
REI	55	.24	125.87	126.11	6932.08	126.0379	.00771	.05718	.003	-.880	.322	.121	.634
PI	55	.58	-.36	.22	-2.77	-.0503	.01792	.13293	.018	-.261	.322	-.648	.634
TWI	55	1730.31	1278.61	451.70	– 3982.17	-72.4031	32.75845	242.94320	59021.399	-2.239	.322	10.533	.634
STI	55	2.10	.00	2.10	10.31	.1874	.04949	.36705	.135	3.668	.322	15.337	.634
Valid N (listwise)	55

Source: Authors’ calculation

Table 4: Descriptive statistics of selected indicators (Bolpur-Sriniketan)

Descriptive Statistics
Indicators	N	Range	Minimum	Maximum	Sum	Mean		Std. Deviation	Variance	Skewness				Kurtosis
	Statistic	Statistic	Statistic	Statistic	Statistic	Statistic	Std. Error	Statistic	Statistic	Statistic	Std. Error			Statistic	Std. Error
R	52	35.00	7.00		42.00		1151.00	22.1346	1.15830	8.35259	69.766	.505	.330	-.363	.650
DI	52	5.76	.91		6.67		121.42	2.3349	.14819	1.06864	1.142	1.853	.330	4.957	.650
Rn	52	.34	.29		.63		25.07	.4820	.01111	.08012	.006	-.226	.330	.063	.650
S	52	7.19	1.51		8.70		172.97	3.3263	.17708	1.27694	1.631	2.215	.330	7.209	.650
Sa	52	286.97	42.21		329.18		9839.52	189.2215	11.43571	82.46409	6800.327	-.218	.330	-1.040	.650
Fdir	52	97.68	1.06		98.74		1498.56	28.8184	3.74705	27.02039	730.102	1.005	.330	.047	.650
Fdis	52	17.31	-1.81		15.50		154.40	2.9692	.51929	3.74465	14.022	2.026	.330	3.931	.650
Flen	52	.20	.00		.20		3.17	.0610	.00743	.05361	.003	.929	.330	-.003	.650
Fs	52	.79	.00		.79		21.22	.4081	.01845	.13302	.018	.080	.330	1.755	.650
DD	52	4.39	.00		4.39		59.44	1.1431	.13695	.98755	.975	1.590	.330	2.826	.650
If	52	2.00	.00		2.00		24.55	.4721	.06191	.44647	.199	1.477	.330	1.957	.650
NDVI	52	.85	-.51		.34		5.62	.1082	.02043	.14731	.022	-1.785	.330	4.916	.650
MNDWI	52	1.01	-.30		.71		-2.74	-.0526	.02049	.14776	.022	2.681	.330	13.630	.650
NDBI	52	.85	-.52		.33		.43	.0083	.02123	.15311	.023	-.142	.330	1.787	.650
BSI	52	.38	-.02		.36		12.11	.2329	.01409	.10160	.010	-1.026	.330	.161	.650
LST	52	5.60	22.39		27.99		1252.20	24.0807	.17325	1.24932	1.561	.928	.330	.288	.650
REI	52	1.38	128.46		129.84		6725.51	129.3367	.04889	.35257	.124	-.948	.330	.253	.650
SPI	52	.64	-.64		.00		-11.13	-.2141	.01478	.10659	.011	-1.283	.330	3.772	.650
TWI	52	3757.01	-3534.84		222.17		740.61	– 110.3963	70.15301	505.88055	255915.131	-6.305	.330	43.153	.650
STI	52	2.79	.00		2.79		8.51	.1636	.05962	.42994	.185	5.057	.330	28.706	.650
Valid N (listwise)	52

Source: Authors’ calculation

Table 5: Rotated component matrix of factor analysis (Rampurhat-I)

Rotated Component Matrixa
Factors/Indicators	Component
	1	2	3	4	5	6	7
R	.142		.120	.847			-.152
DI	-.454		.159	-.564	-.221	.181
Rn		.164			.829	.174
S	.219				.219	.718	.192
SA		-.110					.928
Fdir	.166		.225			-.730	.172
Fdis	.158	-.154	-.738	.270		.126
Fl		.706		.402	-.219
Fs		.818		-.105	-.211		-.141
DD		.811			.449
If		.911			.291
NDVI	-.838	-.143				.211
MNDWI	-.689	.182		-.204		-.360	-.203
NDBI	.940			.181			.138
BSI	.903					.166
LST	.925	-.112
REI	-.141	.522	.114	-.176	-.590	.297
SPI		-.108	.841	.137	-.207
TWI	.170	-.119	.532	-.385	.260		-.218
STI		.149	.737	.385	.128		.130
Extraction Method: Principal Component Analysis Rotation Method: Varimax with Kaiser Normalization.a a. Rotation converged in 6 iterations

Source: Authors’ calculation

Table 6: Total variance explained in PCA (Rampurhat-I)

Total Variance Explained
Component	Initial Eigenvalues			Extraction Sums of Squared Loadings			Rotation Sums of Squared Loadings
	Total	% of Variance	Cumulative %	Total	% of Variance	Cumulative %	Total	% of Variance	Cumulative %
1	4.458	22.292	22.292	4.458	22.292	22.292	4.119	20.594	20.594
2	3.177	15.883	38.175	3.177	15.883	38.175	3.126	15.631	36.225
3	2.236	11.179	49.354	2.236	11.179	49.354	2.220	11.100	47.325
4	1.658	8.291	57.645	1.658	8.291	57.645	1.731	8.656	55.982
5	1.496	7.478	65.123	1.496	7.478	65.123	1.654	8.272	64.254
6	1.346	6.728	71.852	1.346	6.728	71.852	1.430	7.148	71.402
7	1.036	5.181	77.032	1.036	5.181	77.032	1.126	5.631	77.032
8	.900	4.501	81.533
9	.733	3.666	85.199
10	.584	2.921	88.120
11	.531	2.653	90.772
12	.395	1.973	92.745
13	.369	1.843	94.588
14	.331	1.653	96.240
15	.258	1.289	97.530
16	.199	.995	98.525
17	.119	.593	99.118
18	.100	.502	99.620
19	.045	.224	99.844
20	.031	.156	100.000
Extraction Method: Principal Component Analysis

Source: Authors’ calculation

Table 7: Rotated component matrix of factor analysis (Bolpur-Sriniketan)

Rotated Component Matrixa
	Component
	1	2	3	4	5	6	7
R			-.135	-.912
DI		.110	-.179	.873	-.201
Rn	.103	.171	-.153		-.256		-.589
S	-.154	-.774	-.157	-.254	-.110
Sa	.154					.875
Fdir		.131		.152	-.109		.767
Fdis		.241	-.444	-.178	-.460	-.213	.287
Flen	.159	.121	.145	.169	-.563	-.378	-.161
Fs	-.230	-.108	.299	.439	.449		.202
DD		.121	.923		-.101
If		.123	.951				.148
NDVI	.794	.435				.130
MNDWI	.370	-.832		.107			-.124
NDBI	-.747	.565
BSI	.924	.185
LST	-.937	.128	-.114
REI	-.482				.179	.544
SPI	.349	.765					-.205
TWI				-.217	.634	-.188	.219
STI	.131	.137			.692		-.153
Extraction Method: Principal Component Analysis Rotation Method: Varimax with Kaiser Normalizationa a. Rotation converged in 8 iterations

Source: Authors’ calculation

Table 8: Total variance explained in PCA (Bolpur-Sriniketan)

Total Variance Explained

Component	Initial Eigenvalues			Extraction Sums of Squared Loadings				Rotation Sums of Squared Loadings
	Total	% of Variance	Cumulative %	Total	% of Variance	Cumulative %	Total		% of Variance	Cumulative %
1	3.693	18.466	18.466	3.693	18.466	18.466	3.585		17.924	17.924
2	2.844	14.222	32.688	2.844	14.222	32.688	2.636		13.179	31.103
3	2.403	12.016	44.704	2.403	12.016	44.704	2.203		11.013	42.116
4	1.925	9.626	54.331	1.925	9.626	54.331	2.023		10.113	52.229
5	1.650	8.251	62.582	1.650	8.251	62.582	1.806		9.030	61.259
6	1.249	6.245	68.827	1.249	6.245	68.827	1.342		6.708	67.967
7	1.099	5.497	74.324	1.099	5.497	74.324	1.272		6.358	74.324
8	.932	4.662	78.987
9	.891	4.454	83.440
10	.739	3.694	87.135
11	.610	3.049	90.183
12	.513	2.566	92.749
13	.436	2.180	94.929
14	.394	1.970	96.899
15	.271	1.357	98.256
16	.126	.630	98.886
17	.089	.447	99.333
18	.065	.327	99.660
19	.050	.250	99.910
20	.018	.090	100.000
Extraction Method: Principal Component Analysis

Source: Authors’ calculation

Table 9: Composite factor scores and mean composite factor scores of extracted points (Rampurhat-I)

id	Points	F1	Latitude	Longitude	F2	F3	F4	F5	F6	F7	Composite	Mean
1	1	0.77406	24.101	87.752	-1.21214	0.97682	-0.93895	-0.06792	-3.01023	-0.24282	-3.72	-0.53
2	2	-0.47439	24.101	87.780	-1.32805	-0.09565	-0.69625	-0.76695	-0.24901	1.58132	-2.03	-0.29
3	3	0.49281	24.118	87.615	-1.24917	-1.28309	1.18634	-0.66349	-0.72736	0.06729	-2.18	-0.31
4	4	1.8538	24.118	87.643	-1.00488	-1.06308	0.55428	0.37137	0.99357	-0.17611	1.53	0.22
5	5	-1.00742	24.118	87.670	-0.82416	-0.79322	1.28624	0.50287	-0.54184	-0.42969	-1.81	-0.26
6	6	-0.84846	24.118	87.698	-0.70063	-0.39432	0.62924	1.00475	0.68615	0.7555	1.13	0.16
7	7	-0.77626	24.118	87.725	1.2437	-0.34053	-0.10682	1.75701	-0.33954	0.87238	2.31	0.33
8	8	-1.37742	24.118	87.752	-0.70016	-1.18518	-0.08001	0.50089	0.94445	0.45863	-1.44	-0.21
9	9	-0.03958	24.118	87.780	-0.46305	-0.94424	-0.85596	-0.10738	0.28436	1.12653	-1.00	-0.14
10	10	1.26466	24.136	87.615	0.51282	-0.62436	0.02387	2.27166	0.70956	-0.89107	3.27	0.47
11	11	-0.78834	24.136	87.643	-0.22445	1.15059	0.79685	1.60177	0.22075	-0.99464	1.76	0.25
12	12	0.17742	24.136	87.670	-0.57314	-1.70104	1.40807	1.00207	-0.19787	0.31167	0.43	0.06
13	13	-0.85067	24.136	87.698	0.16333	-1.03575	1.05659	-0.00474	-0.43035	-1.05409	-2.16	-0.31
14	14	1.8311	24.136	87.725	1.34635	0.2019	-1.58574	-0.75105	0.10673	-0.92911	0.22	0.03
15	15	-0.78392	24.136	87.752	1.08324	-0.68604	-0.26601	1.32476	0.06101	0.54935	1.28	0.18
16	16	-0.81365	24.136	87.780	0.09552	-0.88303	-0.43561	-0.22163	-0.2814	-0.52432	-3.06	-0.44
17	17	0.16842	24.136	87.807	-0.25798	1.80056	-1.10127	-0.07143	0.22014	1.75998	2.52	0.36
18	18	-0.55488	24.153	87.588	-1.58316	0.46033	0.90658	1.96992	-1.98931	1.11397	0.32	0.05
19	19	0.23363	24.153	87.615	-0.48372	0.31328	0.17713	0.59536	-1.21241	0.01227	-0.36	-0.05
20	20	1.24813	24.153	87.643	-1.23831	2.17677	1.02079	0.42496	-1.63861	0.16411	2.16	0.31
21	21	0.57921	24.153	87.698	-0.94944	-0.76998	1.43565	-0.5522	1.97394	-0.90899	0.81	0.12
22	22	1.10351	24.153	87.725	0.82789	-0.26081	-0.31579	-0.11202	-0.63708	-0.531	0.07	0.01
23	23	-0.77838	24.153	87.752	-0.16458	-0.79177	0.587	-0.53138	0.40306	-0.03968	-1.32	-0.19
24	24	0.51257	24.153	87.780	0.23308	0.24429	-0.59851	-0.93298	-0.97607	-1.17853	-2.70	-0.39
25	25	-0.84055	24.153	87.807	-0.09022	-0.44682	-0.29921	0.48082	1.19527	0.28775	0.29	0.04
26	26	-0.72868	24.171	87.725	0.62804	-0.62645	2.6039	-2.65133	-0.05263	0.50415	-0.32	-0.05
27	27	-0.3291	24.171	87.752	0.80212	1.89639	-0.3741	-0.46968	1.00296	0.1219	2.65	0.38
28	28	0.97541	24.171	87.780	1.48527	-1.24622	-1.18717	0.00956	-0.84042	1.22486	0.42	0.06
29	29	-0.32129	24.171	87.807	-1.09418	0.30887	-0.84336	-0.06463	1.44522	0.84292	0.27	0.04
30	30	-1.46617	24.171	87.835	-1.10505	-0.33869	-0.79089	0.6278	0.97049	-0.3139	-2.42	-0.35
31	31	2.92247	24.188	87.698	-1.08729	1.28911	0.5228	-0.24077	2.83861	-0.86047	5.38	0.77
32	32	-0.09858	24.188	87.725	0.70382	1.28279	2.07574	-1.11563	-0.30215	0.20465	2.75	0.39
33	33	-0.86383	24.188	87.752	2.19209	0.59834	-0.84632	0.02182	-0.01772	-1.31449	-0.23	-0.03
34	34	0.79631	24.188	87.780	0.22819	0.44277	-0.41725	-0.99191	0.02769	0.51728	0.60	0.09
35	35	1.51687	24.188	87.807	1.57755	-1.11325	-1.09714	2.07818	-0.06979	0.13677	3.03	0.43
36	36	-0.45348	24.206	87.698	2.3617	2.34764	2.76993	1.14865	0.38383	0.08299	8.64	1.23
37	37	0.00747	24.206	87.725	1.10799	0.81761	1.15465	-1.18245	-0.90138	-1.39875	-0.39	-0.06
38	38	0.38897	24.206	87.752	1.19987	-1.05513	0.10639	0.18502	1.1487	0.70246	2.68	0.38
39	39	-1.09083	24.206	87.780	-0.37762	-0.24038	-0.40616	-0.45304	0.99379	-1.4682	-3.04	-0.43
40	40	-0.12378	24.206	87.807	0.93884	-1.04088	-1.34035	-0.2575	-0.67484	-1.82176	-4.32	-0.62
41	41	0.12315	24.223	87.643	-1.11062	0.04551	-0.82966	0.02107	-0.64423	-0.33492	-2.73	-0.39
42	42	0.40838	24.223	87.670	1.34779	0.04378	0.37932	0.53381	-0.68456	1.91956	3.95	0.56
43	43	0.6314	24.223	87.698	1.41916	-0.34549	0.25314	0.06281	0.72708	0.68022	3.43	0.49
44	44	-0.97063	24.223	87.725	1.01465	-0.254	0.24788	-0.93954	-0.20119	-1.30423	-2.41	-0.34
45	45	0.61786	24.223	87.752	0.44214	-1.06817	-0.54354	-0.84375	-1.41085	1.20896	-1.60	-0.23
46	46	-1.11052	24.223	87.780	-0.53435	0.68742	-0.75632	-0.47939	1.39515	-0.10693	-0.90	-0.13
47	47	-1.04057	24.223	87.807	-0.58819	1.12768	-1.53233	-0.41597	1.16198	0.94239	-0.35	-0.05
48	48	-0.19834	24.240	87.643	-1.4795	0.41325	0.11803	1.61726	-0.93838	-1.75352	-2.22	-0.32
49	49	1.61077	24.240	87.670	-0.97471	0.23683	-0.22353	0.05947	0.51729	-0.81166	0.41	0.06
50	50	1.76849	24.240	87.698	-0.8594	-0.22994	-0.14847	-1.13762	-0.42923	-0.80659	-1.84	-0.26
51	51	0.38597	24.240	87.725	-0.28254	-1.296	0.70264	-1.30858	-0.18921	1.71312	-0.27	-0.04
52	52	0.08941	24.240	87.752	0.02134	2.07488	-0.43265	-0.29648	0.42329	2.36271	4.24	0.61
53	53	-1.13277	24.240	87.780	-0.33011	-0.43355	-0.4861	-2.17661	-0.76776	-0.08279	-5.41	-0.77
54	54	-1.51799	24.240	87.807	-0.32404	1.02665	-1.18123	0.49316	-0.3272	-1.52062	-3.35	-0.48
55	55	-1.10172	24.258	87.780	0.21836	0.62298	-1.28632	-0.8588	-0.15247	-0.42681	-2.98	-0.43

Source: Authors’ calculation

Table 10: Composite factor scores and mean composite factor scores of extracted points (Bolpur-Sriniketan)

id	Points	Latitude	Longitude	F1	F2	F3	F4	F5	F6	F7	Composite	Mean
1	1	23.569	87.804	-1.30455	0.09856	-0.75834	-0.34103	0.13165	1.93122	1.45094	1.208	0.173
2	2	23.569	87.831	-1.13569	-5.83748	-0.14576	0.25575	0.41182	-0.584	-0.42872	-7.464	-1.066
3	3	23.591	87.777	-1.94225	0.79139	0.24612	-0.3285	0.61196	1.10021	-1.97008	-1.491	-0.213
4	4	23.591	87.804	-2.02024	-0.08983	1.18406	-0.2832	0.30453	-0.13602	0.19098	-0.850	-0.121
5	5	23.591	87.831	1.021	0.95119	-0.4195	0.07369	3.89275	0.00532	-1.3458	4.179	0.597
6	6	23.612	87.640	-0.26664	-0.03321	-0.94911	-1.59815	0.70246	1.33373	-0.39569	-1.207	-0.172
7	7	23.612	87.695	-2.43465	0.11031	-0.63572	-0.25417	0.29773	0.83902	-0.74492	-2.822	-0.403
8	8	23.612	87.722	-2.12103	-0.1053	-0.94595	1.02638	0.99777	-0.01966	0.06043	-1.107	-0.158
9	9	23.612	87.749	0.70699	0.0851	0.14377	1.02555	-0.08348	1.52258	0.72402	4.125	0.589
10	10	23.612	87.777	-1.08549	1.06247	0.30728	-1.61483	0.4888	0.18809	-0.20124	-0.855	-0.122
11	11	23.612	87.804	0.98653	-0.6848	0.02788	-1.21888	0.35166	0.75038	1.45805	1.671	0.239
12	12	23.633	87.640	0.40815	-0.98046	-0.45408	-1.21634	0.03645	1.5013	0.02952	-0.675	-0.096
13	13	23.633	87.668	-1.43113	0.52536	-0.39844	-1.04436	0.51397	-0.34212	0.25169	-1.925	-0.275
14	14	23.633	87.695	0.63543	0.54954	-0.33252	1.10257	1.1806	0.30549	0.39201	3.833	0.548
15	15	23.633	87.722	0.58646	0.08456	1.72747	0.27802	1.8668	1.07962	-0.31448	5.308	0.758
16	16	23.633	87.749	0.46223	0.52042	0.95403	0.59301	-0.13493	1.31829	-1.56681	2.146	0.307
17	17	23.633	87.777	0.819	0.80209	-0.73443	1.81595	0.08786	0.26197	0.23057	3.283	0.469
18	18	23.633	87.804	-0.12466	-0.15446	1.27355	0.30804	0.77527	0.79869	0.05728	2.934	0.419
19	19	23.655	87.640	-0.84412	-1.30501	-0.08323	0.28885	-0.1382	-0.42182	-0.8417	-3.345	-0.478
20	20	23.655	87.668	0.22917	0.5027	1.53047	-0.17363	-0.32125	-0.26645	0.35694	1.858	0.265
21	21	23.655	87.695	-1.20504	0.50304	0.12747	-0.02411	-0.08791	-1.0184	-1.46182	-3.167	-0.452
22	22	23.655	87.722	-0.48835	1.28508	-1.56404	0.46415	-1.05712	-0.4636	1.36746	-0.456	-0.065
23	23	23.655	87.749	1.02546	0.39519	2.49404	0.00557	-0.8972	-0.05683	-0.12386	2.842	0.406
24	24	23.655	87.777	1.10116	-0.31646	-0.50148	0.41556	-0.30226	1.08801	-0.42497	1.060	0.151
25	25	23.676	87.613	0.55512	-0.23972	-0.88532	1.0016	0.58972	0.13792	1.51534	2.675	0.382
26	26	23.676	87.640	0.5653	0.25973	1.57884	0.20953	-0.58373	-0.1492	-0.18693	1.694	0.242
27	27	23.676	87.668	-0.38702	0.00404	-0.03975	1.44147	-3.53476	1.41456	-1.75935	-2.861	-0.409
28	28	23.676	87.695	0.23022	0.14192	-1.23967	0.06883	-0.98243	0.50522	-0.02067	-1.297	-0.185
29	29	23.676	87.722	-0.31283	0.14013	-1.49707	-0.7286	-1.11029	0.92381	2.54009	-0.045	-0.006
30	30	23.676	87.749	-0.87886	1.42686	-1.83346	-1.56876	-0.79637	-1.28449	-0.18901	-5.124	-0.732
31	31	23.676	87.777	1.23997	-0.09352	0.6335	-0.7349	1.06559	0.02136	0.30756	2.440	0.349
32	32	23.697	87.613	-0.04087	0.73334	-0.05372	0.23141	-0.11188	0.71379	0.73455	2.207	0.315
33	33	23.697	87.640	-0.84434	0.92824	0.15957	2.17126	-0.37356	-0.68214	-1.08389	0.275	0.039
34	34	23.697	87.668	0.36428	-0.24643	0.74853	1.00886	-0.1621	-0.04043	1.23072	2.903	0.415
35	35	23.697	87.695	-0.81719	-0.27213	0.17721	-0.14548	-0.17519	-0.14951	0.46035	-0.922	-0.132
36	36	23.697	87.722	0.48722	0.65798	-0.43309	-0.99426	-0.10051	-1.50545	0.13706	-1.751	-0.250
37	37	23.697	87.749	0.99613	-0.2756	0.38586	-1.0934	0.1516	0.86245	-0.44035	0.587	0.084
38	38	23.697	87.777	1.48369	-0.50033	-0.53065	-2.19559	-1.00098	1.2669	-0.79659	-2.274	-0.325
39	39	23.719	87.613	1.11924	-0.36093	0.07759	-0.91376	-0.5414	0.18297	0.32886	-0.107	-0.015
40	40	23.719	87.640	0.70481	-0.58921	0.27665	-1.00943	-0.48445	-1.48381	1.39521	-1.190	-0.170
41	41	23.719	87.668	0.04932	0.25	-0.38977	1.51769	0.74005	-1.06996	1.42843	2.526	0.361
42	42	23.719	87.695	0.52076	-0.74761	-1.178	2.84632	-0.34053	-0.22926	0.15136	1.023	0.146
43	43	23.719	87.722	-1.79783	0.41715	2.96461	-0.28639	-0.0918	-1.45927	1.67987	1.426	0.204
44	44	23.719	87.749	0.9225	0.18369	-0.383	0.26366	0.27841	-1.24847	-0.95425	-0.937	-0.134
45	45	23.740	87.586	0.71841	0.2921	-0.8221	0.13861	1.01111	-1.484	0.21061	0.065	0.009
46	46	23.740	87.613	0.906	0.12445	-0.00066	-0.70658	-0.8981	-2.48302	-1.39184	-4.450	-0.636
47	47	23.740	87.640	0.73897	0.44634	-0.66649	-0.47302	-1.20754	-0.53823	-1.43926	-3.139	-0.448
48	48	23.740	87.668	0.42479	-0.89235	-0.40788	-0.55551	-0.07061	-1.51552	-1.3329	-4.350	-0.621
49	49	23.740	87.695	0.14681	0.11697	1.99357	0.40505	-1.25461	0.72871	0.32794	2.464	0.352
50	50	23.740	87.749	0.02993	-0.37525	0.37551	-0.31073	-0.5865	-1.50614	1.12472	-1.248	-0.178
51	51	23.761	87.668	0.934	-0.29033	-0.61723	0.68916	0.52068	-0.17823	-0.58117	0.477	0.068
52	52	23.761	87.749	0.36371	0.00047	-0.4871	0.16709	0.42045	-0.46561	-0.14626	-0.147	-0.021

Source: Authors’ calculation

Table 11: Correlation among the indicators of gully erosion (Rampurhat-I)

Indicators		R		DI	Rn	S	As		Fdir	Fdi s	Fl	Fs		Dd		If	ND VI		MN DW I	ND BI	BI	LS T	REI	SPI	TW I
R	Pearson Correlation	1		– .36 3**	– .01 5	– .05 8	– .12 7		.03 6	.10 4	.17 9	– .06 7		– .02 9		– .06 7	– .16 6		– .25 9	.25 6	.10 2	.24 9	– .20 3	.22 6	– .20 8
	Sig. (2- tailed)			.00 6	.91 5	.67 4	.35 7		.79 4	.44 8	.19 2	.62 7		.83 3		.62 9	.22 5		.05 6	.06 0	.45 8	.06 6	.13 8	.09 7	.12 8
	N	55		55	55	55	55		55	55	55	55		55		55	55		55	55	55	55	55	55	55
DI	Pearson Correlation	-.36 3**		1	– .08 9	– .05 2	– .00 6		– .14 8	– .33 0*	– .29 0*	.07 9		– .11 9		– .09 9	.40 0**		.32 8*	– .48 1**	– .30 0*	– .43 3**	.28 2*	.02 5	– .03 2
	Sig. (2- tailed)	.00 6			.51 7	.70 4	.96 4		.28 2	.01 4	.03 1	.56 6		.38 6		.47 3	.00 2		.01 4	.00 0	.02 6	.00 1	.03 7	.85 8	.81 7
	N		55	55	55	55		55	55	55	55	55		55		55	55		55	55	55	55	55	55	55
Rn	Pearson Correlation		– .01 5	– .08 9	1	.19 7		.01 7	– .01 5	– .05 2	.01 2	– .04 2		.40 2**		.31 2*	.01 3		– .04 9	– .00 4	.14 3	.04 9	– .23 8	– .14 7	.11 8
	Sig. (2- tailed)		.91 5	.51 7		.14 9		.90 0	.91 6	.70 4	.93 0	.76 2		.00 2		.02 0	.92 3		.72 5	.97 8	.29 7	.72 4	.08 0	.28 4	.39 0
	N		55	55	55	55		55	55	55	55	55		55		55	55		55	55	55	55	55	55	55
S	Pearson Correlation		– .05 8	– .05 2	.19 7	1		.10 7	– .20 8	.13 7	.00 9	– .18 7		.06 4		– .03 2	– .06 4		– .30 8*	.18 8	.25 4	.19 4	.06 1	– .05 4	.08 8
	Sig. (2- tailed)		.67 4	.70 4	.14 9			.43 5	.12 7	.31 7	.95 0	.17 2		.64 5		.81 9	.64 1		.02 2	.17 0	.06 1	.15 6	.65 6	.69 5	.52 5
	N		55	55	55	55		55	55	55	55	55		55		55	55		55	55	55	55	55	55	55
As	Pearson Correlation		– .12 7	– .00 6	.01 7	.10 7		1	.01 2	.07 7	– .08 8	– .14 8		– .08 3		– .13 4	.03 6		– .19 0	.11 3	.04 7	– .03 7	– .13 4	– .04 1	– .10 2
	Sig. (2- tailed)		.35 7	.96 4	.90 0	.43 5			.93 3	.57 7	.52 1	.28 2		.54 7		.32 9	.79 3		.16 4	.41 2	.73 3	.78 8	.32 9	.76 9	.45 8
	N		55	55	55	55		55	55	55	55	55		55		55	55		55	55	55	55	55	55	55
Fdir	Pearson Correlation		.03 6	– .14 8	– .01 5	– .20 8		.01 2	1	– .10 6	.03 1	– .14 2		.02 2		– .05 1	– .19 3		.06 8	.19 9	.04 3	.17 9	– .18 9	.13 2	.18 9
	Sig. (2- tailed)		.79 4	.28 2	.91 6	.12 7		.93 3		.43 9	.82 5	.30 3		.87 3		.71 4	.15 8		.62 3	.14 6	.75 4	.19 2	.16 8	.33 5	.16 7
	N		55	55	55	55		55	55	55	55	55		55		55	55		55	55	55	55	55	55	55
Fdis	Pearson Correlation		.10 4	– .33 0*	– .05 2	.13 7		.07 7	– .10 6	1	– .01 0	– .09 6		– .20 1		– .19 6	– .01 9		– .20 9	.16 5	.17 6	.10 1	– .12 1	– .49 4**	– .30 5*
	Sig. (2- tailed)		.44 8	.01 4	.70 4	.31 7		.57 7	.43 9		.94 0	.48 5		.14 1		.15 1	.89 3		.12 6	.22 8	.19 8	.46 2	.38 0	.00 0	.02 4
	N		55	55	55	55		55	55	55	55	55		55		55	55		55	55	55	55	55	55	55
Fl	Pearson Correlation		.17 9	– .29 0*	.01 2	.00 9		– .08 8	.03 1	– .01 0	1	.51 2**		.36 4**		.45 7**	– .04 5		.14 2	.06 4	– .13 4	– .11 1	.45 6**	.00 7	– .19 3
	Sig. (2- tailed)		.19 2	.03 1	.93 0	.95 0		.52 1	.82 5	.94 0		.00 0		.00 6		.00 0	.74 6		.30 1	.64 0	.33 1	.41 9	.00 0	.96 2	.15 9
	N		55	55	55	55		55	55	55	55	55		55		55	55		55	55	55	55	55	55	55
Fs	Pearson Correlation		– .06 7	.07 9	– .04 2	– .18 7		– .14 8	– .14 2	– .09 6	.51 2**	1		.49 8**		.70 3**	.01 7		.15 4	– .07 9	– .03 2	– .14 5	.45 5**	– .08 9	– .05 9
	Sig. (2- tailed)		.62 7	.56 6	.76 2	.17 2		.28 2	.30 3	.48 5	.00 0			.00 0		.00 0	.90 1		.26 2	.56 7	.81 5	.29 0	.00 0	.52 0	.66 8
	N		55	55	55	55		55	55	55	55	55		55		55	55		55	55	55	55	55	55	55
Dd	Pearson Correlation		– .02 9	– .11 9	.40 2**	.06 4		– .08 3	.02 2	– .20 1	.36 4**	.49 8**		1		.92 6**	– .17 3		.09 2	.05 2	– .04 3	– .06 1	.14 2	– .08 5	.01 0
	Sig. (2- tailed)		.83 3	.38 6	.00 2	.64 5		.54 7	.87 3	.14 1	.00 6	.00 0				.00 0	.20 6		.50 3	.70 4	.75 4	.66 0	.29 9	.53 5	.94 1
	N		55	55	55	55		55	55	55	55	55		55		55	55		55	55	55	55	55	55	55
If	Pearson Correlation		– .06 7	– .09 9	.31 2*	– .03 2		– .13 4	– .05 1	– .19 6	.45 7**	.70 3**		.92 6**		1	– .15 1		.13 8	.01 0	– .04 3	– .08 2	.27 2*	– .07 0	.03 1
	Sig. (2- tailed)		.62 9	.47 3	.02 0	.81 9		.32 9	.71 4	.15 1	.00 0	.00 0		0		.00 .27 2			.31 6	.93 9	.75 5	.55 1	.04 5	.61 0	.82 4
	N		55	55	55	55		55	55	55	55	55		55		55	55		55	55	55	55	55	55	55
ND VI	Pearson Correlation		– .16 6	.40 0**	.01 3	– .06 4		.03 6	– .19 3	– .01 9	– .04 5	.01 7		– .17 3		– .15 1	1		.32 9*	– .82 1**	– .62 0**	– .77 7**	.10 8	– .08 5	– .09 9
	Sig. (2- tailed)		.22 5	.00 2	.92 3	.64 1		.79 3	.15 8	.89 3	.74 6	.90 1		.20 6		.27 2			.01 4	.00 0	.00 0	.00 0	.43 2	.53 7	.47 0
	N		55	55	55	55		55	55	55	55	55		55		55	55		55	55	55	55	55	55	55
MN DWI	Pearson Correlation		– .25 9	.32 8*	– .04 9	– .30 8*		– .19 0	.06 8	– .20 9	.14 2	.15 4		.09 2		.13 8	.32 9*		1	– .68 6**	– .69 7**	– .59 9**	.18 5	– .08 4	– .01 5
	Sig. (2- tailed)		.05 6	.01 4	.72 5	.02 2		.16 4	.62 3	.12 6	.30 1	.26 2		.50 3		.31 6	.01 4			.00 0	.00 0	.00 0	.17 7	.54 0	.91 2
	N		55	55	55	55		55	55	55	55	55		55		55	55		55	55	55	55	55	55	55
ND BI	Pearson Correlation		.25 6	– .48 1**	– .00 4	.18 8		.11 3	.19 9	.16 5	.06 4	– .07 9		.05 2		.01 0	– .82 1**		– .68 6**	1	.81 5**	.83 6**	– .18 2	.06 9	.01 5
	Sig. (2- tailed)		.06 0	.00 0	.97 8	.17 0		.41 2	.14 6	.22 8	.64 0	.56 7		.70 4		.93 9	.00 0		.00 0		.00 0	.00 0	.18 3	.61 9	.91 1
	N		55	55	55	55		55	55	55	55	55		55		55	55		55	55	55	55	55	55	55
BI	Pearson Correlation		.10 2	– .30 0*	.14 3	.25 4		.04 7	.04 3	.17 6	– .13 4	– .03 2		– .04 3		– .04 3	– .62 0**		– .69 7**	.81 5**	1	.80 9**	– .14 3	– .00 1	.12 3
	Sig. (2- tailed)		.45 8	.02 6	.29 7	.06 1		.73 3	.75 4	.19 8	.33 1	.81 5		.75 4		.75 5	.00 0		.00 0	.00 0		.00 0	.29 9	.99 4	.37 3
	N		55	55	55	55		55	55	55	55	55		55		55	55		55	55	55	55	55	55	55
LST	Pearson Correlation		.24 9	– .43 3**	.04 9	.19 4		– .03 7	.17 9	.10 1	– .11 1		– .14 5	– .06 1	– .08 2	– .77 7**	– .59 9**		.83 6**	.80 9**	1	– .19 0	.09 7	.15 9	.096
	Sig. (2- tailed)		.06 6	.00 1	.72 4	.15 6		.78 8	.19 2	.46 2	.41 9		.29 0	.66 0	.55 1	.00 0	.00 0		.00 0	.00 0		.16 5	.48 2	.24 7	.484
	N		55	55	55	55		55	55	55	55		55	55	55	55	55		55	55	55	55	55	55	55
REI	Pearson Correlation		– .20 3	.28 2*	– .23 8	.06 1		– .13 4	– .18 9	– .12 1	.45 6**		.45 5**	.14 2	.27 2*	.10 8	.18 5		– .18 2	– .14 3	– .19 0	1	.09 0	– .10 3	.032
	Sig. (2- tailed)		.13 8	.03 7	.08 0	.65 6		.32 9	.16 8	.38 0	.00 0		.00 0	.29 9	.04 5	.43 2	.17 7		.18 3	.29 9	.16 5		.51 3	.45 2	.816
	N		55	55	55	55		55	55	55	55		55	55	55	55	55		55	55	55	55	55	55	55
SPI	Pearson Correlation		.22 6	.02 5	– .14 7	– .05 4		– .04 1	.13 2	– .49 4**	.00 7		– .08 9	– .08 5	– .07 0	– .08 5	– .08 4		.06 9	– .00 1	.09 7	.09 0	1	.37 1**	.502 **
	Sig. (2- tailed)		.09 7	.85 8	.28 4	.69 5		.76 9	.33 5	.00 0	.96 2		.52 0	.53 5	.61 0	.53 7	.54 0		.61 9	.99 4	.48 2	.51 3		.00 5	.000
	N		55	55	55	55		55	55	55	55		55	55	55	55	55		55	55	55	55	55	55	55
TWI	Pearson Correlation		– .20 8	– .03 2	.11 8	.08 8		– .10 2	.18 9	– .30 5*	– .19 3		– .05 9	.01 0	.03 1	– .09 9	– .01 5		.01 5	.12 3	.15 9	– .10 3	.37 1**	1	.192
	Sig. (2- tailed)		.12 8	.81 7	.39 0	.52 5		.45 8	.16 7	.02 4	.15 9		.66 8	.94 1	.82 4	.47 0	.91 2		.91 1	.37 3	.24 7	.45 2	.00 5		.160
	N		55	55	55	55		55	55	55	55		55	55	55	55	55		55	55	55	55	55	55	55
STI	Pearson Correlation		.31 0*	– .10 1	.11 7	.11 0		– .03 7	.21 7	– .35 8**	.18 1		– .01 7	.22 4	.19 1	– .14 8	– .16 3		.17 3	– .03 7	.09 6	.03 2	.50 2**	.19 2	1
	Sig. (2- tailed)		.02 1	.46 2	.39 3	.42 3		.79 1	.11 2	7	.00		.18 7	.90 0	.10 1	.16 2	.27 9		.23 5	.20 7	.78 8	.48 4	.81 6	.00 0	.16 0
	N		55	55	55	55		55	55	55	55		55	55	55	55	55		55	55	55	55	55	55	55
**. Correlation is significant at 0.01 level (2-tailed)
*. Correlation is significant at 0.05 level (2-tailed)

Source: Authors’ calculation

Table 12: Correlation among the indicators of gully erosion (Bolpur-Sriniketan)

Correlations
Indicators		R	DI	Rn	S	As	Fdir	Fdis	Fl	Fs	Dd	If	ND VI	MN DW I	ND BI	BI	LS T	REI	SPI	TW I	ST I
R	Pearson Correlatio n	1	– .78 1**	.07 5	.19 2	– .03 6	– .14 7	.16 6	.02 2	– .34 3*	– .14 9	– .21 7	.08 0	– .05 2	.04 2	.07 1	– .01 5	– .04 5	.01 5	.10 9	– .05 7
	Sig. (2- tailed)		.00 0	.59 7	.17 3	.80 2	.29 9	.24 0	.87 9	.01 3	.29 3	.12 2	.57 3	.71 6	.76 8	.61 7	.91 8	.74 9	.91 7	.44 1	.68 6
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
DI	Pearson Correlatio n	– .78 1**	1	– .00 5	– .20 4	.04 8	.12 2	.06 4	.20 2	.08 2	– .10 4	– .08 9	.09 8	.00 8	– .05 8	.10 9	– .06 4	– .06 5	.08 6	– .17 3	– .09 2
	Sig. (2- tailed)	.00 0		.97 3	.14 7	.73 4	.38 8	.65 3	.15 2	.56 2	.46 4	.52 8	.49 1	.95 3	.68 4	.44 1	.65 2	.64 6	.54 3	.22 0	.51 4
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
Rn	Pearson Correlatio n	.07 5	– .00 5	1	– .05 0	.09 6	– .13 6	.08 5	.21 3	– .20 9	– .12 7	– .19 0	.11 8	– .05 6	– .06 0	.09 3	– .13 4	– .09 4	.19 9	– .08 6	– .07 7
	Sig. (2- tailed)	.59 7	.97 3		.72 3	.49 6	.33 5	.54 9	.13 0	.13 6	.36 8	.17 7	.40 3	.69 1	.67 5	.51 1	.34 3	.50 9	.15 8	.54 6	.58 8
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
S	Pearson Correlatio n	.19 2	– .20 4	– .05 0	1	– .07 7	– .13 8	– .06 6	– .15 1	– .16 9	– .16 5	– .25 5	– .34 0*	.47 1**	– .25 9	– .31 0*	.06 9	.02 9	– .57 9**	– .02 2	– .10 4
	Sig. (2- tailed)	.17 3	.14 7	.72 3		.58 7	.32 9	.64 1	.28 5	.23 0	.24 3	.06 8	.01 4	.00 0	.06 4	.02 5	.62 5	.84 0	.00 0	.87 7	.46 1
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
As	Pearson Correlatio n	– .03 6	.04 8	.09 6	– .07 7	1	.06 6	– .18 9	– .17 0	.01 5	.02 4	.05 3	.23 0	– .10 9	– .06 6	.08 5	– .11 2	.20 3	.14 6	– .08 1	.00 2
	Sig. (2- tailed)	.80 2	.73 4	.49 6	.58 7		.64 4	.18 0	.22 9	.91 8	.86 5	.70 9	.10 2	.44 2	.64 3	.54 7	.43 1	.14 8	.30 1	.57 0	.98 8
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
Fdir	Pearson Correlatio n	– .14 7	.12 2	– .13 6	– .13 8	.06 6	1	.09 9	.08 9	.21 4	.11 2	.14 4	.07 0	– .18 7	.04 5	.00 5	.02 5	.04 7	– .00 3	.03 1	– .02 4
	Sig. (2- tailed)	.29 9	.38 8	.33 5	.32 9	.64 4		.48 4	.53 0	.12 7	.42 9	.30 9	.62 4	.18 5	.75 3	.97 0	.85 9	.74 2	.98 1	.82 9	.86 6
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
Fdis	Pearson Correlatio n	.16 6	.06 4	.08 5	– .06 6	– .18 9	.09 9	1	.13 7	– .36 2**	– .21 7	– .25 9	.15 9	– .16 7	.09 6	– .02 7	.14 1	– .07 1	.08 7	– .12 7	– .23 0
	Sig. (2- tailed)	.24 0	.65 3	.54 9	.64 1	.18 0	.48 4		.33 1	.00 8	.12 3	.06 4	.26 1	.23 8	.49 9	.85 0	.31 9	.61 9	.53 8	.36 8	.10 2
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
Fl	Pearson Correlatio n	.02 2	.20 2	.21 3	– .15 1	– .17 0	.08 9	.13 7	1	– .10 2	.08 2	.07 1	.12 3	.00 1	– .03 0	.15 5	– .15 5	– .25 9	.11 6	– .33 5*	– .18 7
	Sig. (2- tailed)	.87 9	.15 2	.13 0	.28 5	.22 9	.53 0	.33 1		.47 3	.56 4	.61 6	.38 6	.99 6	.83 3	.27 2	.27 4	.06 4	.41 2	.01 5	.18 4
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
Fs	Pearson Correlatio n	– .34 3*	.08 2	– .20 9	– .16 9	.01 5	.21 4	– .36 2**	– .10 2	1	.04 3	.34 9*	– .31 4*	.02 6	.08 2	– .19 0	.12 0	.19 7	– .22 3	.13 3	.18 2
	Sig. (2- tailed)	.01 3	.56 2	.13 6	.23 0	.91 8	.12 7	.00 8	.47 3		.76 1	.01 1	.02 3	.85 2	.56 1	.17 6	.39 8	.16 2	.11 2	.34 6	.19 6
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
Dd	Pearson Correlatio n	– .14 9	– .10 4	– .12 7	– .16 5	.02 4	.11 2	– .21 7	.08 2	.04 3	1	.92 2**	.05 4	– .12 1	.12 2	.02 5	– .09 8	.02 3	.13 9	– .00 9	– .02 9
	Sig. (2- tailed)	.29 3	.46 4	.36 8	.24 3	.86 5	.42 9	.12 3	.56 4	.76 1		.00 0	.70 5	.39 4	.38 9	.86 1	.48 8	.87 3	.32 7	.94 8	.83 8
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
If	Pearson Correlatio n	– .21 7	– .08 9	– .19 0	– .25 5	.05 3	.14 4	– .25 9	.07 1	.34 9*	.92 2**	1	– .04 6	– .15 4	.17 5	– .05 0	– .03 1	.07 5	.05 9	.03 3	.02 2
	Sig. (2- tailed)	.12 2	.52 8	.17 7	.06 8	.70 9	.30 9	.06 4	.61 6	.01 1	.00 0		.74 9	.27 5	.21 6	.72 4	.82 8	.59 7	.67 7	.81 7	.87 6
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
ND VI	Pearson Correlatio n	.08 0	.09 8	.11 8	– .34 0*	.23 0	.07 0	.15 9	.12 3	– .31 4*	.05 4	– .04 6	1	– .15 8	– .31 4*	.83 5**	– .64 1**	– .29 9*	.51 9**	– .05 5	.11 1
	Sig. (2- tailed)	.57 3	.49 1	.40 3	.01 4	.10 2	.62 4	.26 1	.38 6	.02 3	.70 5	.74 9		.26 3	.02 3	.00 0	.00 0	.03 1	.00 0	.70 1	.43 3
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
MN DWI	Pearson Correlatio n	– .05 2	.00 8	– .05 6	.47 1**	– .10 9	– .18 7	– .16 7	.00 1	.02 6	– .12 1	– .15 4	– .15 8	1	– .80 7**	.18 1	– .44 7**	– .12 8	– .42 4**	– .00 9	.04 5
	Sig. (2- tailed)	.71 6	.95 3	.69 1	.00 0	.44 2	.18 5	.23 8	.99 6	.85 2	.39 4	.27 5	.26 3		.00 0	.20 0	.00 1	.36 7	.00 2	.95 2	.75 0
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
NDB I	Pearson Correlatio n	.04 2	– .05 8	– .06 0	– .25 9	– .06 6	.04 5	.09 6	– .03 0	.08 2	.12 2	.17 5	– .31 4*	– .80 7**	1	– .55 4**	.75 0**	.25 2	.11 9	.03 1	– .09 0
	Sig. (2- tailed)	.76 8	.68 4	.67 5	.06 4	.64 3	.75 3	.49 9	.83 3	.56 1	.38 9	.21 6	.02 3	.00 0		.00 0	.00 0	.07 1	.40 3	.82 8	.52 7
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
BI	Pearson Correlatio n	.07 1	.10 9	.09 3	– .31 0*	.08 5	.00 5	– .02 7	.15 5	– .19 0	.02 5	– .05 0	.83 5**	.18 1	– .55 4**	1	– .82 3**	– .38 8**	.39 3**	.01 1	.13 6
	Sig. (2- tailed)	.61 7	.44 1	.51 1	.02 5	.54 7	.97 0	.85 0	.27 2	.17 6	.86 1	.72 4	.00 0	.20 0	.00 0		.00 0	.00 4	.00 4	.93 8	.33 6
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
LST	Pearson Correlatio n	– .01 5	– .06 4	– .13 4	.06 9	– .11 2	.02 5	.14 1	– .15 5	.12 0	– .09 8	– .03 1	– .64 1**	– .44 7**	.75 0**	– .82 3**	1	.45 0**	– .26 7	.05 2	– .08 2
	Sig. (2- tailed)	.91 8	.65 2	.34 3	.62 5	.43 1	.85 9	.31 9	.27 4	.39 8	.48 8	.82 8	.00 0	.00 1	.00 0	.00 0		.00 1	.05 6	.71 3	.56 1
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
REI	Pearson Correlatio n	– .04 5	– .06 5	– .09 4	.02 9	.20 3	.04 7	– .07 1	– .25 9	.19 7	.02 3	.07 5	– .29 9*	– .12 8	.25 2	– .38 8**	.45 0**	1	– .11 6	– .00 6	.15 4
	Sig. (2- tailed)	.74 9	.64 6	.50 9	.84 0	.14 8	.74 2	.61 9	.06 4	.16 2	.87 3	.59 7	.03 1	.36 7	.07 1	.00 4	.00 1		.41 3	.96 7	.27 5
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
SPI	Pearson Correlatio n	.01 5	.08 6	.19 9	– .57 9**	.14 6	– .00 3	.08 7	.11 6	– .22 3	.13 9	.05 9	.51 9**	– .42 4**	.11 9	.39 3**	– .26 7	– .11 6	1	– .06 4	.19 1
	Sig. (2- tailed)	.91 7	.54 3	.15 8	.00 0	.30 1	.98 1	.53 8	.41 2	.11 2	.32 7	.67 7	.00 0	.00 2	.40 3	.00 4	.05 6	.41 3		.65 1	.17 4
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
TWI	Pearson Correlatio n	.10 9	– .17 3	– .08 6	– .02 2	– .08 1	.03 1	– .12 7	– .33 5*	.13 3	– .00 9	.03 3	– .05 5	– .00 9	.03 1	.01 1	.05 2	– .00 6	– .06 4	1	.16 8
	Sig. (2- tailed)	.44 1	.22 0	.54 6	.87 7	.57 0	.82 9	.36 8	.01 5	.34 6	.94 8	.81 7	.70 1	.95 2	.82 8	.93 8	.71 3	.96 7	.65 1		.23 5
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
STI	Pearson Correlatio n	– .05 7	– .09 2	– .07 7	– .10 4	.00 2	– .02 4	– .23 0	– .18 7	.18 2	– .02 9	.02 2	.11 1	.04 5	– .09 0	.13 6	– .08 2	.15 4	.19 1	.16 8	1
	Sig. (2- tailed)	.68 6	.51 4	.58 8	.46 1	.98 8	.86 6	.10 2	.18 4	.19 6	.83 8	.87 6	.43 3	.75 0	.52 7	.33 6	.56 1	.27 5	.17 4	.23 5
	N	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52	52
**. Correlation is significant at 0.01 level (2-tailed) *. Correlation is significant at 0.05 level (2-tailed)

Source: Authors’ calculation

Table 13: Covariances are among the indicators of gully erosion (Rampurhat-I)

Covariances
Variabl es	DI	Rn	S	As	Fdir	Fdis	Fl	Fs	Dd	If	NDV I	MND WI	NDBI	BI	LST	REI	SPI	TWI	STI
R	–	–	–	–	7.83	2.658	0.812	–	–	–	–	-0.15	0.219	21.521	2.115	–	0.228	-382.79	0.864
	0.669	0.01**	0.886	73.833				0.294	0.012	0.044	0.148					0.088
DI	0.059	–	–	-0.115	-1.026	-0.268	–	0.011	–	–	0.011	0.006*	–	–	–	0.004	0.001	-1.887	-0.009
	**	0.002	0.026				0.042*	*	0.002	0.002		*	0.013	2.025**	0.117	**	*
													*		*
Rn	–	0.008	0.036	0.119	-0.037	-0.016	0.001	–	0.002	0.002	0*	0	0	0.357	0.005	–	–	2.576	0.004
	0.002							0.002		**						0.001	0.002
S	–	0.036	4.059	16.629	–	0.928	0.01	–	0.007	–	–	-0.047	0.043	14.239	0.437	0.007	–	42.894	0.082
	0.026				11.989			0.218		0.006	0.015		*				0.014
As	–	0.119	16.62	5905.6	25.376	19.80	-4.072	-6.56	–	-0.9	0.325	-1.111	0.981	100.59	-3.18	-0.59	–	-1906.955	-1.033
	0.115		9	34		9			0.343					8			0.414
Fdir	–	–	–	25.376	817.04	-10.21	0.524	-2.34	0.034	–	–	0.147	0.642	34.279	5.704	–	0.503	1312.732	2.275
	1.026	0.037	11.98		4					0.126	0.646					0.308
			9
Fdis	–	–	0.928	19.809	-10.21	11.25	-0.021	–	–	–	–	-0.053	0.063	16.423	0.379	–	-0.22	–	-0.44*
	0.268	0.016				4		0.186	0.036	0.058	0.007					0.023		248.291**
		*
Fl	–	0.001	0.01	-4.072	0.524	-0.021	0.359	0.177	0.012	0.024	–	0.006	0.004	-2.223	–	0.016	0.001	-28.055	0.04
	0.042	*							**	**	0.003				0.074		**
											**
Fs	0.011	–	–	-6.56	-2.34	-0.186	0.177	0.334	0.015	0.036	0.001	0.007	–	-0.518	–	0.015	–	-8.318	-0.004
		0.002	0.218					**		**	**		0.005		0.094		0.007
																	**
Dd	– 0.002	0.002	0.007 **	-0.343	0.034	-0.036	0.012	0.015 **	0.003 **	0.004	– 0.001 **	0	0	-0.065	– 0.004	0	– 0.001	0.133	0.004
If	– 0.002	0.002	– 0.006 *	-0.9	-0.126	-0.058	0.024	0.036 **	0.004 **	0.008 **	– 0.002	0.001	0	-0.105	– 0.008	0.001	– 0.001 *	0.652	0.006
NDVI	0.011	0**	– 0.015	0.325	-0.646	-0.007	-0.003	0.001	– 0.001	– 0.002	0.014	0.003	– 0.011 *	– 2.016**	– 0.102 **	0.001 **	– 0.001	-2.827	-0.006
MND WI	0.006	0*	– 0.047	-1.111*	0.147	-0.053	0.006	0.007	0	0.001	0.003	0.006*	– 0.006	– 1.472**	– 0.051 **	0.001 **	– 0.001	-0.282	-0.005
NDBI	– 0.013	0**	0.043	0.981	0.642	0.063	0.004	– 0.005	0	0	– 0.011	– 0.006* *	0.013 **	2.559	0.106 **	– 0.001 **	0.001	0.422	0.007
BI	– 2.025	0.357 *	14.23 9	100.59 8	34.279	16.42 3	-2.223	– 0.518	– 0.065	– 0.105	– 2.016	– 1.472* *	2.559 **	771.71 4**	25.12 5	– 0.226 **	– 0.004	827	-0.378
LST	– 0.117	0.005 **	0.437	-3.18	5.704	0.379	-0.074	– 0.094	– 0.004	– 0.008	– 0.102	– 0.051* *	0.106 **	25.125* *	1.248 **	– 0.012	0.014	43.097	0.04
REI	0.004	– 0.001 *	0.007	-0.59	-0.308	-0.023	0.016	0.015 **	0**	0.001	0.001 *	0.001	– 0.001	-0.226	– 0.012	0.003	0.001	-1.437	0.001
SPI	0.001	– 0.002	– 0.014	-0.414	0.503	-0.22	0.001* *	– 0.007	– 0.001	– 0.001	– 0.001	-0.001	0.001	-0.004	0.014	0.001	0.018	11.968	0.024**
TWI	– 1.887	2.576	42.89 4	– 1906.9 55	1312.7 32	– 248.2 91	– 28.05 5*	– 8.318	0.133	0.652	– 2.827	-0.282	0.422	827	43.09 7	– 1.437	11.96 8	59021.39 9**	17.117
STI	– 0.009 *	0.004	0.082	-1.033	2.275	-0.44	0.04**	– 0.004	0.004	0.006	– 0.006	-0.005	0.007	-0.378	0.04	0.001	0.024	17.117**	0.135
**. Significant at 0.01 level (2-tailed) *. Significant at 0.05 level (2-tailed)

Source: Authors’ calculation

Table 14: Covariances are among the indicators of gully erosion (Bolpur-Sriniketan)

Covariances
Variabl es	R	DI	Rn	S	As	Fdir	Fdis	Fl	Fs	Dd	If	NDV I	MND WI	NDB I	BI	LST	REI	SPI	TWI	STI
	1	– .781* *	0.07 5	0.192	-0.036	-0.147	0.166	0.02 2	-.343*	– 0.149	– 0.217	0.08	-0.052	0.042	0.071	– 0.015	– 0.045	0.015	0.109	– 0.057
R	69.76	–	0.05	2.047	-24.506	-33.16	5.189	0.01	–	–	-0.81	0.098	-0.064	0.054	0.06	–	–	0.013	461.351	–
	6	6.97*							0.381	1.226						0.153	0.134			0.206
		*							*
DI	-6.97**	1.142	0	–	4.246	3.527	0.255	0.01	0.012	–	–	0.015	0.001	–	0.012	–	–	0.01	-93.481	–
				0.278				2		0.109	0.043			0.009		0.085	0.025			0.042
Rn	0.05	0	0.00	–	0.637	-0.295	0.026	0.00	–	-0.01	–	0.001	-0.001	–	0.001	–	–	0.002	-3.476	–
			6	0.005				1	0.002		0.007			0.001		0.013	0.003			0.003
S	2.047	–	–	1.631	-8.12	-4.768	–	-0.01	–	–	–	–	0.089**	–	-0.04*	0.111	0.013	–	-14.166	–
		0.278	0.00				0.316		0.029	0.208	0.145	0.064		0.051				0.079		0.057
			5									*						**
As	–	4.246	0.63	-8.12	6800.3	146.3	–	-0.75	0.161	1.961	1.954	2.789	-1.328	–	0.716	-11.5	5.912	1.285	–	0.073
	24.50		7		27	82	58.28							0.832					3364.62
	6						9												4
Fdir	-33.16	3.527	– 0.29 5	– 4.768	146.38 2	730.1 02	10.03 4	0.12 9	0.771	2.993	1.737	0.277	-0.746	0.185	0.015	0.849	0.446	-0.01	420.095	– 0.278
Fdis	5.189	0.255	0.02 6	– 0.316	-58.289	10.03 4	14.02 2	0.02 8	– 0.181 **	– 0.801	– 0.432	0.088	-0.092	0.055	-0.01	0.659	– 0.093	0.035	-241.48	-0.37
Fl	0.01	0.012	0.00 1	-0.01	-0.75	0.129	0.028	0.00 3	– 0.001	0.004	0.002	0.001	0	0	0.001	-0.01	– 0.005	0.001	-9.094*	– 0.004
Fs	– 0.381*	0.012	– 0.00 2	– 0.029	0.161	0.771	– 0.181 **	– 0.00 1	0.018	0.006	0.021 *	– 0.006 *	0.001	0.002	– 0.003	0.02	0.009	– 0.003	8.968	0.01
Dd	-1.226	– 0.109	– 0.01	– 0.208	1.961	2.993	– 0.801	0.00 4	0.006	0.975	0.407 **	0.008	-0.018	0.018	0.003	– 0.121	0.008	0.015	-4.616	– 0.012
If	-0.81	– 0.043	– 0.00 7	– 0.145	1.954	1.737	– 0.432	0.00 2	0.021 *	0.407 **	0.199	– 0.003	-0.01	0.012	– 0.002	– 0.017	0.012	0.003	7.421	0.004
NDVI	0.098	0.015	0.00 1	– 0.064 *	2.789	0.277	0.088	0.00 1	– 0.006 *	0.008	– 0.003	0.022	-0.003	– 0.007 *	0.012	– 0.118 **	– 0.016 *	0.008 **	-4.069	0.007
MNDW I	-0.064	0.001	– 0.00 1	0.089 **	-1.328	-0.746	– 0.092	0	0.001	– 0.018	-0.01	– 0.003	0.022	– 0.018 **	0.003	– 0.082 **	– 0.007	– 0.007 **	-0.641	0.003
NDBI	0.054	– 0.009	– 0.00 1	– 0.051	-0.832	0.185	0.055	0	0.002	0.018	0.012	– 0.007 *	– 0.018**	0.023	– 0.009 **	0.144 **	0.014	0.002	2.392	– 0.006
BI	0.06	0.012	0.00 1	-0.04*	0.716	0.015	-0.01	0.00 1	– 0.003	0.003	– 0.002	0.012 **	0.003	– 0.009 **	0.01	– 0.105 **	– 0.014 **	0.004 **	0.568	0.006
LST	-0.153	– 0.085	– 0.01 3	0.111	-11.5	0.849	0.659	-0.01	0.02	– 0.121	– 0.017	– 0.118 **	– 0.082**	0.144 **	– 0.105 **	1.561	0.198 **	– 0.036	33.015	– 0.044
REI	-0.134	– 0.025	– 0.00 3	0.013	5.912	0.446	– 0.093	– 0.00 5	0.009	0.008	0.012	– 0.016 *	-0.007	0.014	– 0.014 **	0.198 **	0.124	– 0.004	-1.045	0.023
SPI	0.013	0.01	0.00 2	– 0.079 **	1.285	-0.01	0.035	0.00 1	– 0.003	0.015	0.003	0.008 **	– 0.007**	0.002	0.004 **	– 0.036	– 0.004	0.011	-3.467	0.009
TWI	461.3 51	– 93.48 1	– 3.47 6	– 14.16 6	– 3364.6 24	420.0 95	– 241.4 8	– 9.09 4*	8.968	– 4.616	7.421	– 4.069	-0.641	2.392	0.568	33.01 5	– 1.045	– 3.467	255915. 13	36.48 6
STI	-0.206	– 0.042	– 0.00 3	– 0.057	0.073	-0.278	-0.37	– 0.00 4	0.01	– 0.012	0.004	0.007	0.003	– 0.006	0.006	– 0.044	0.023	0.009	36.486	0.185
**. Significant at 0.01 level (2-tailed) *. Significant at 0.05 level (2-tailed)

Source: Authors’ calculation