Geomorphic character and dynamics of gully morphology, erosion and management in laterite Terrain: few observations from Dwarka – Brahmani Interfluve, Eastern India

ABSTRACT The most intensified form of soil erosion is the gully which is an extreme figure of land degradation in India, and alongside gully erosion signifies instability in the fragile landscape. The present study encompasses the lateritic interfluve and badlands of Dwarka – Brahmani River Basin (eastern India). The geomorphic research reveals a variable range of annual erosion rates (8.12–24.01 kg m−2 y−1) at watershed scale (i.e., three catchments of permanent gullies) using models and field measured data. It is found that the main cause of gully formation is too much runoff water at a certain location of slope – a threshold condition that may be brought about by external factors (land cover) or internal factors (slope). It is found that in the event-based rainfall range of 42 mm to 137.2 mm the gullies can yield runoff of 40.02 mm to 118.0 mm in excess moisture condition of monsoon. Intense rainfall is the primary trigger, but the local conditions such as slope morphometry, land use, barren soil cover and soil-plant characteristics control the triggering of gully erosion. The potential erosion map of area depicts annual erosion rate beyond the soil tolerance limit (T-value – 1.0 kg m−2 y−1). Finally it is suggested that reduction of runoff discharge, channel grade, vegetative measures and structural control of gully headcut erosion and sedimentation are the key procedures of erosion protection in the laterite terrain.

The population density of India increases from 382 person km −2 (2011) to 464 person km −2 (2020). Land degradation is greatly aggravated in recent decades because of increasing population of India at the rate of about 1.8 percent requiring marginal areas to be brought under the plough to meet the growing food demand (Aulakh & Sidhu, 2015). It is estimated that 120.72 million ha area is affected by various forms of land degradation and desertification in India with water erosion being chief contributor (68.4%) . It was found that ravines along the banks of the Yamuna, Chambal, Mahi, Tapti and Krishna Rivers, revealed soil losses exceeding 18 t ha −1 y −1 and the erosion rates on the alluvial Indo-Gangetic Plains of Punjab, Haryana, Uttar Pradesh, Bihar and West Bengal were 5 to 10 t ha −1 y −1 (Singh et al., 1992). In this country major rainfed crops suffer an annual production loss of 13.4 Mt due to water erosion which amounts to a loss of Rs. 205.32 billion in monetary terms (Bawa, 2017). The reservoirs of India are losing about 1.3 billion m 3 of storage capacity each year due to soil erosion and siltation. That should be alarming enough for everyone as at today's rates creation of 1.3 billion m 3 storage capacity would cost Rs. 1448 crores (Thakkar & Bhattacharyya, 2006). Now, the major point of research interest is the laterite terrain of West Bengal, eastern India (known as Rarh Bengal, i.e., the land of red soil), which is severely dissected by the dense network of rills and gullies (Ghosh & Guchhait, 2017), developing miniature forms of badland topography. The lateritic Rarh region and plateau fringe (districts of Purulia, Bankura, Paschim Barddhaman and Paschim Medinipur) show lower T value (soli-tolerance level) ranging from 2.5 to 5.0 Mg ha −1 y −1 (Lenka et al., 2014;Mondal & Sharda, 2011). Here, realistic assessment of gully erosion risk or soil loss rate thus constitutes the first step for understanding the ground reality of erosion and raising awareness among governmental and other stakeholders in a given region to adopt appropriate strategies for erosion protection. Before taking any erosion protection measures, the estimation of channel and slope erosion (Strahler, 1964) at plot to basin scale is the fundamental step towards achieving soil conservation and sustainable development. In this regard, this study can give light on the hydro-geomorphic aspects of gully erosion assessment and erosion modelling using minimal data inputs and measured plots at basin scale to take fundamental steps of erosion management in the lateritic badlands. Therefore, three major objectives of the study are set forth as follows - (1) To understand the gully morphology, erosion and threshold in laterite terrain, (2) To estimate channel and slope erosion using model and field experimental database; and (3) To suggest erosion protection strategies for the study area.

