Debris flow susceptibility assessment of Leh Valley, Ladakh, based on concepts of connectivity, propagation and evidence-based probability

The Leh Valley which lies within the Trans Himalayan state of Ladakh, India, is known to be affected almost annually by debris flows ranging from minor to catastrophic scale events. The effect has been getting magnified due to increased urbanization and rapid growth in tourism industry. Though these flows are triggered by intense and abnormal rainfall events the conditioning factor has always been the topography and sediment availability. A lucid acknowledgement of the terrain condition and the degree of vulnerability of such events is required. For this a detail investigation of sediment availability, topographic conditions and their relation with known events becomes crucial. This study utilizes index of connectivity (IC) model to understand the sediment source-sink relationship and farther applied Flow-R model to simulate the probable scenario of events through predefined algorithms. We then use the Weights of evidence (WOE) method to compute the statistical probability of debris flow occurrence. This paper demonstrates the application of these three independent techniques and their implementation in a highly rugged terrain of Ladakh which is a region of frequent debris flows onslaught. The IC and Flow-R models are found to be counter supportive and effective in delineating areas which could be affected by flows that will solely originate in upstream areas where high angle channels directly connected to sediment sources are present. WOE-based model determines the probability of the rare and extensive flows that results from downward integration of other drainage networks in an open fan area.