Geographical settings of study area
The geomorphic unit of study is recognized as an interfluve in between Brahmani (north) and Dwarka (south) rivers (encompassed by 24° 20´ N to 23° 40´ N, and 87° 26´ E to 88° 21´ E) ( Figure 1). Dwarka River Basin (2,978 km 2 ) is a sub-basin of Mayurakshi River Basin and Brahmani River (1,139 km 2 ) is a sub-basin of Dwarka River. The selected study area of Dwarka -Brahmani interfluve (about 176 km 2 , encompassed by 24°08´N to 24°14´ N and 87°38´ E to 87°44´ E), covers Shikaripara block (Dumka, Jharkhand), Rampurhat I and Nalhati I blocks (Birbhum, West Bengal) ( Figure 2). Geologically, the present research work deals with the contiguous unit between Rajmahal Basalt Traps (RBT) (Early Cretaceous origin) and the Bengal Basin which exhibits shallow Quaternary alluvium deposits and palaeogenesis of the deep weathering profiles under intense tropical wet -dry palaeoclimate on the basaltic surface to form hard ferruginous crust, i.e., Ferricrete (Palaeogene -Early Pleistocene) (Ghosh & Guchhait, 2015). The climate of this region has been identified as sub-humid and sub-tropical monsoon type, receiving mean annual rainfall of 1300 to 1437 mm. The peak monsoon and cyclonic rainfall intensity of 21.51 mm hr −1 (minimum) to 25.51 mm hr −1 (maximum) which is the most powerful climate factor to develop this lateritic badlands. Generally, the thin solum is loamyskeletal and hypothermic (sandy loam to sandy clay loam texture) in nature developing on the barren lateritic wastelands with sparse bushy vegetation and grass. The dark reddish to brown coloured sandy clay loam of 0-16 cm (A horizon, maximum grass root zone) is developed over the fragmented secondary Annual soil erosion is taking place at the rate of 16.35 t ha −1 y −1 . Narayana and Babu (1983) 2 Indo-Gangetic Plains of Punjab, Haryana, Uttar Pradesh, Bihar and West Bengal are affected by erosion rate of 5 to 10 t ha −1 y −1 . Singh et al. (1992) 3 About 20% of India's existing reservoirs will have lost 50% of their previous storage capacity due to soil loss and siltation. Kothyari (1996) 4 3.975 million ha of wastelands are severely affected by gullies and ravines. Yadav and Bhushan (2012) 5 Due to siltation, India is losing about 1.3 billion m 3 of storage capacity each year and to create this storage capacity India will require Rs. 1448 crores. Thakkar and Bhattacharyya (2006) 6 Paddy is the most affected among all crops in terms of both productions 4.3 million tonne and monetary loss of Rs. 24.4 billion. Sharda et al. (2010) 7 The Lower Gangetic Plain and eastern part of Chotanagpur Plateau has soil loss tolerance level of 2.5 to 12.5 t ha −1 y −1 . Mondal and Sharda (2011) 8 India suffers an annual loss of 13.3 million tonne in production of cereals, oilseeds and pulses due to water erosion.  9 About 69.5% area of India has soil loss tolerance limit of <10 t ha −1 y −1 .
Sharda and Dogra (2013) 10 About 5.4 million tone of fertilizer worth US $ 245 million is washed away by water erosion. Lenka et al. (2014) 11 Erosion escalates the siltation rate of reservoirs in India -Maithon (1.076 mm y −1 ), Panchet (0.631 mm y −1 ), Tilaiya (2.792 mm y −1 ), Tenughat (0.716 mm y −1 ), Durgapur barrage (0.042 v), Kangsabati (0.752 mm y −1 ) and Massanjore (0.557 mm y −1 ) Central Water Commission (2015) 12 The soil pool loses 110 Mt Carbon into the atmosphere due to soil erosion. It is projected that 1% increase in rainfall intensity may increase the rainfall Erosivity by 2-6%. Annual loss due to soil degradation ranges from Rs. 89 to 232 billion. Bawa (2017) laterites. These ferruginous soils has weak fine crumb and granular structure (slightly hard, friable and slightly sticky), 2-5 mm size of manganese nodules, > 2 mm size of ferruginous nodules with goethite cortex, 30 to 40% travels and pebbles, excessive drained surface and p H of 5.4-5.7. The natural vegetation of the study area belongs to the tropical moist and dry deciduous type with few  evergreen types. The observed natural vegetation species are: Babul (Acacia nilotica), Bel (Aegle marmelos), Behara (Terminalia belerica), Sal (Shorea robusta), Mahua (Madhuca indica), Khair (Acacia catechu), Khajur (Phoenix sylvestris) and Jamun (Syzygium cumini) etc. Though one upon a time the most of the region was covered under thick forest, mainly Sal (Shorea robusta), due to encroachment of stone crushers, mining and agriculture the forests are fragmented and vanished from some places. Now, the natural green vegetation and grassland cover an area of 58.83 km 2 . Stone quarries, roads, built-up area and barren fallow land (including basalt exposure and non-arable land) cover an area of 79.56 km 2 which acts as main source region of surface runoff in the peak monsoon ( Figure 3).

Secondary data collection
The key sources of main secondary data are regional soil report, geology report and other physical environmental report published by NBSS and LUP (National Bureau of Soil Service and Land Use Planning), Census of India, district gazetteer, official websites of IMD (Indian Meterological Department) Pune and Kolkata, Irrigation and Waterways Dept. of Govt. of West Bengal (IWD), Geological Survey of India (GSI), related e-books and e-journals. The topographical sheets of Survey of India (72 P/12/ NE, R.F. 1:25,000 and 72 P/12, R.F. 1:50,000), District Resource Map of Geological Survey of India, District Planning Map of NATMO (National Atlas Thematic Mapping Organization) and Block map of Census of India are most important sources of spatial information. Landsat TM (30 m resolution) images are downloaded from the website of Global Land Cove Facility (GLCF) and and ASTER (Advanced Spaceborne Thermal Emission and Reflection Radiometer, 30 m resolution) elevation data are downloaded from the websites of GLCF and Consortium for Spatial Information (CGIAR-CSI). The spatial information is stored in Geographic Information System (GIS) and the thematic maps are prepared using GIS software (ArcGIS 9.2, Erdas Image 9.1 and MapInfo Professional 11.5).
In this case we have gathered the daily, monthly and annual rainfall data from three IWD (Irrigation and Waterways Department, Government of Wes Bengal) raingauge stations at Nalhati (24°17´25´´N, 87°49´44´´E), Rampurhat (24°10´13´´N, 87°46´50 ´´E) and Mollarpur (24°04´35´´N,87°42´36´´E) which are situated at eastern part of study area, having areal distance of 18 to 25 km. The calculated mean annual rainfall for this region is 1510 mm in 2016 (maximum intensity of erosive rain is 25.21 mm hr −1 ) and the per day rainfall amount is 17.48 mm, considering total rainfall and rainy days in a year. From the database we have selected the rainfall events of greater than 40 mm per day and these events vary from 42 mm to 137.2 mm of rainfall in 2016 for the SCS-CN (Soil Conservation Service -Curve Number) modelling.

Experimental design of erosion plot and sedimentation measurement at field
In this study the water erosion of a basin will be recognized as taking two basically different forms: (a) slope erosion (dominated by rill and inter-rill erosion at hillslope) and (b) channel erosion (cutting-away of bed and banks of clearly marked channel) (Strahler, 1964). The spatial scale to study erosion processes is selected here as plot-scale (10 to 100 m 2 ) and field scale (100 to 10,000 m 2 ). Since it is only possible to take measurements at specific points in the landscape and it is important that these area representative of the catchment as high erosion prone zone (where maximum erosion is observed). From the field survey it is observed that except channel erosion gully head slope (average 7°34 ´) is the key pathways of sediment transport to the main gully channel. In this lateritic terrain the high erosion risk catchment of gully (selecting three catchments of gullies) is firstly selected and it has well defined basin area (about 1,09,250 to 2,16,050 m 2 ) and dense network of gullies (7.57 to 8.33 km km −2 ) ( Figure 4). Firstly, eighteen gully heads of three basins (selected randomly within seventeen basins at study area, based on high drainage density of greater than 7.5 km km −2 ) were identified and then eighteen gully head slope elements (considering 2 m width of slope strip to incorporate soil -land use parameters) were selected. The steepness of hillslope was measured using Leica Sprinter 150 m digital levelling instruments (accuracy -± 0.7 mm of the 250 m distance) and other parameters of erosion model were estimated in the recurrent field survey (2016)(2017) and guided values of erosion model (Morgan, 2005;Renard et al., 1997;Wischmeier & Smith, 1978). The steepness of slope elements varies from 3°45´to 11°06´, whereas slope length varies from 22.1 m to 106.8 m. Maintaining a certain distance (1.5 to 2 m) from active gully head, eighteen check dams (used as sedimentation pits) were developed (denoting Dam 1 to Dam 18) at the base (i.e., gully floor) of representative slope elements to trap eroded sediments coming from upslope in a year (2016-2017) ( Figure 5).
Development of check dams served two essential purposes -(1) soil conservation approach (measuring sedimentation and checking gully erosion), and (2) soil erosion estimation approach (comparing the model-based erosion rates against observed erosion rate at dam sites). Following the shape of gully channel it was decided to built v-shape design using sand, cement and laterite boulders of irregular shape and the gap between gully head base and dam was used as sedimentation pits to collect eroded materials. The dams had height range of 40-55 cm and width range of 92-190 cm. The dams were developed in January, 2016 with the help of local manpower and resources. Then after one year of observation the sedimentation was measured in January, 2017 and the mass volume was measured as multiplying the area of sedimentation behind dam and mean depth of sedimentation at eighteen dam sites. The bulk density of eroded materials was calculated at laboratory (mean bulk density of materials is 1.717 gm cm −3 ) and the mass weight of sedimentation materials was measured by multiplying the volume of mass (unit in kg) by bulk density.

Gully morphometry
In this study the levelling survey was conducted to assess the morphology (width, depth, shape, basin area and slope) in three gully catchments using Leica Sprinter 150 m (height accuracy ± 1.5 mm and distance accuracy ± 1 mm) and Garmin GPS (horizontal accuracy ±3 m). To understand the changing transverse profile of gully along the channel from head to mouth, several cross-sectional surveys were taken and then, plotting successive cross-profiles along the downstream gully floor it can be learn that shape of gully can reflect the stages of development and dominancy of processes over underlying structure. Based on field studies the morphological distinction among gullies were examined and various classification of gullies were prepared -(1) physical based classification, (2) hydro-geomorphic classification and (3) geometric shape classification ( Figure 6) .

Estimation gully erosion
The role of terrain has a definite role in the development of badlands. The susceptibility of soil to water erosion (i.e., K factor of RUSLE) is a crucial factor to initiate gully head in the laterite terrain. The K factor is a composite function of percentage of sand, silt, clay and organic matter, including soil cohesion and permeability. In this case the relation between gully density (km km −2 ) and Erodibility Index (K) are assessed in seventeen gully catchments. To estimate channel erosion a simple approach of transverse cross-profile and GPS-based transect survey (from gully head to outlet) was performed in three sample catchments, viz. gully catchment 1, 2 and 3 ( Figure 7). The crossprofile wise erosion estimation includes both channel and side-wall erosion in a gully. To estimate the potential sites of bank failure a survey (3D area measurement tool) was performed using Leica DistoS910 Laser Distance Meter (height accuracy ± 0.05 mm, distance accuracy ± 1.0 mm). The main task is to calculate the eroded area between two transect profiles and then to estimate volume, multiplying the area with length between two transect ( Figure 8). The area of profile is measured using the trapezoidal equation or cone. To estimate loss of earth materials between two transects can be estimated using the following equation (Figure 9).

Loss of Earth materials kg ð Þ ¼ Bulk density of earth materials=volume of each segment
(1) The bulk density of laterite profile (up to 1.45 m) was estimated in the laboratory taking different samples from different depths and the calculated average bulk density is 2.205 gm cm −3 .
Garmin 76csx receiver and Google Earth Pro are used to assess the retreat rate of permanent gully head in between 2010 and 2018 (thirteen GPS points). The geographic locations are recorded from the GPS receiver and these locations are plotted in Google Map to estimate the recession of head as metre per year. The volumetric retreat rate of gully head is estimated by Poesen (2018). It is found that gully headcut retreat rate (GHR) is a function of contributing area of runoff (A R ) and normal rainy day per year (RDN) (Poesen, 2018;Vanmaercke et al., 2015). The expression of equation is mentioned as follows: GHR ¼ 0:001 : AR 0:52 : RDN 4:97 (2)

Geomorphic threshold model
Channel initiation by hydro-geomorphic processes has been viewed as a threshold phenomenon related to size of contributing area (A) and its slope (S) (Schumm, 1979;Torri & Poesen, 2014). There are wide ranges of threshold conditions or limits (viz., thresholds of hydraulic, rainfall, topography, lithology and land use -land cover control etc.) which are responsible for the initiation of gullies in different environments (Begin & Schumm, 1984;Dong et al., 2013;Moeyersons, 2003;Montgomery & Dietrich, 2004;Morgan & Mngomezulu, 2003;Patton & Schumm, 1975;Poesen et al., 2003;Samni et al., 2009;Torri & Poesen, 2014;Valentin et al., 2005;Vandaele et al., 1996). It is hypothesized that the gullies over the laterites develop when the geomorphic thresholds (may be extrinsic or intrinsic) are transgressed due to either a decrease in the resistance of the materials (i.e., erodibility) or an increase in the erosivity of the runoff or both. In terms of identifying the geomorphic thresholds in gully initiation, the present experimental work includes the 118 gully heads (both valley-floor and valley-side gullies). The relation between critical valley slope and drainage basin area (S = a A −b , where a = coefficient and b = exponent of relative area) is used as a predictive model to locate those areas of instability within alluvial valleys where gullies will form (Ghosh & Guchhait, 2017).
A threshold line is drawn through the lower limit of scatter of the points and this line represents, for a given area, a critical value for valley slope above which entrenchment of the laterite should occur.

SCS-CN based rainfall-runoff model
The Soil Conservation Service Curve Number (SCS-CN) method is used for quick and accurate estimation of surface runoff in any storm event in the un-gauged watersheds of India (Gajbhiye et al., 2014;Mishra & Singh, 2003;Patil et al., 2008;Srivastava & Imtiyaz, 2016). This model was used in three selected catchments of gullies using the database of daily rainfall event of 2017. The SCS-CN rainfall-runoff model is based on the water balance equation and two fundamental hypotheses. This model includes different parameters, viz., direct surface runoff (Q), total rainfall (P), the amount of actual infiltration (F), initial abstraction (I a ) and the potential maximum retention (S). Thus, SCS-CN method consists of following equation (Chow et al., 1988;Gajbhiye et al., 2014;Mishra et al., 2006;Srivastava & Imtiyaz, 2016): In practice, S is derived from a mapping equation expressed in terms of the curve number (CN) (Mishra et al., 2006): The CN is a dimensionless value (30 to 100) to reflect the runoff yield for a specific event-based rainfall in respect of land use -land cover (derived from supervised classification of Landsat TM images), soil type, hydrologic condition and antecedent moisture condition (AMC). To derive the average CN-values for AMC II mathematically from the rainfall-runoff data of a gauged watershed, it is suggested S-computation using the following equation (Mishra et al., 2006): For dry conditions (AMC I) and wet condition (AMC III), the equivalent curve numbers can be computed by (Chow et al., 1988) Now, it is essential to assess the role of catchment runoff in the flow erosion in the lateritic badlands. In this case critical tangential shear stress (T s , pascal) is estimated at thirty one sites of permanent gully heads. The expression of the equation is mentioned as follow (Leonard & Richard, 2014;Sidorchuk, 2020): where, g is acceleration due to gravity, p is water mass density, d is flow depth and S f is friction slope.
To estimate transport capacity of runoff (H, kg m −2 ), Morgan (2005) has developed an equation at slope scale. It is a function of soil cohesion (Z), runoff (QR), sine of slope (S) and ground canopy cover (GC, percentage). The equation is expressed as follows (Morgan, 2005):