Introduction
Debris flows occur when a mixture of loose material such as rocks, water, soil and any other sort of debris including plant remains are set in motion under the influence of gravity. They are usually found though not restricted to piedmont zones and steep terrains such as alluvial fans and canyons. Availability of a liquefying agent is an important conditioning factor for a frictionless flow to occur provided debris and loose materials are present in surplus along a sloping terrain (Takahashi 1981). Such flows with bulk densities ranging from 1.5 to 2.2 g/m 3 and velocities ranging from 2 to 12 m/s are extremely buoyant and are capable of loading up to 95% of the entire volume with clasts (Rodine and Johnson 1976;Wang et al. 2018). Debris flow snouts loaded with boulders can exert impact pressures up to 60 kN in an area of just 15 cm 2 (Okuda 1977;Takahashi 1981) and represents a formidable geological agent of denudation and aggradation (Eaton et al. 2003;Eyles and Kocsis 1988). Debris flows rapidly shapes landforms, acts as major channel scouring agent and supplies clasts in fluvial systems (Wohl and Pearthree 1991). The spontanous, erosive nature and association of such flows with drainage networks possess threats especially to mountain dwellers. During 1950During -2011 report suggests 77,779 people died and several settlements and towns were buried by debris flows (Badoux et al. 2016). They were commonly triggered by storm surges, earthquake, lake outbursts, snowmelts, etc. (Wohl and Pearthree 1991). Destruction of infrastructures such as roads, rail lines, water supply networks, telephone lines and agricultural lands introduces secondary long-term socioeconomic impacts (Wohl and Pearthree 1991). Rapid urbanization is a leading cause of occupation of such hazard zones which calls attention for formulating accurate prevention/ mitigation plans.
In the Himalaya, the rugged terrain produces massive debris through gully and glacial erosions. Higher relief above the threshold for orographic precipitation again triggers cloud bursts in steep catchments as a consequence of rapid rise of air masses causing high intensity rainfalls (Dimri et al. 2017). These events leads to high recurrence of landslides and debris flows in the downstream areas due to of the geometry of the catchment (Bookhagen and Burbank 2010;Dobhal et al. 2013). Ladakh is one of the worst affected regions in Himalaya of which one of such cases is the 2010 flood that caused series of debris flows and claimed at least 600 lives and led to destructions of 71 villages, the major hospital, radio station, and public water supply network (Juyal 2010;Ziegler et al. 2016). The disaster was triggered by a short but intense storm that produced 356 mm of rainfall in just 2 h (Gupta et al. 2012;Hobley et al. 2012;Thayyen et al. 2013). Choglamsar was the most affected region of this catastrophic event where some of the survivors are still recovering from post-traumatic stress disorders (Tabassum and Kanth 2012;Ishikawa et al. 2013;Stolle et al. 2015). Similarly, a series of cloudburst and glacial lake outburst triggered disasters have been reported in the years 2006 in Leh and Phyang, 2009in Chushot, 2014at Gya, 2015in Leh, 2017 in Khaltse and 2018 in Saboo and Shey. Though there was no record of casualties, the 2015 flash flood damaged infrastructures and properties worth USD 12 million in Leh (https:// thewi re. in/ envir onment/ ladakh-floods-timel ine-disas ter). It is obvious that Leh and its surroundings are extremely vulnerable and it is increasing as a result of unchecked urbanization (Stolle et al. 2015;Ziegler et al. 2016). With such a growing urbanization and projected warming climate scenario, the higher Himalayan regions especially Ladakh is becoming increasingly vulnerable to such disasters. This necessitates a systematic identification and analysis of the terrain in light of debris flows occurrence probability.
This study, using high resolution digital elevation model (DEM), geomorphic indices and field-based evidences of debris flows attempts to prepare a debris flow susceptibility map of Leh valet, Ladakh. The index of connectivity (IC) is used to assess the degree of linkage between hillslopes and downstream drainage systems. To delineate probable debris source areas we applied Flow-R algorithm and simulated the propagation model of typical alpine debris flows based on spreading algorithms and friction laws. These models were used to compare with field evidences to assess the presence of certain geomorphic features that could cause debris flows (Horton et al. 2013). Then, a susceptibility model based on the evidence of such flows identified in the field is built using GIS-based weights of evidence (WoE) analysis of selective topographic factors antecedent to debris flows. The paper also examines the possibility of assessing debris flows vulnerability of Leh Valley by applying these three independent techniques and assesses their limitation in defining hazard zones. As the Flow-R model has never been used in this area and the equations were tested in the Central Alps, we applied this theory in conjunction with the IC model which has been proven to be effective in the Himalayan Region.

Study area
The study area lies within the Trans Himalayan belt where Ladakh Batholith is the dominant bedrock flanked by Karakorum Range in the North and Zansakar ranges in the south west. The Indus River runs along the contact between the Zanskar range and the Ladakh Batholith termed as the Indus Tsangpo Suture Zone (ITSZ). The Leh valley and altogether Ladakh is composed of complex landforms produced by glacial, fluvial, lacustrine, aeolian and mass wasting process (Juyal et al. 2014). The tributaries of the Indus traverse the Ladakh Batholith forming magnificent valleys and alluvial fans hosting most of the quaternary deposits of the region. They were formed as a result of both mass wasting and sheet flood alluvial processes. Aeolian sand ramps, relict lake sediments and older terraces of Indus and tributary rivers also represents prominent quaternary elements of the landscape (Pant et al. 2005;. The Pleistocene glaciation has a profound influence on the landforms apart from the current surface processes. Extensive U-shaped landforms are the key features of all valleys and remnant glaciers are still present at much higher altitudes. Such glaciers serve as lifelines of inhabitants of major towns such as Leh and surrounding valleys. As a consequence of glacial retreat these valleys are scattered with moraines. The ridges are sharp crested forming cirques. Eventually with the onset of fluvial activity after the Post-Pleistocene deglaciation, the valleys are steadily getting transformed into alluvial plains draining to the Indus. The amphitheater like valleys found in the upstream of these valleys are produced due to high erosion processes (Sant et al. 2011a, b). Several glacial boulders were studied and dated to be produced during recent glacial cycles in the Ladakh Range (Owen et al. 2006). In the last four decades there has been 14% decrease in the glaciated area in Ladakh Range (0.3% −1 year) with maximum shrinkage between 1991 and 2002 (0.6 −1 year) . Glaciers and glacial melts have been acting as an important source of artificial irrigation and agriculture . In Leh and surroundings, areas transitional to glaciers and fluvial systems are connected by fans traversed by perennial and ephemeral streams. This is the location where frequent mass wasting processes occurs. Debris flow is an important agent that occasionally supplies sediment to the alluvial systems. These amphitheater valleys due to their windward orientation and orography experience large number of debris flow disasters (Juyal 2010). The production of massive amount of debris particularly through physical weathering is a leading cause of debris flows during extreme precipitation events.
According to a Land use Land cover analysis through supervised classification of Lnadsat-8 satellite image (taken on August 14, 2015; 7 bands) the study area is divisible into five different units. These are water, settlements with farm lands, snow cover, hillslope and barren area (Fig. 1). Approximately 74.7% of the total area comprises of hillslopes. These are highly rugged granitic terrain (Ladakh Batholith) where mass wasting is a predominant. 16.7% of the area is composed of gently sloping barren lands. These area is least traversed by perennial drainage systems and hence making is dry (Fig. 1). 8% of the total areas is cultivable and utilized for settlement and farming activities. Almost all tributaries that flow toward the Indus River here are used for agricultural purposes. There are no large water bodies apart from the Indus fluvial system comprising of just 0.21% of the study area which is 2943 km 2 (Fig. 1). The snow cover area can change drastically depending on the season of observation. (https:// issuu. com/ tjprc/ docs/ 2ijhr mraug 20192).
This high altitude desert lies 3000 m above the sea level with temperatures dropping below − 20 °C during winters and rises up to ~ 20 °C during summers at present. Leh receives approximately 115 mm of annual precipitation (Banerjee and Dimri 2019). A study of variability of precipitation in Ladakh during 1901-2000 indicates rise in winter precipitation whereas the summer precipitation is decreasing (Shafiq et al. 2016). The main town of Leh has a population of 27,513 and converge main highways such as the Leh-Manali (NH 3), Leh-Srinagar (NH 1D) highways and Leh-Nubra road. Towards the east through the Shakti Village a highway that ultimately leads to Pangong Tso is Located (Fig. 2). The desertic high altitude landscape and unique culture has made Ladakh one of the most sought after high altitude tourist destination leading to rapid infrastructural developments. The rate of expansion of household area is exponential leading to rapid loss of agricultural and barren lands (Dame et al. 2019).

Methods
A 10 × 10 m high-resolution Digital Elevation Model (DEM) of the area is generated from the Cartosat stereo pair satellite images where ground control points (GCPs) were mapped using differential global positioning system (DGPS). The DEM is used to run the Index of Connectivity (IC) and Flow-R. IC and Flow-R are used to delineate areas that are connected as source and sink for debris flows and to simulate the extension of flow events by considering probable source areas. The Flow-R model is particularly applied for detecting triggering areas. Different thematic layers consisting of topographic data diagnostic to debris flows occurrence are also extracted from the DEM and then used for probabilisticbased prediction of debris flow occurrence using weights of evidences (WOEs) method.

Field survey
An extensive field survey was done to detect sediments sources (connected/disconnected), triggering areas, collect GCPs to rectify the DEM and to prepare inventory so as to use it as training data for weight of evidence (WOE) method. The field survey was planned to understand how sediments are stored in the upslope and to see whether they are actively in the process of erosion or dormant in different valleys as well as to understand the overall geomorphology of the fans around Leh and surroundings ( Fig. 3A-D). The relationship of the past debris flows and with the active channels were observed (Fig. 3E). We also observed the relationship of slopes and gullies with main channel to assess the agreement of field data with the constructed models (Fig. 3F, G). Recent and other past flows in the downstream areas were also mapped to use as training data for constructing susceptibility models and for validation (Fig. 3H). The damages caused by recent debris flows were also mapped to compare with the developed topographic models and their derivatives (Fig. 3I). The field survey was done on 2017 and 2019. We mapped several past debris flows that are exposed on the surface irrespective of the time as evidences. For validation data we Fig. 3 Field photographs: A widespread debris sources at Nang fan depocenter, B incised gullies in Nang fan depocenter presently not active, C landscape view of Saboo fan from apex, D landscape view showing sand ramp at Shey, E intersection of past debris flows with stream at Shey, F debris supply gully to main stream at Phyang where renovation of pipeline damaged from gully erosion is in process, G channel intersection with debris sources at Phyang, H hyperconcentrated flows at bassse of Thiksey fan, I damage of culvert due to a recent debris flow at Shaboo Nala randomly selected 10% of the debris flow evidences which was not used in the WOE calculation.

The DEM
The high resolution DEM was created using Cartosat stereo pair satellite images (Path/ Row: 0518/0241 to 0521/0243) and DGPS (Differential Global Positioning System; Model-Leica Viva TS12 Robotic TSV) field survey. We use a stereo derived elevation model of Cartosat-1 where two images taken from different angles are used to calculate the elevation. Cartosat-1 is equipped with double panchromatic camera (2.5 m spatial resolution) one oriented + 25° (fore) and other at −5° (rear) for simultaneous capture of stereo data (https:// www. isro. gov. in/ Space craft/ carto sat-1). The DEM is generated using Leica Photogrammetric Suite (LPS) available in Erdas imagine 14. A block file (Projection: WGS 1984 UTM Zone 43) was created. Stereo data and RPC (Rational Polynomial Function) file was added to it (Biswas et al. 2016). Raw satellite images contain inherent errors and distortions that arise during measurement. For this the RPC file is used to approximate the interior and exterior orientation parameters and then corrects the systematic noise. It is a third order coefficient that gives the mathematical relationship of an image with the sensor geometry (Biswas et al. 2016). Thus, the RPC provides prior constrains for block adjustment of high resolution DEM (Grodecki and Dial 2003). To further augment the orthorectification process, GCPs (ground control points) measured in the field using DGPS, (Differential Global Positioning System) having millimeter scale precision were added to the block file. The x, y and z values of selective GCPs common to both stereo pairs were added using the point measurement tool for each pair while maintaining an even distribution. This is followed by auto tie points generation and finally block triangulation was performed. Before the final extraction, repeated refinement of the process was done to obtain a minimum spatial RMSE (Root mean square error) of 2 pixels. A DEM of 10 m pixel resolution was prepared by mosaic of all the individual DEMs and filling the random voids using spatial analyst tool (ArcGIS 10.1). The vertical accuracy of the DEM was tested by calculating the RMSE (Root Mean Square Error) of various points with respect to the precise DGPS data using the equation- where 'e' is the model error (difference of elevation between the DGPS values and the DEM generated of selected coordinated) and 'n' is the number of observations. We obtain a vertical RMSE of approximately 6 m. The overall vertical accuracy of ASTER (Advanced Spaceborne Thermal Emission and Reflection Radiometer) DEM shows an RMS error of 12.62 m and SRTM (Shuttle Radar Topography Mission) shows an RMS error of 17.76 m in Western Himalaya (Mukherjee et al. 2013). The vertical RMSE also depends on the relative variation of a pixel with respect to the neighboring pixels. As the topography of Ladakh is highly rugged characterized by high density ridges, valleys and steep slopes, a high vertical RMSE is expected. The resolution of a DEM does not always imply the predictive potential of a geomorphic event. Performance of different DEM resolutions (50, 25, 10, 5 and 2 m) in simulating debris flows propagation in Swiss Alps was tested by Horton et al. (2013). The study demonstrated that 10 m resolution as the best indicating good correspondence between source area and flow tracks. Systematic analysis of grid size effect on prediction and simulation of geomorphic and hydrologic process again suggest a 10 m grid resolution gives the perfect balance of smooth data handling and good resolution (Horton et al. 2013). So, the DEM was in optimal configuration for assessing the debris flows susceptibility of Leh and its surroundings. The DEM was used to extract topographic information, channel long profiles and model the debris flows susceptibility.

Index of connectivity (IC)
Index of connectivity (IC) or sediment connectivity is a spatial parameter established to understand the degree of linkage between sediment sources and downstream affected areas (Cavalli et al. 2013). IC depends on coupling between hill slopes where sediments are being produced and valleys where sediments are transported. It is crucial to model the potential efficiency of sediment transfer from source to sink in a catchment scale for hydrological hazard monitoring. We used a method developed by Cavalli et al. (2013) that runs on ArcGIS 10.1 and 10.3 (ESRI 2012) modified after Borselli et al. (2008). This model is specially designed for mountain catchments to tackle debris flows. As study of the IC in the Himalayan river basin demonstrate that the model can be highly beneficial in assessing the hillslope-channel linkages (Mishra et al. 2019). It is based on two aspects: (1) sediment delivery, (2) sediment coupling-decoupling between source and sink. Index of connectivity is defined as where D up and D dn are upslope and downslope components of connectivity. D up is the sediment routing potential of the upslope area expressed as where W is the average weighing factor of the upslope contributing area. Here, we use roughness index (RI) as the weighing factor as suggested by Cavalli et al. (2013), S is the average slope gradient of the upslope contributing area (m/m) and A is the upslope contributing area in m 2 . The D dn quantifies the length of flow path that sediment particle has to travel from the source to the nearest sink defined as where d i is the length of flow path of a particle along the maximum downward slope along the ith cell, W i is the weighing factor and S i is the slope gradient of the ith cell (Borselli et al. 2008;Cavalli et al. 2013).

Flow-R model
Flow-R is a MATLAB-based debris flow susceptibility mapping tool developed at the University of Lausanne by Horton et al. (2013) particularly for regional susceptibility mapping that only requires DEM and its derivatives. The simulation model is exclusively based on the terrain parameters and no climatic variables were considered. Though one can statistically calculate the probability of the occurrence of events or asses the source sink linkage, a realistic simulation is always called for. This is where the application of Flow-R becomes important. It could be considered to be one of the most effective ways to see a realistic simulation and spreading area of debris flows based on exact parameters such as probable source area, travel angle, velocity, and friction. However, it has never been tested in the region. The index of connectivity is comparatively well studied in this context. So, the application of these two models can be complimentary. The Flow-R model has been tested at Sichuan Province (China), Umyeon Mountain (South Korea), Swiss Alps (Switzerland), Frence Alps (France), Scandinavian Mountains (Norway), Italian Central Alps (Italy), Karakoram Highway (Pakistan) and is proved to be effective in modeling debris flows susceptibility (Hürlimann et al. 2003;Horton et al. 2010Horton et al. , 2013Kappes et al. 2011;Fischer et al. 2014;Kang and Lee 2018). It is used to identify (1) potential source areas and (2) their propagation extent. Certain criteria/thresholds were used to identify potential source areas such as slope gradient and upslope contributing area and then choose the algorithms to assess the flow propagation. A threshold of 15° slope is selected (Rickenmann and Zimmermann 1993;Takahashi 1981). Upslope contributing area of a channel should be large enough to supply enough sediments and water for extreme flows to occur. Here, 0.01 km 2 is used as the minimum contributing area which is designed for rare extreme flows (Rickenmann and Zimmermann 1993;Heinimann 1998). The source area identification is purely based on the topographic aspect rather than physical verification. It is assumed that areas with upslope catchment larger than 0.01 km 2 and slopes higher than 15° contain sufficient debris. Based on this criteria all pixels are classified as favorable when initiation is probable, excluded if not or ignored when no decision is feasible. The program contains two options, one for extreme events and one for other rare events. In this case we use extreme events (Rickenmann and Zimmermann 1993;Horton et al. 2013).

Weights of evidence method
Weights of evidence (WoE) method is a Bayesian data analysis approach used for establishing hypothesis using statistical relationship of certain evidences with occurred events (Bonham-Carter 1994). Here, two types of data sets are used: (1) training data and (2) factors. Training data consist of locations of known occurred events. Factors are variables directly or indirectly related to occurrence of that event. In this study training data consist of debris flows inventory prepared using combination of flows mapped in the field and extracted from a previous work (Stolle et al. 2015) (Fig. 1). We use Topographic Wetness Index (TWI), Stream Power Index (SPI), Stream Transport Capacity Index (STI), slope, plan curvature, stream density and aspect as factors (Fig. 4). In the absence of proper evidential data, prior probability (P{D}) is calculated to assess probability of occurring events due to unknown factors (Eq. 5).
where N{D} is the total number of pixels containing the debris flows and N{A} is the total number of pixels of the study area.
When the relationship of causal factors and evidences are quantifiable this can be modified to calculate conditional probability. Conditional probability is more reliable as it is based on the presence or absence of different factors (Fig. 4). The equation used to (5) P{D} = N{D}∕N{A} calculate the presence of debris flow D in the presence of causative factor F in the area A is as follows where Npix is the number of pixels.
Four possible combinations of probability are calculated for each evidential layers using the numbers of pixels: Npix1-when flows occur in the presence of a conditioning factor, Npix2the absence of it, Npix3-absence of flows but factor is present, Npix4-absence of both flows and a particular factor (Fig. 5). The W+ and W− which are positive and negative weights are then calculated and the degree of correlation of a causative factor with the flows is quantified using the formula (Bonham-Carter 1994). where F is the factor and D is the debris flow Positive weight (W+) quantifies the correlation between a causative factor and the flows whereas negative weight (W−) indicates the absence of correlation. The weight contrast factor is then calculated as Following are the parameters that are used in of weights of evidence method: Topographic wetness index (TWI): It is a steady-state index used to assess topographic control of wetness conditions of catchments. It is a function of upslope contributing area and of slope with respect to a point, often used to assess water saturated areas and soil moisture patterns (Chen and Yu 2011;Grabs et al. 2009). TWI is mathematically defined as where A s is the local upslope area draining to a particular point and β is the local slope in radians. Areas with similar TWI will respond similarly during rainfall or any other hydrological event (Qin et al. 2011). Higher TWI is suggestive of higher soil saturation such as in landslide bodies whereas lower TWI suggest high run offs.
Stream power index (SPI) It is a measure of potential erosive power of a flow on a given topographic surface. It is also based on slope and contributing area. where A s is the upslope contributing area and β is the slope of the grid cell.
Sediment transport capacity index (STI) It is a parameter that quantifies maximum amount of sediment that a flow can carry. This is also based on slope and contributing area, defined as where A s is the upslope contributing area and β is the slope of the grid cell (Moore and Burch 1986).
Slope and curvature Slope is the rate of change of elevation from a pixel to another neighboring cell toward the downward direction of maximum rate of change of elevation in the DEM.
where x 2 -x 1 represents the distance and y 2 -y 1 represents the difference in elevation between the two neighboring cells.
Curvature is the second derivative of slope which tells whether a terrain is concave up or concave down. Curvature is generally uses as a factor that controls the formation of gullies and also determines whether a flow will accelerate or decelerate (Kang and Lee 2018). It can also influence runoff erosion (Conoscenti et al. 2013).
Drainage density It is the ratio of the sum of drainage lengths per unit area.
where m is the total length of all streams and channels of an area and m 2 is the area. Drainage density depends on multiple factors such as infiltration capacity, sediment texture, geology, slope, and rainfall. Areas with high drainage density are susceptible to higher water supply and sediment flux. Aspect of slope It is the orientation of slope measured clockwise in degrees from 0 to 360 where 0 is north facing, 90 is east facing, 180 is south facing and 270 is west facing. As rainfall and intensity of sunlight differ in different directions of slope, aspect can be an influencing factor on erosion intensity and regolith thickness.
Based on this method weights are assigned for each class of evidential layer and finally the overall weights of the layers (Table 1).

Results
Seven individual valleys that supports Leh and the surrounding settlements were selected for debris flows susceptibility analysis namely Shakti, Nang, Thiksey, Stakna, Saboo, Leh and Phyang (Figs. 1, 2 and 3). These valleys were abandoned by the glaciers in the last glaciations event. Presently, they are influenced by paraglacial processes. These valleys are triangular in shape, bounded in its sides by Ladakh ranges (Ladakh Batholith) and converses into cirques toward their upstreams (Sant et al. 2011a, b). These valleys host surplus amount of clasts derived from glacial activities and other mass wasting processes. The tributaries of the Indus River which are fed by melts from the small-sized glaciers are the primary sources of water in the valleys (Schmidt et al. 2020). With rapid urbanization several barren as well as former agricultural lands were utilized for housing and settlement activities (Dame et al. 2019). Such zones that are frequently traversed by floods and mass wasting processes imposes menaces to developmental activities. Several debris flows extension and debris sources were identified in the field survey. Most of the debris flows in Phyang were all along the channels ( Fig. 3F and G). These are incised by the channel and tributaries whereas upstream of the valley contains extensive glacial debris that are directly intersected by the main channels (Figs. 2, 3G). The comparatively smaller valley that lies between Phyang and Leh town also contains debris flows and the downstream is heavily occupied by for military purposes. The Leh main valley and Saboo is heavily affected by the recent 2010 flows extending from the upstream to the trunk channel (Fig. 2). Here, debris flows were evidently traversing the settlements and present infrastructures. Though in other portions of the study area several older flow events were identified but their spatial extend were not comprehensible. In the upstream portions, many higher energy flows containing a larger and a higher number of clasts were identified that are spatially restricted to the highly sloping regions of the valleys. (Fig. 3B, E, G). In the downstream portions the flows are extensive but the sediments are much finer indicating reduction in energy and increase in liquidity thus making them more fluid (Fig. 3H). A number of debris sources were identified in upstream portions and sloping gullies (Figs. 3A-C, 12). Majority of the sources are glacial moraines in Phyang, Stakna and Nang whereas several past debris flows that could supply debris in the channels were also identified in Nang and Saboo (Fig. 3E, G, I). The result of the three method used to analyze the intensity of such problems are described in the in the subsequent sections.

Index of connectivity (IC) and Flow-R
The resulting IC map shows the degree of linkage of channelized fluvial system with the hill slopes where maximum sediments are produced as a result of mass wasting and glacial activity. According to the result several source areas that are heavily connected to the main target stream are delineated. An example of Shakti valley is given in Fig. 6. Here, most of the area along the main stream falls within the low connectivity class which may be due to lower relief and presence of flat terrains. The highest connectivity zone is found in the northeastern part of the catchment (Fig. 6). The result of Flow-R shows almost similar result where the north eastern part of the catchment is highly susceptible to debris flow and as such this region will likely produce debris during abnormal precipitation event and reach the target channel. Hence, the results of IC and Flow-R have fair correlation. Source areas that were identified in the field are located mostly in the downstream areas of high connectivity and high debris flow propagation zones in the remaining areas. Likewise, upstream of Phyang valley (Fig. 7A) indicates a similar condition whereas in main Leh area/valley (Fig. 7B) such observation is not deducible as major source materials were not identified in the field. Similarly Saboo, Stagmo and Nang valleys (Fig. 8A, B) also shows similar scenario like that of Phyang where high connectivity and the result of Flow-R simulation is supported by availability of debris. However in Shakti valley again no major debris was identified in the field (Fig. 9B). The overall scenario indicates that the simulation and field evidences are counter supportive. This again suggests connectivity of debris lying upstream into the main channel plays a crucial role in debris reaching the downstream areas and would determine the possibility of finally reaching to the main fan areas. Here, it is important to note that the connectivity does not necessarily mean stream interconnectivity, but rather hill slope and alluvial system connectivity. Thus, IC map shows that Phyang and Leh have higher IC whereas Saboo, Stagmo. Nang and Shakti valleys are low in connectivity (Figs. 7,8,9). The degree of connectivity seems to follow a pattern of systematic reduction towards the east. The probable source area and propagation mechanism/extent computed using Flow-R is based on topographic variables and no field-based source area is usedn. All gullies and channels steeper than 15° and upslope contributing area greater than 1 hectare are considered potential source areas. The result of Flow-R model agrees well with the channel connectivity model. Here according to the models Shakti, Nang and Phyang valleys are potentially vulnerable areas since the ratio of Flow-R-based susceptible zones to catchment area is high (Figs. 6,7,8,9). This would imply higher frequency of debris flow events.

WOEs based susceptibility analysis
The debris flow susceptibility map obtained using Weight of Evidence method is divided into four classes-Very low, Low, Moderate and High (Fig. 10). According to this analysis aspect has the highest weightage amongst the conditioning factors, followed by SPI, TWI, slope, STC and lastly stream density (Table 1). Percentages of area which comes under classes very low, low, moderate and high are 50.25%, 35.10%, 12.80% and 2.00%, respectively (Table 1). For validation, we randomly selected 10% of the debris flow evidences which was not used in the WOE calculation. Within the validation areas maximum pixels falls in the highest susceptibility which is 4 th class (41% area) and the 38% of this area falls within the second most susceptible zone the 3 rd class and the area progressively decreases with decrease in susceptibility. So, 79% of the flows used for validation fall within the two highest susceptible zones of the model. The prediction rate curve indicates the quality of the hazard map and goodness of the weights given to the factors in predicting the validation data (Fig. 11).

Discussion
The fundamental precursor of debris flows that has been occurring in Ladakh, especially in Leh valley is the topographic condition but the triggering factor has always been climate extremes or failure of hydrological structures (Ziegler et al. 2016;Schmidt et al. 2020). Normally intense and persistent rainfall events would favor the rise in sediment pore water pressure that could lower the overall strength of sediment and debris that exist in threshold conditions. This is likely to happen during storm event that produces anomalous amount of water that infiltrate the loosed sediments (during convective storms) and results in slope failure. The changing climatic condition has been known to cause such extreme climatic events and rapid mobilization of sediments downstream. Occupation of such hazard zone has enlarged the risk while the frequency of such events is also known to be increasing (Ziegler et al. 2016). Although the triggering factor is well understood to be mostly metrological it is clear that the long-term factor is the topography and the availability and production of sediments upstream where such threshold conditions of slope failures exists. (Hobley et al. 2012;Thayyen et al. 2012;Ziegler et al. 2016). It is observed from the field evidences that most of the plain areas of Leh Valley and the surrounding fans are stuffed with lobes of past debris flows. Their vertical sections indicate lateral and vertical variation of rheology evident from the sediment architecture which was rapidly released from upslope during these events. At steeper and higher altitude, they are mostly clast supported whereas the downstream parts are composed of finer sediments. The vertical and lateral variability indicate the importance of both magnitudes of the events and topography in determining their rheology. Higher magnitude implies higher energy and deposit larger clasts. These are identifiable in the vertical sections which are manifested as density stratification. High matrix content in all these events are typical examples of water saturation of around 50% volume creating a dense slurry which enables it to carry all the clasts without much traction and clast collisions unlike that of a typical fluvial system. At the apex of a Stagmo fan several incised channels were identified in the field containing sand bodies that indicate the finer fractions ware usually picked up here and the clasts are progressively added downstream transforming into debris flows. The fine and sandy matrix act as both lubricating and transporting medium. These clasts are then deposited downstream where the topography does not allow further progression. Subsequent to deposition it seems to suffer dewatering where only the fine slurry is mobilized further downstream. The aridity of the region does not promote vegetation and consequently enormous debris is produced through physical weathering. During high precipitation events these sediments are evacuated as debris flows and leads to either dissection of the fan or entrenchment. Streams generally perform gradual reworking of these sediments in longer timescale. The Fig. 11 Prediction rate curve used for validation of the weight of evidence-based susceptibility model likelihood of a debris flow event to occur depends on parameters such as slope, the effective catchment area and their shape (Chang and Chao 2006).
The result of this study demonstrates that channel connectivity and propagation model agrees well and the result highlights areas that are well connected can function as potential source areas. The maximum area of Phyang Valley is composed of high connectivity zones (Fig. 7A). Upstream region of the catchment is highly susceptible to debris flows. The filed survey also validates presence of older debris flows materials that can still act as source area. As the Phyang main channel profile indicates high angle and uniform gradient the sediments produced might be routed efficiently (Fig. 12A). In Leh town the upslope areas has high connectivity to the main stream (Fig. 6B). The Flow-R model suggests higher susceptibility in the NE part of the catchment. As the connectivity is high they can reach the main channel that traverses the settlement areas. Leh has comparatively broader catchment and lower channel gradient (Fig. 12B). This might result in production of higher debris with lower routing efficiency as the downstream areas have low gradient (Fig. 12B). This would then result in sediment bulking and it is dangerous as they can be mobilized when the threshold of water infiltration is achieved and would ultimately lead to high magnitude hydrological disaster. Broader catchment area can potentially collect more water in Fig. 12 Channel long profiles of main (longest) streams of all valleys studied the catchment drastically integrating the magnitude of any disaster. It is strongly advisable to avoid any construction activity near the channels or block the natural flow/connectivity of channel networks. Saboo valley on the other hand has slightly lower connected areas however the upslope area is fairly susceptible to debris flow (Fig. 8A). The channel long profile of Saboo Valley posses a typical concave upwards geometry indicating the efficiency of the drainage system (Fig. 12C). Hence, the main channel still possesses a debris flow risk and the fact, that in 2010 Debris flow event in Leh valley, the Saboo Nala was one of most active debris flow stream, provides credence to this inference. Stagmo and Nang have moderate level of connectivity (Figs. 8B, 9A). The upslope of Stagmo has high susceptibility (on basis of Flow-R) (Fig. 8B). The channel profile of Stagmo Valley has convex upward geometry that might have resulted from sediment bulking due to low flow efficiency (Fig. 12D). The same is true for Thiksey fan (Fig. 12E). Poorer routing and sediment overwhelming might be the reason for this geometry which in general is common for dryland rivers. The highly elongated Nang Valley is safe in the downstream whereas upstream regions is not. This is inferred from the Flow-R and IC model that suggest potentials for production of enormous debris and high connectivity to the main channel (Fig. 9A). The typical concave stream profile of Nang fan again suggests the downstream portion is less likely to reach debris flows as the energy might die out (Fig. 12F). The overall Shakti valley has lesser connected areas (Fig. 9B). Although Shakti has overall lower connected region the simulation suggests the likelihood of extreme flow occurrence. This might be a result of presence of highly rugged areas traversed by number of steep sloping streams. However, upslope of Shakti valley is highly susceptible to flows (Fig. 9B). The lower part of the catchment is comparatively safer. Moreover, the channel profile having convex upwards suggest sediment jamming and poorer routing efficiency (Fig. 12G). The neighboring valleys Kharoo and Igoo seem to possess high connectivity with the slopes and Flow-R simulation hazardous zones in comparison to its total area. This observation is reliable as the area was heavily affected during the 2010 flood.
In an overall scenario the result also suggests that susceptibility of debris flows and simulation of Flow-R only considers the upslope steep sopping gullies. However, debris flow events that got transformed into a wider scale in the flat downstream areas cannot be addressed through this model. Such events occur when debris flows get integrated with the main channel flows during anomalously high precipitation events. Here, WoE method is used to tackle this problem that employs past flows as evidences and uses their relationship with controlling factors to generate a probability model. The evidences here are the flows identified in open areas of the main bodies of the fans (Fig. 1). The flows which are restricted to piedmont zones as debris flow fan/cone rarely reach the open fan areas. The only way for debris fans to be transformed into larger fan scale flow is integration with stream networks of the fan body. Occurrence of debris flows however depends on many variables such as presence of debris upslope, propagation routes, climate, topography and rheology of the materials. In this study the evidences were collected from the lower and flatter part of the valleys and fans and not those areas in the steep slopes and valley walls where debris were produced.. So, the WoE-based susceptibility model here only provides susceptibility of deposition and extension of channel modified flows toward the trunk Indus channels and in the open fan area. This will happen when discharge increases leading to integration of debris fans materials produced from slopes. This does not imply areas upstream have low probability of debris flow occurrence. The areas upslope to the highest susceptibility zones can function as source areas and propagation routes that supplies debris to the main channel network of fan.
Altogether this study provides the most likely scenario of debris flows hazards. Susceptibility based on Flow-R in combination with IC provides the probable outcome of debris flows processes when the source area is triggered by intense storm event. These event if simultaneously occurs in a wider scale would lead to convergence of debris from all gullies. An integration effect would be observed in the downstream area. The WOE model will account for these types of flows that affect the flatter downstream areas in the form of hyper concentrated flows. Though the conditioning factors may be determined from terrain analysis the ultimate triggering mechanism will be water availability through atmospheric or snow or permafrost melting. The models presented above can be used as a decisive data for urban planning providing a caution towards debris flows, a common natural disaster in mountains. Improvement in the inventory of debris flows and source areas is suggested as this would increase the prediction accuracy of models. We also encourage identification of sediment availability and rheology analysis complimentary with topography to further assess the threshold precipitation for triggering debris flows for each valley. Debris flow gullies and channels are subject to avulsion which implies occurrence in one gulley does not make it highly susceptible or non-occurrence does not imply low susceptibility. This makes inter gulley comparison for susceptibility on the basis of past occurrence unfeasible. So, we highly advocate avoiding constructions and clearing any civil or private infrastructures within 100 m distance of any ephemeral streams or gullies regardless of their hazard histories. On the other hand since all watersheds are fed by glaciers, frequent monitoring of snow cover area and glacial dam is suggested.
As the result of all three methods portray unique components and projects different scenarios, it is hard to integrate them into a single unified index. Albeit we attempted to give scales for each valleys by taking the percentages of areas falling under the highest two classes of the index of connectivity, percentage of areas falling under the Flow-R susceptibility zones and percentage area falling under the two highest susceptibility class constructed using WoE. An overall rank is assessed by taking the average of these three values (Table 2). According to this method Phyang (44.35) happens to get the highest score followed by Leh (34.91). These could mean both valleys are can produce debris, connected to the main channel and high probability for future occurrences. The remaining valleys have similar scores with Shakti valley having the lowest (23.39). However Shakti valley can produce maximum sediment since 25.22% of the total valley area is under the Flow-R simulation. The likelihood of them transforming into wider scale event seems low due to poorer connectivity as compared to Phyang and Leh.

Conclusion
The Ladakh region produces surplus debris through physical erosion owing to arid climate and absence of vegetation. Inefficient sediment routing due to water limitation leads to sediment overwhelmed hillslopes. Abnormal atmospheric interaction then causes spontaneous rainfall effluxes, which finally results in mass failures such as debris and hyperconcentrated flows. Such flows are known to threaten the survival of the communities residing in Leh and surroundings and jeopardizes civil infrastructure. Our analysis suggests that if the sediment mantled hill slopes in headwater region have high IC with downstream area then there is a large likelihood that even a small debris flow event initiated in headwater region can turn into a disastrous event in the downstream. Therefore a detailed investigation of the topography, in context of debris flows, is needed in town planning. Different models are computed using high resolution DEM and produce probable zones of future occurrence, sediment source and sink linkage and propagation mechanisms. The index of connectivity and Flow-R simulation are counter supportive and adequate in delineating high risk zones in highly rugged and purely alpine settings. However, they seem to possess limitations in the relatively flat and extensive distal parts of alluvial fan system. The susceptibility model constructed based on weights of evidence method fills up this gap. The results of this study can be used as a preliminary data though we recommend further development of the susceptibility model in addition to a detail inventory of past debris flows.