Estimating slope erosion by RSULE Model
The Revised Universal Soil Loss Equation (RSULE) (Renard et al., 1997) is an empirical equation for predicting long-term average soil erosion (rill and interrill erosion) from agricultural field under specific cropping and management practice. Smith (1999), Sovrin (2003) have successfully applied RUSLE to assess erosion rate in different environmental settings and they have found the suitability and effectiveness of RUSLE in comparison to other models. The applied version of RUSLE is mentioned as follows (Renard et al., 1997) ( Table 2 and Figure 10).
where, • A p is the computed soil loss per unit area (t ha −1 y −1 or kg m −2 y −1 ); it can transformed into SI unit • R, the rainfall and runoff factor, is the number of rainfall erosion index units, i.e., EI 30 • K, the soil erodibility factor • L, the slope-length factor • S, the slope-steepness factor • C, the cover and management factor • P, the support practice factor

Model validation and effectiveness coefficient
The model efficiency coefficient (MEC), firstly proposed by Nash and Sutcliffe (1970), is now increasingly used an alternative to the correlation coefficient to express the performance of model (Morgan, 2011). Generally, a MEC value of greater than 0.5 is considered that the model performs satisfactorily in the region, and one should not expect values to exceed 0.7 (Morgan, 2005(Morgan, , 2011. where X obs is the observed value, X pred is the value predicted by the model andX' obs is the mean of a set of observed values. A model "effectiveness coefficient" (E c ) was defined by Nearing (2000) for studies undertaken on large numbers of prediction versus measured data comparisons. The relative difference (R diff ) between predicated and measured values are calculated and then a particular set of conditions that 95% of the values for differences in erosion (fall within a certain range) is calculated.
Relative difference values (Y-axis) are plotted against measured values (X-axis) to get a trend in the scatters.
It is defined E c as the fraction of simulation model predictions for which a model is effective in predicting the measured erosion, using the acceptance criteria. Using the 95% occurrence intervals from the replicated erosion data would result in a value, E c (a = 0.05) . The value of E c (a = 0.05) signifies that the percentage of the difference between measured and predicted soil loss fell within the expected range of difference for two measured data points within the same population (Nearing, 2000). In this case Chisquare test and t-test are also applied to compare predicated and observed erosion rates.

Categorization of gully in laterite terrain
A simple, physically based classification distinguishes gullies as (a) ephemeral gullies (i.e., these are impermanent channels that are obliterated periodically by cultivation, e.g., deep tillage or land-levelling operation, or natural processes, i.e., deposition) and (b) permanent gullies (i.e., these are deeply incised channels that have cross sections permanently recognizable without flowing water and have identifiable banks) ( Figure 11). The present classification of gullies was based on field observations and aerial view, and it includes the gully planform, gully side morphology and forms of the longitudinal and transverse profiles. Based on 118 gully head samples of study area two classification schemes have been formulated -(1) gully classes based hydro-geometry and (2) gully classed based on transverse shape.

Gully classes based on hydro-geometry
The simplest classification system is based gully depth (signifying the extent of incision and stage of gully development): < 1-2 m depth is recognized as small gully, 2-5 m depth is recognized as medium gully and > 5 m depth is identified as large gully. Flow rate measurement is an effective way to classify gullyincreasing flow convergence and high level of erosion. So, it is identified that < 0.3 m 3 s −1 discharge rate is associated with small gully, 0.3-2.0 m 3 s −1 is associated with medium gully and > 2.0 m 3 s −1 is recognized as large gully in the study area. The annual runoff yield of < 100-350 mm is associated with small gully, whereas 350-600 mm and > 600 mm runoff yield are associated with medium and large gully. From the result (Table  4.1) it is found that 38.72% gullies are associated with small gullies which have < 2400 m 2 basin area, < 3-5 m width and < 1-2 m depth. About 48.25% of gullies are assonated with medium gullies which have basin area of 2400-7200 m 2 , 5-10 m width and 2-5 m depth. Out of 118 samples, 13.03% of gullies are recognized as large gully which have basin area of > 7200 m 2 , > 10 m width and > 5 m depth (Table 3).

Gully classes based on shape
Three principal transverse gully shapes are identified in the study area and these shapes have deep relation with erosion dominancy, stage of development and erodible soil layers. Gullies are classified according to shape of cross-section ( Figure 12): • U-shaped Gullies -These gullies are formed where both topsoil and subsoil have the same OM is the percentage of organic matter in soil, M is the particle size parameter, s is the soil structure code and p is permeability code recommended by (Wischmeier & Smith, 1978) Bayramov L is the slope length (m) and θ is slope steepness in percent.
Rahaman et al. (2015); resistance against erosion. Since the subsoil is eroded as easily as the topsoil, nearly vertical walls are developed on each side of the gully. • V-shaped Gullies -These gullies develop where the subsoil has more resistance than topsoil against erosion. This is the most common gully form. • Trapezoidal Gullies -These are formed where the gully bottom is made of more resistant material than the topsoil and subsoil because the erosion rate along the gully bank is greater than along the bottom.
In the study area, out of total 118 gullies 49.27% of gullies are recognized as U-shaped gullies where vertical and sidewall erosion are both operated, having meandering course of main channel (Table 4). In these gullies the secondary duricrust and mottle zone are evenly eroded due to similar resistance power. About 39.55% of gullies are identified as V-shaped gullies   where valley incision is more active than lateral expansion. These active gullies are the initial stage of permanent gully development on the Laterites and these are formed at high angle slope base of hillsides. Only 11.18% of gullies are recognized as trapezoidal gullies where the primary hard laterites are exposed due to deeply incised erosion and valley deepening is stop and valley widening is increased.

Significance of transverse profiles in gullies
Gully catchment 1 has drainage area (DA) of 1,09,250 m 2 , having average drainage density (DD) of 1.8 km km −2 , 28 m of relative relief (RN) and average slope (AS) of 6° towards south-west. High pediment slope with loose ferruginous materials and soils are observed and most of the land is appeared as barren laterite cover with few patches of grasslands. Initially two v-shape gullies are formed, reflecting active phase of erosion and headward migration. At a distance these two gullies are joined to form a single incised channel. The gully becomes wide (dominancy of deposition) with vertical banks due to bank failure and development of bank-side rills and gullies (developing wide U-shape gully floor).
Gully catchment 2 has (BA -1,18,325 m 2 , DD -1.35 km km −2 , RR -19 m and AS -4°) has loose secondary laterite cover at the surface but the hard primary laterite is remained at basal part. Upper catchment is mostly covered under barren laterite and grassland and lower part is under thick vegetation cover due to Acacia plantation site. In the middle part few laterite morum quarries are active. In this land use and land cover many sub-parallel rills are merged to form single gully at upper catchment. The successive profiles reflect that initially two v-shape gullies are formed and at a distance these two gullies are joined to form a single incised channel, having U-shape form. Similarly in this case at initial part vertical erosion is dominated, but with distance the lateral erosion is more dominated.
In the upper part of gully catchment 3 (BA -2,16,050 m 2 , DD -1.9 km km −2 , RR -14 m and AS -3°30ʹ) natural vegetation and plantation patches are observed in the upper catchment, but most of land is covered under grassland and lateritic barren cover. Re-vegetation and excessive deposition are observed in the gully floor (that phenomenon gives an idea about erosion management technique, discussed later). Initially a wide v-shape gully is formed with very steep vertical headcut and with increasing the distance the U-shape valley is formed with increasing wideness. Alongside the banks become very steep and these are eroded by sub-parallel rills and bank gullies.  At the lower section of catchment due to exposure of hard laterite and creamy-white kaolinite layer the vertical incision is stopped, but valley widening is quite active during the rainstorms. It forms a trapezoidal shape of gully floor (Figure 13).

Geomorphic threshold in gully initiation
A geomorphic threshold is a point or period of time that separates different modes of operation within part of a landscape system and the concept of geomorphic threshold is useful identifying those conditions at which a landform is incipiently unstable (Bull, 1980;Schumm, 1980). Based on the data of slopes (S) and drainage areas (A) of 118 gully-head catchments an empirical power regression can be used as geomorphic intrinsic threshold condition for gully initiation on this lateritic terrain. The upstream slopes above gully heads are negatively correlated (r = -0.55) with upstream drainage areas which are used as surrogate for the volume of runoff yield in the study area (Equation (3)). A significant line is fitted through the lower-most scatter points for the study sites which are incised to form gully heads. This empirical straight line (S = 17.419 A −0.2517 , with R 2 of 0.52) represents an approximation to critical slope -area threshold relationship for gully incision (Figure 14). Any site (may be un-trenched or trenched by gullies) lying above this critical line is much prone to gully erosion on this laterite terrain. It is derived that mean critical threshold slope for the initiation of gullies is 2.34° in this region. The high value of a (i.e., 17.419) signifies the initiation of gullies by high volume of overland flow and landslide at micro scale in the study sites (Morgan & Mngomezulu, 2003). Most importantly the constant b is variously interpreted as relative area exponent or relative shear stress indicator (Begin & Schumm, 1984;Morgan & Mngomezulu, 2003). The negative value of b (i.e., -0.2517) and in general consideration b > 0.2 is considered to identify the dominancy of overland flow erosion over sub-surface processes in the study area (Dong et al., 2013;Morgan & Mngomezulu, 2003;Samni et al., 2009;Vandaele et al., 1996;Vandekerckhove et al., 1998). To test b value of regression, the null hypothesis (H O ) is that there is no Figure 13. Final instrument-based head to outlet transverse profiles of gullies, showing downstream changes of profile shapes (narrow V-shape to wide U-shape). significant correlation between the two variables. For 116 degree of freedom (N -2) the tabulated t value is 3.29 in 0.01 significance level (two-tailed) but our calculated t value (7.09) much greater that tabulated t. Thus H O is rejected and alternative hypothesis is accepted, which reflects high importance of b value and it favours a significant inter-relation between S and A in the geomorphic system of gully erosion.

Estimating runoff variability
In gully catchment 1 the values of weighted CN II and CN III are 85.88 and 93.45 respectively (  Figure 15). R c (runoff coefficient) varies from 0.35 to 0.71 in AMC II condition and it increases from 0.62 to 0.86 in AMC III condition respectively. In gully catchment 2 the derived S varies from 42.97 mm (AMC II) to 18.35 mm (AMC III). In the same rainfall events, Q c ranges from 14.62 mm to 97.41 mm in AMC II condition and it varies from 25.94 mm to 118.03 mm in AMC III condition (Table 5.6). Therefore, R c of AMC   In each rainfall event the runoff values of catchment 3 is less than other two catchments due to differences in land use -land cover area. It is understood that for AMC III condition the CN value is increased and S is decreased, because in AMC III condition the soils hold maximum moisture previously due to heavy rainfall of monsoon and the situation promotes more runoff than storage in a rainfall event. If the rainfall amount is increased gradually (during cyclonic condition), the R c is also increased consecutively and it will be high runoff event which is the sign of high vulnerability of flow erosion on the bare slope. The SCS-CN analyses reflects that in prolong rainfall event of tropical depression more than 86% of rainfall can be transformed into direct runoff, as R c is reaching up to 0.86. Alongside, it gives more hydraulic energy to gully initiation and gully head migration in AMC III condition. In short, it can be said the as that CN values of AMC II and AMC III condition are quite high (due to minimum areal coverage of natural vegetation and grassland), the runoff potential and erosion potential of the catchments will remained high in the extreme rainfall event or any torrential rainfall event, if any conservation measure or the transformation of existing land use -land cover practice is not taken.
In a observed record of event based rainfall period (17 th to 19 August 2017) the total rainfall is recorded as 147 mm at the study area and the situation is very resemble with antecedent moisture condition III type (i.e., high moisture condition in the laterite surface). Using SCS-CN model the calculated catchment runoff (sample of 17 gully catchment) varies from 70 to 115 mm in that rainy period. To understand the role of overland flow velocity, an estimation of flow shear stress (T s ) is calculated at 31 sites of gull headcuts within 17 catchments and it varies from 1.78 to 4.56 pascal (Equation (9)), generating mean flow depth of 1.68 cm. The linear regression, between runoff (Q R ) and shear stress (T s ), reveals that with increasing catchment runoff the flow shear stress at gully headcuts is escalating in the study area, having a positive trend relation -Q R = 0.0333 T s + 0.424 (R 2 = 0.342) and product moment correlation coefficient (r) of 0.584 ( Figure  16). In those 31 sites of gully headcuts, the potential transport capacity of overland flow (H) varies from 6.92 to 15.74 kg m −2 (Equation (10)) and it depends on slope steepness of laterite terrain, soil cohesion and ground canopy cover as roughness factor. The linear regression (between H and Q R ) reveals that with increasing runoff (during high moisture condition of monsoon rainfall) the flow transport capacity of eroded material is escalating at the upstream slope of gully headcuts, having a strong positive correlation of 0.837. The empirical relation is expressed as Q R = 2.556 H -12.287 (R 2 = 0.701) (Figure 16).

Channel erosion assessment
In the gully catchment 1, the estimation of land loss due to gully erosion is based on the eight transverse cross-sectional areas at different interval of distance along the gully floor. The average volume of crosssections ranges in between 1504.66 m3 to 9973.37 m 3 (measured bulk density of laterite profile is 2.205 gm cm −3 ), whereas the amount of land loss (Equation (1)) varies from 3317.78 tonne to 21,991.28 tonne (Table 7). It is calculated that all total 79,188.16 tonne of lateritic land is permanently eroded due to gully erosion in this catchment 1. Similarly, in the gully catchment 2 ten consecutive cross-sections reveal that the average amount of land loss varies from 251.64 tonne to 3985.19 tonne. In this catchment all total 21,534.58 tonne of lateritic land is lost due to gully erosion till now. In the gully catchment 3 the analysis of six cross-sections reveals that on an average 6730.43 tonne to 1,61,356.46 tonne of land is eroded along the gully channel and all total 9,48,501.18 tonne of land mass is permanently lost by this gully.
Analyzing the bank morphology of gully (3D geometry tool of Leica Disto S910) it is found that average width of tension cracks in bank sidewall varies from 0.05 to 0.18 m and the gradient of wall varies from 70° 21´ to 85° 15´ in the selected sites. The ultimate results show that the average weight of bank materials (i.e., potential loss of land) varies from 457.23 kg to 1913.72 kg in the bank of gullies. This amount of land mass is vulnerable to next bank failure in the gully catchments. The key processes of bank failure are the weathering and tension crack formation, rill erosion, mass movement and undercutting by ephemeral flow. The GPS-based locational analysis of permanent points (sample of 13 points) reveals that in between 2000 and 2018 (10 years interval) the mean rate of headcut retreat varies from 0.231 m yr −1 to 0.681 m yr −1 in the gully catchments ( Figure 17).
A similar experiment was performed in the 116 gully heads using GPS locations of 2018 and the retrieve data reveals that average rate of gully head migration varies from 0.07 m yr −1 to 1.14 m yr −1 in this laterite terrain.
The role of slope and drainage area on gully initiation is discussed in topographic threshold model (Torri & Poesen, 2014). Also the role of bedrock direction and joint orientation as we all seismic activity and the mechanism by which this happens in controlling and density of gully needs to be better understood (Poesen, 2018). To assess the relation between terrain susceptibility to water erosion (i.e., K factor of RUSLE) and gully density, a linear regression is plotted at catchment scale, taking sample of 17 catchments. It is learn that the erodibility index of laterite terrain varies from 0.23 to 0.61 (mostly sandy clay loam soils) and percentage of coarse materials (i.e., ferruginous coarse sands and Fe-nodules) increase the K value of laterite surface. In the sample catchments the gully density (GD) varies from 1.12 to 2.67 km km −2 , having relative relief of 11-19 m. The established empirical equation is framed as GD = 3.5276 K + 0.2638 (R 2 = 0.557) and the product moment correlation coefficient is 0.745 ( Figure 18). The result reveals that with increasing terrain susceptibility or erodibility the density of gully or gully head development is increasing in the study area.
Although several attempts have been made to develop process-based models for predicting gully erosion rates in a range of environments, there are still no reliable (i.e., validated) models available allowing one to predict gully erosion rates at various temporal and spatial scales (Poesen, 2018). Several regional empirical relations have been proposed to predict gully headcut retreat rates (GHR) and a global analysis of actively retreating gully heads ad their environmental characteristics has been made recently (Poesen, 2011;Vanmaercke et al., 2015). Using the Equation (2) it is found that GHR varies from 48.29 m 3 y −1 to 171.25 m 3 y −1 in the 31 sample point of gully headcuts (Equation (2)). The linear regression reveals that with increasing barren lateritic surface area (B A ) of catchment the GHR is escalating, following the trend as B A = 102.1 GHR + 53.39 (R 2 = 0.544) (Figure 18). There is a strong positive correlation between the two variables, i.e., r value of 0.738. It is proved that barren surface area generates excessive runoff during rainfall and escalates the flow convergence at the vulnerable concave sites of downslope to develop headcut.

Observed sedimentation and erosion rate
It was calculated that in eighteen dam sites the estimated weight of trapped sediments (i.e., mostly ferruginous nodules and coarse sands) varies to a great extent due to activeness of water erosion, slope angle and overland flow length, ranging from 566 kg to 3581 kg. The observed annual erosion rate (O) of  three sample catchments was finally measured as -(a) 10.50 to 24.27 kg m −2 y −1 (gully catchment 1), (b) 8.12 to 20.82 kg m −2 y −1 (gully catchment 2) and (c) 11.87 to 20.82 kg m −2 y −1 (gully catchment 3) respectively (Table 8). The average observed rate is near about 16.27 kg m −2 y −1 which is much greater than the soil Figure 18. Bi-variate regression analysis: (a) with increasing terrain susceptibility to water erosion, i.e., Erodibility Index, the density of gully is rising in the laterite terrain and (b) with increasing barren surface area of catchment the gully headcut retreat is escalating in sample catchments. loss tolerance T-value of this region (i.e., 1.0 kg m −2 y −1 ). Therefore, it can be said that the lateritic badlands of study area have high erosion risk (rendering organic rich top-soil development and increasing Fecrusting, badlands area and degradation of biomass) and the region needs immediate protective measures to check erosion and land degradation at basin scale. Now, it is decided to compare the predicted data of erosion model (RUSLE) with the observed data at field scale.

Predicted erosion rate and comparison
The mean P is estimated as 1510 mm in 2016-17 and I m is calculated as 25.52 mm hr −1 for this climatic region. The analysis has assigned the Rainfall erosivity factor (R) of RUSLE modelling, i.e., 654 for this region. The mean K-factor of laterite terrain is estimated by soil texture and organic matter content of sample soils and the average K values of the catchments varies from 0.23 to 0.28 (Table 9). The length of slope elements or erosion plots varies from 22.1 m to 106.8 m (length in between gully head and water divide), having 55 to 75% of bare lateritic stony surface with development of rills. The steepness of hillslope varies from 3° 45´ to 11° 06´, having average slope of 7° 14´ 30´´ in the sample sites. The C-factor is estimated as weighted value in respect of land use condition in three gully catchments and it varies in each slope elements -(1) 0.61 to 0.91 (gully catchment 1), (2) 0.65 to 0.83 (gully catchment 2) and (3) 0.68 to 0.82 (gully catchment 3). The most important phenomenon is that the study area is not protected under any erosive control measures, except few patches of Accacia plantation. Therefore, in each slope element the P-factor is regarded as 0.1 for RUSLE modelling. Based on the above estimation of inputs, multiplied R, K, LS, C and P factors are taken to get potential or predicated values of annual soil erosion rate (A P ). A P of three gully catchments are estimated as (Equation (11)) -(1) 13.22 to 20.87 kg m −2 y −1 (gully catchment 1), (2) 7.86 to 19.71 kg m −2 y −1 (gully catchment 2) and (3) 16.06 to 24.47 kg m −2 y −1 (gully catchment 3). It is obtained from database that A P of hillslope yields maximum erosion value due to high LS-factor (>1.50). It is found that if the slope is recognized as short length and high steepness, it has high potential for erosion (at dam sites 1, 2, 6, 11 and 14). Based on eighteen dam sites, the average A P is 16.63 kg m −2 y −1 which is beyond the soil tolerance T-value limit (1.0 kg m −2 y −1 ), showing high risk of erosion.
At 0.05 level of significance and 17 (n-1) degree of freedom, the Chi-square (ϗ 2 ) test statistic sets forth the null hypothesis (H o , O -A P or S P = 0) which states that there is no difference between certain characteristics of a population, i.e., difference between predicted and observed value is zero and good correlation. The alternate hypothesis (H 1 , O -A P or S P ≠ 0) reflects significant difference. The value of ϗ 2 statistic is assigned as 27.59 at 0.05 significance level with 17 degree of freedom. The ϗ 2 statistic value of RUSLE is estimated respectively as 10.43 which is much lower than the tabulated ϗ 2 value at 0.05 level. Therefore, it is concluded that H o is accepted and H 1 is rejected. So, there is no significant difference between observed and predicted values in the study. Now, applying model efficiency coefficient (MEC) into the relation between observed and predicted data (Equation (12)), it is found that MEC value of RUSLE is 0.48 (Table 10). The MEC > 0.50-0.70 signifies good and satisfactory performance of model in reference to observed erosion results (Morgan, 2011).
For model evaluation the effectiveness coefficient (E C ) of erosion model is applied to get R diff value (relative difference) which varies from ). It is found from the regression analysis (R diff = m log 10 O E + b) that 55.55% of RUSLE results is placed in overpredicted zone ( Figure 19). The confidence interval of observed erosion rate is calculated as 14.15 to 18.39 kg m −2 y −1 . If the large number of predicted values is fallen within this confidence interval, then E C yields high value (>0.5), signifying the good performance of the model (Nearing, 2000). Calculated E C of RUSLE modelling is 0.61. Therefore, it can be concluded that at 0.05 significance of confidence interval RUSLE model can provide satisfactory results in this region.
It is finally estimated that the predicted values of RUSLE is statistically inter-related with the observed values (A P = 5.90 + 0.659 O E ), having good coefficient of determination (R 2 ) of 0.521 and notable slope (b) value of trend line, i.e., 0.659 ( Figure 12). The t-test statistic of b value is 2.120 at 0.05 significance level with 16 (n -2) degree of freedom (H o : b = 0, Y does not depend on X; H 1 -Y depends on X) and the estimated value of t-test statistic is 2.99 in RUSLE modelling. This analysis reflects that test statistic of RUSLE b-value is greater than the tabulated t-value and it means high dependence of predicted values (H 1 accepted) on the observed values (i.e., RUSLE predicted values resemblance with observed erosion rates). Therefore, based on the seventeen subcatchments of gullies (a part of study area) and RUSLE modelling (considering 118 gully head slope) we have developed an erosion map which depicts the potential annual rate of soil loss due to rill and interrill erosion in the lateritic region. The erosion map ( Figure 20) shows that the western and eastern part is very much susceptible to soil erosion (greater than 9.4 kg m −2 yr −1 ) due to high LS factor and bare soil cover, but the erosion rate (less than 9.4 kg m −2 yr −1 ) is much lower in the central part, because this part is covered with Acacia plantation, Sal forest, aerodrome pavement and relatively low LS factor. The whole region is under very high erosion risk, because the erosion rate is beyond the acceptable T value limit (i.e., 1 Kg m −2 yr −1 ).

Dominant erosion processes
In the profile of secondary laterite there distinct upper zones are identified -(1) dismantled loose laterite and crust, (2) mottled sandy clay zone and (3) kaolinite pallid zone (Figure 21a). At the gully head the formation of grooves and tension cracks in the laterite curst layer develop several tunnels or pipes as seepage lines into the mottled zone, but the due to presence of impermeable kaolinte clay layer the pipes are restricted only in the upper two zones. During rainstorm the tunnel erosion expands the pipes and the overhanging mass of laterite is destabilized over the pallid zone. The roof of gully head has been collapsed and slumped in the gully floor, enhancing the upward migration of head. The gully head and gully wall collapse are a composite and cyclical process resulting from downslope creep, tension crack development, crack saturation by overland flow, head or wall collapse followed by debris erosion which facilities the next failure. Scour of the bed and bank toe increases the bank's height and slope angle, decreasing its stability with respect to mass failure under gravity ( Figure  21(b)).
Three general badland initiation patterns are distinguished in this region. The first two patterns correspond to the expansion of hillslope gullies that can themselves be initiated at mid-slope sections, caused by within-slope conditions or at the slope bottom, through a combination of within-slope and basal conditions. The third pattern corresponds to the disruption of non-channelized hillslope by mass movements that open a bare soil or rock scar to tropical weathering and water erosion. It is observed that mid-slope gullies initiate within hillslopes as "discontinuous gullies" where the local erosive power of runoff overcomes soil-vegetation resistance. Differently, slope-base gullies are initiated by basal scour of the hillslope, through a combination of headwall or sidewall recession and incision of the master network. It is assumed that not only climate fluctuations (strengthening monsoon at late Holocene) or land use changes may play a decisive role, but tectonically induced base-level lowering (movement along Rajmahal basement fault) or local inclination or other discontinuities in drainage network evolution are also decisive factors.
The glimpses of rills and gullies along with bare laterite erodible surface are the sign of active and intensive water erosion. In the study area three broad (1) Above gully head many sub-parallel rills are converged downstream from the water divide. These sites of convergence need immediate vegetative cover.
(2) High slope (> 5°) and long stretch of semiconvex runoff slope (> 70 m) are much prone to dep incision in the torrential rain. (3) Bank failure (> 70°) due to mass wasting, pipe flow, rill erosion and undercutting by gully channel are one of key problems in gully expansion.

Erosion control strategies
In this section the principles and methodology of erosion and sediment control are borrowed from the pioneer works of Gray and Sotir (1996), Morgan (2005), Stokes et al. (2008), Norris et al. (2008), Toy et al. (2013) and Osman (2014). Almost in all areas where inter-rill and rill erosion re concerns include concentrated-flow areas that collect overland flow in a defined channel network and deliver it to a gully and then a watershed outlet. Based on total analysis it can be suggested here that in each gully catchment the focus of erosion control should centred on five aspects -(1) reduction of discharge rate, (2) reduction (1) Reduction of Discharge Rate -A change in land use can change the runoff rate significantly.
Maintaining the good growth of biomass in upper catchment is a key strategy. Water retention basins on the site can be used to decrease the discharge rate.
(2) Reduction of Grade -Channel grade can be reduced by installing a series of drop structure (concrete and rock chutes). Drop structures must be designed and constructed carefully to avoid excessive erosion in the plunge pool.
(3) Control of Headcut Erosion -Gully erosion frequently occurs as headcuts that move upstream in a concentrated-flow area. Drop structures such as those used to reduce the energy grade along a channel can used to protect an overfall. Another way to control headcuts in a gully involves a gully plug, which is basically a dam placed across the gully just below headcut. (4) Downstream Sediment Control -Gravel bags and hay or straw bales are sometimes used as flow barriers. Gravel bags work better than straw bales, and if places sufficiently high and have a wide base, they can be more stable during high rainstorm. In the study area the upper catchment should be protected using gravels bags (placing these bags as upper convex shape along the slope) to decrease the effect of flow convergence (Figure 22). Small sized and convex-shaped check dams (across the gully) just below the gully headcut are found to be an efficient structure to check the flow erosivity and to increase sedimentation during the high rainstorm events. It reduces channel depth and slope. The experimental study reveals that the sixteen dams can trap a significant mass of eroded sediments which varies from 566 kg to 3583 kg respectively. In these trapped sediment the vegetation was growing naturally due to power of natural resilience and these biomass had stopped the active migration of gully heads. (5) Vegetative Measures -Grasses can be used (Lygeum spartum, Brachypodium restusum and Stipa tenacissima) in combination with deeper rooted trees or shrubs (Acacia auriculiformis, Ziziphus mauritiana, Atriplex halimus and Salsola genistoides) along gully walls. On steep slopes, shrubs, e.g., Salsola genistoides would be useful. Brachypodium retusum and reed species, e.g., Juncus acutus could be planted to vegetate drainage lines whereas for stabilizing gully floors a combination of Vetiver grasses and deep rooted tress, e.g., Indian Jujube (Ziziphus mauritiana) should be considered. It is found in the field that if vegetation gets chance to grown in the gully floor, then the whole gully system is started to stabilize upstream ( Figure 23).

Conclusion
In summary, it can be said at last that the present research work has fulfilled the objectives with mentioning the region as a high potential erosion risk at basin scale using measured data and models (16.27 to 18.63 kg m −2 y −1 ). It is found that soil thickness of 0.47 to 1.41 cm y −1 (mean 0.95 cm y −1 ) is permanently lost from the surface of gully catchments. It is estimated that the water erosion will require 127 to 223 years (average 176 years) to erode the mean soil thickness of 1500 mm in this region. It is also found that 52.51% of gullies are affected by overland flow erosion (Slope steepness, S -1.2° to 5.2° and Drainage area, A -2129.1 to 10,513.9 m 2 ) while 27.96% belongs to landslide erosion (S -5.2° to 9.5° and A -457.1 to 5702.5 m 2 ). The topographical character of laterite terrain directs the flow convergence and location of channel initiation and growth of gully network because slope steepness is a prime intrinsic factor in the erosion system, catalyzing the overland flow energy down the slope during thunderstorms. The most challenging task is to grow new plants in the infertile and heavily eroded surface of laterites where progressive expansion of rills and gullies, surface crusting, water crisis in lean period (November to April) and bareness are the key issues. Before applying any vegetative measures we have to understand the root morphology of plants, criteria of re-vegetation, plant selection and design of plantation. To decrease the amount of overland flow and flow convergence no part of upper gully catchment (above gully headcut) should be not be left barren. The convex slope of barren lateritic land signifies active erosion phases, so this land should be protected through grass plantation using flow barriers. As gully control measures are very expensive, therefore, prevention is always better than cure. One gully is never exactly the same as another one, even in the same watershed. So each gully needs separate conservation treatments through identifying the most vulnerable sites of active erosion.
Erosion not only lowers soil quality on-site (hampering pedogenesis processes), but causes also significant sediment-related problems off-site (increasing siltation in rivers and reservoirs). Given the large number of research papers on this topic and present study of gully erosion on laterite terrain, it might conclude that we know now almost everything about the character of gully erosion and its control so that little new knowledge can be added. But still there are many research gaps of gully erosion studies which are key issues to understand the erosion dynamics in Anthropocene and to assess the role of climate changes and active tectonics in the Rarh Plain of Bengal Basin. The present study reveals that there is a need for more research attention: (a) improved understanding of both natural and anthropogenic soil erosion processes and their complex interactions, (b) innovative techniques and strategies to prevent soil erosion or reduce erosion rates with implementation of soil governance and bioeconomy, (c) better understanding of hydro-geomorphic processes and their interactions with external environmental factors at a range of spatial and temporal scales (d) studying the historical gully erosion and palaeogeomorphology (as past is key to present and future), (e) developing a user-friendly empirical or process-based gully erosion models for the predication of erosion rates at a range of environments and (f) assessing the risk of gully erosion in future with implementation of effective conservation practices in the developing countries, like India.