Unreported mass movement cascades and future Glacial Lake Outburst Flood (GLOF) hazard in Jammu and Kashmir, Western Himalaya

The Warwan sub-basin in Jammu and Kashmir is a remote glaciated region containing several glacial lakes and has recently experienced population growth and infrastructure development. Due to inaccessibility and geomorphic masking, multiple mass movement events, including avalanches and glacial lake outburst floods (GLOFs), have remained largely unreported. This study analyzes three major avalanche events from the past two decades. The September 2005 and September 2020 avalanches originated from glaciers GL-B and GL-A within the same glacier complex, while the March 2020 rock–ice avalanche initiated from the headwall of glacier GL-F in an adjacent valley and terminated before reaching its ablation zone. Runout mapping shows that the September 2020 avalanche descended from the headwall of GL-A and impacted its proglacial lake, triggering a GLOF. In contrast, the September 2005 event terminated before reaching the glacial lake that began forming in 1999 at the terminus of GL-B. Geomorphic analysis indicates persistent sediment influx from meltwater streams of GL-D and GL-E into the lake associated with GL-A, progressively infilling the lake. Continued glacier retreat caused meltwater to accumulate behind the sediment infill, dividing the lake into two disconnected lakes. Impact from the September 2020 avalanche led to complete drainage of stored water and sediment, generating a downstream debris flow. Pre- and post-GLOF imagery reveal breaching and widening of the outflow channel and deposition of a debris fan downstream. The repeated history of mass movements and continued growth of glacial lakes raise serious hazard concerns in the Warwan sub-basin. GLOF modeling identifies Lake-B as a potentially hazardous lake in its present condition, posing significant downstream risk. Settlements such as Youdu and Qaderna, along with bridges, roads, and residential structures, lie within potential GLOF inundation zones. High sediment availability from past mobilization and deposition along the valley increases the potential of debris-flow cascades, posing elevated risks to downstream infrastructure and hydropower, underscoring the need for preparedness and mitigation planning. GLOF modeling suggests that an early warning system closer to the source will be able to provide a good lead time in case of potential GLOFs in the valley.

1 Introduction

Since 2000, accelerated glacier ice loss driven by rising temperatures (Azam et al., 2018; Bhattacharya et al., 2021; Bolch et al., 2012; King et al., 2019; Maurer et al., 2019) has been observed in High Mountain Asia (HMA), with lake-terminating glaciers showing substantially greater mass loss compared to land-terminating glaciers (King et al., 2018; Pronk et al., 2021; Sato et al., 2022). This climate change-induced retreat has transformed deglaciated regions into debris- and bedrock-dominated landscapes, fostering glacial lake formation and expansion due to enhanced calving and warm water-ice interactions (Carrivick and Tweed, 2013; Furian et al., 2022; King et al., 2019; Nie et al., 2018; Rai et al., 2024; Sakai et al., 2000; Sattar et al., 2019; Watson et al., 2020). Together, these complex and diverse glacier responses pose a challenge for accurate assessments of future glacier dynamics and associated cascading hazards (Kargel et al., 2011; Pritchard, 2019; Wiltshire, 2014). Glacial lakes have notably increased in number and volume globally in the last three decades (Shugar et al., 2020). With rapidly increasing downstream population, agricultural activity, infrastructure, and hydroelectric development, these expanding lakes pose increasing threats to communities and infrastructure (Haeberli et al., 2017; Sattar et al., 2025b; Schwanghart et al., 2016).

A GLOF event can release millions of cubic metres of water (Sattar et al., 2025b) impacting downstream infrastructure, settlements, and agricultural land, often resulting in significant loss of life and livestock (Sattar et al., 2025a; Zhang et al., 2025a, 2023b; Zheng et al., 2021b). For instance, the October 2023 South Lhonak GLOF in Sikkim, India, devastated the downstream valley (Sattar et al., 2025b; Zhang et al., 2025b). This event also demonstrated the enormous sediment mobilization typical of GLOFs, capable of altering landscapes far downstream. Interestingly, there have been recent GLOF events in the Himalaya that originated from small lakes but had amplified impacts due to erosion dynamics. For instance, the 2016 Gongbatongsha GLOF event in the Poiqu Basin of China (Sattar et al., 2022), which was amplified by the huge debris eroded along the channel, previously deposited during and post-earthquake landsliding of the 2015 Gorkha earthquake.

In the Union Territories of Jammu and Kashmir (J&K) and Ladakh, only a limited number of destructive GLOFs have been documented (Ahmed et al., 2022; Majeed et al., 2021; Schmidt et al., 2020). Historical examples include a 1971 GLOF 18 km upstream of Nymmo village that killed 16 people (Ikeda et al., 2016), and a 2003 GLOF affecting Domkhar village, which destroyed bridges, water mills, and agricultural land (Ikeda et al., 2016). More recently, the 2014 Gya Glacier GLOF damaged agricultural land and other infrastructures in Gya village (Majeed et al., 2021; Schmidt et al., 2020). All these events highlight the dominant role of debris transport, which depends on flow dynamics, channel morphology, erosion, and the availability of loose sediment for transportation. Water released during a GLOF event can efficiently erode and mobilize debris, transforming floods into debris flows or hyper-concentrated flows (Sattar et al., 2025b).

Although numerous GLOFs and mass movements have been documented at both global and regional scales (Emmer and Cochachin, 2013; Lützow et al., 2023a; Nie et al., 2018; Shrestha et al., 2023; Zhang et al., 2024), many events likely remain undocumented owing to the remoteness of affected areas, the absence of significant downstream community impacts, and the limited temporal frequency of satellite observations (Zheng et al., 2021b). Veh et al. (2019) identified 22 previously unreported GLOFs that have occurred since the 1980s. Zheng et al. (2021a) found 176 unreported GLOFs in the HMA, ∼ 95 % were from moraine-dammed glacial lakes, and ∼ 57 % were mapped in the Himalaya. Such unreported events create significant gaps in regional GLOF inventories, leading to underestimated frequency-magnitude relationships (Veh et al., 2022). Identifying these unreported past events is therefore crucial for reconstructing geomorphic processes and understanding the GLOF frequency in mountain catchments (Carrivick and Tweed, 2013; Richardson and Reynolds, 2000; Veh et al., 2025). These undocumented events offer critical insights into GLOF process chains, sediment dynamics, and landscape evolution (Cook et al., 2018; Sattar et al., 2025b), enabling improved hazard modeling and risk mitigation.

While most Himalayan GLOF research focuses on predictive modeling, reconstruction studies remain scarce despite their importance for understanding multi-step process chains capable of triggering GLOFs. Reconstructing GLOF cascades also improves parametrization for hydrodynamic models. For instance, the reconstruction of the South Lhonak GLOF has highlighted how unstable lateral moraines can initiate a sequence of moraine failure, displacement waves generation, overtopping, moraine erosion, and downstream flooding (Sattar et al., 2025b). This study is motivated by the limited research on unreported GLOFs, avalanches, and other mass movements in the Himalayan region. Here, we report and reconstruct mass movements and evaluate future GLOF hazard in a peri-glacial complex. The specific objectives of the study include modeling past mass movements (avalanches) and quantifying future downstream GLOF exposure in the Warwan sub-basin of J&K. Our findings provide a foundation for future GLOF research and offer critical information for downstream infrastructure and land-use planning towards building hazard-resilient mountain communities.

Figure 1 (a) The Warwan sub-basin, in the Chenab Basin of Jammu and Kashmir (Kishtwar district), showing the distribution of the glaciers, glacial lakes, major settlements, previous avalanche source, modeled GLOF path (also see Fig. 11), and individual buildings; (b) glacier complex (GL-A to GL-E) showing the locations of the 2020 avalanche-triggered GLOF and the 2005 avalanche; (c) glacier showing the 2020 avalanche from the headwall of the glacier. (left top subset: HP = Himachal Pradesh; J and K = Jammu and Kashmir; LDK = Ladakh).

2 Study area and importance of the region

The study area is situated in the Fariabaad catchment within the Warwan sub-basin of the Chenab River, located in the Kishtwar district of J&K, India (Fig. 1). The catchment covers an area of 1404 km², with an average elevation of 4308 m above sea level. It contains 272 glaciers with a glaciated area of 315 km² (Pfeffer et al., 2014). The average surface elevation change (2000–2018) across all glaciers in the catchment was −0.33 ± 0.28 m yr⁻¹, with a mean mass balance of −0.28 ± 0.24 m we yr⁻¹ (Shean et al., 2020). This study focuses on a glacier complex comprising 5 glaciers (GL-A to GL-E) (Fig. 1a) and an adjacent glacier (GL-F) located opposite the glacier complex (Fig. 1c). The individual glacier area for GL-A to GL-F are 2.2 ± 0.014, 4.23 ± 0.022, 3.49 ± 0.02  3.47 ± 0.025, 2.16 ± 0.016, 20.5±1.22 km², respectively (Pfeffer et al., 2014). In the glacier complex, Lake-A and Lake-B are proglacial lakes associated with GL-A and GL-B, respectively. The meltwater stream originating from this glacier complex forms the Fariabaad Naala, which merges into the Marusudar River near Yordu village. In the Fariabaad catchment, the number of glacial lakes has increased to twelve, rising from seven lakes in 1990 (Zhang et al., 2023a). Most of the glacial lakes in the catchment are isolated and not connected to any glacier.

The nearby settlements (viz Qaderna, Yurod, Anyar, Rinaie, and Metwan village) are located ∼ 50–55 km downstream from the glacier complex, while major infrastructure such as the Bursar (800 MW) and Pakal Dul (1000 MW) hydropower sites are situated ∼ 90 km downstream. Additionally, several new settlements and agricultural activities have emerged along the riverbank over the past two decades, as identified through high-resolution imagery from Google Earth. Geologically, the region is seismically active with many past earthquakes, clustered in the vicinity of the Kishtwar (Pandey et al., 2017).

3 Data and methods

3.1 Mapping glaciers, glacial lakes, and geomorphic processes

The evolution of glaciers and associated glacial lakes was investigated using cloud- and snow-free Landsat 7 and 8 imagery from 1999 to 2015 and PlanetScope imagery from 2016 to 2024 (Table S1 in the Supplement). These imagery, captured during post-monsoon and peak ablation periods, were used to delineate glacier and glacial lake boundaries using Normalized Differential Water Index (NDWI) and standard FCC (False Colour Composite). Glacier boundaries were further refined through cross-validation with high-resolution imagery from Google Earth and the Randolph Glacier Inventory (RGI) Version 6.0 (Pfeffer et al., 2014). Time-series Landsat and PlanetScope imageries were specifically used to map sediment deposits, crevasses, ice calving, avalanche release areas, and avalanche runouts. Calving processes were identified based on glacier ice loss at the terminus and floating ice on the lake surface. Past GLOF events were identified based on the pre- and post-GLOF moraine changes and the debris deposited downstream of the lake. Associated geomorphic changes from these events were documented using post-event PlanetScope imagery. Mapping uncertainty was calculated using a half-pixel buffer, while lake area uncertainty was estimated using a full-pixel buffer (Granshaw and Fountain, 2006; Racoviteanu et al., 2015).

Figure 2 Schematic summary of modeled scenarios (SCs) for glacial lake outburst flood process chains from Lake-B in the Warwan sub-basin, illustrating mass-movement triggering and propagation, impact-wave formation, lake overtopping, moraine erosion, and downstream flood routing; three scenarios were simulated, involving avalanche volumes of 25 × 10⁶, 12.5 × 10⁶, and 6.25 × 10⁶ m³.

3.2 Avalanche reconstruction

To reconstruct past avalanche events, we used the Rapid Mass Movement Simulation (RAMMS) model, a well-established numerical tool designed for mass movement hazard simulation. RAMMS operates based on the Voellmy-Salm rheological model, a finite volume method that solves two-dimensional depth-averaged mass flow equations for avalanche-like phenomena (Christen et al., 2010). This model accounts for complex interactions between the avalanche material and terrain, making it highly suitable for alpine and glaciated environments (Abhinav and Sattar, 2025; Sattar et al., 2021, 2022). The dynamics of avalanche flow within RAMMS are governed by several key input parameters: the release area, release depth, material density, friction parameters, and the terrain over which it flows. The release areas were identified by detecting changes from pre- and post-event PlanetScope and Landsat imagery. ALOS PALSAR DEM (resampled to 12.5 m spatial resolution) was used to represent terrain for the reconstruction.

To accurately represent avalanche initiation, we conducted simulations with variable release depths/volumes, adjusting these iteratively until modeled runouts matched observed avalanche runouts. Scenario-based modeling initiated with a smaller release volume (R1) and progressively increased volumes (R2–R4) to evaluate the influence of flow magnitude on runout behavior. The material density (ρ) was assumed to be 1000 kg m⁻³, reflecting a mixed composition of snow, debris, and meltwater, consistent with previous avalanche modeling studies (Sattar et al., 2021). The rheological behavior of the avalanche material in RAMMS is controlled by two friction parameters: Coulomb-type basal friction (μ), representing basal sliding resistance, and Turbulent friction (ξ), representing velocity-dependent resistance due to internal turbulence. We set μ=0.12 and ξ=1000 m s⁻², values widely adopted in high-mountain avalanche and debris-flow simulations (Frey et al., 2018; Sattar et al., 2022; Schneider et al., 2014; Somos-Valenzuela et al., 2016). The RAMMS outputs provide spatially distributed information on avalanche runout distances, flow depth, and flow velocity (Fig. 2).

3.3 GLOF process-chain modeling and downstream exposure

To simulate the GLOF process chains from Lake-B we first estimate the future volume of Lake-B by mapping its potential maximum extent. This is achieved by reconstructing the glacier bed using a spatially distributed glacier ice thickness product combined with a DEM (see results Fig. 9a). In the absence of direct ice thickness measurements and subglacial topography data, we utilized ensemble ice thickness data (Farinotti et al., 2019). Raster-based calculations were performed to derive a spatially distributed glacier-bed surface by subtracting the ice thickness from the glacier surface (DEM). The spatially distributed depths of the frontal overdeepening on the glacier bed were mapped following Sattar et al. (2021) (Fig. 9). We followed a contour-based bathymetry generation approach where the depth contours of the frontal overdeepening are extrapolated towards the frontal end of the lake (Fig. 9). This approach was used previously in western Himalayan lake-bathymetry reconstructions (Sattar et al., 2021, 2023). Lake volumes were further compared with volume-area scaling relationships (Cook and Quincey, 2015; Emmer and Vilímek, 2014; Evans, 1986; Fujita et al., 2013; Huggel et al., 2002; Kapitsa et al., 2017; Loriaux and Casassa, 2013; Muñoz et al., 2020; Wang et al., 2012) (Table S2).

To assess potential GLOF impact from Lake-B, we employed a combination of r.avaflow and HEC-RAS hydrodynamic models (v5.0.7) (More et al., 2024; Mergili et al., 2017). Both r.avaflow and HEC-RAS have demonstrated effectiveness in modeling outburst floods in steep Himalayan terrains (Anacona et al., 2015; More et al., 2024; Klimeš et al., 2014; Sattar et al., 2019; Wang et al., 2018). ALOS PALSAR DEM was also used as the terrain for the GLOF simulations.

We defined three scenarios (SCs) for potential GLOF process chain modeling using r.avaflow (see https://www.landslidemodels.org/r.avaflow/, last access: 10 July 2025). Each scenario includes the following processes- (1) avalanche triggering and propagation, (2) avalanche impact on glacial lake and impact wave generation, (3) lake overtopping, (4) moraine erosion, and (5) downstream flood propagation (Fig. 2). The inputs to the GLOF process-chain model including internal friction angle, basal friction angle, and erosion coefficient are adopted from Sattar et al. (2023). The density of the modeled avalanche was assumed to be 2300 kg m⁻³ assuming it to be a mixture of rock and ice, triggering a potential GLOF event from Lake-B. The simulations for processes from process-1 to process-4 above, were performed using r.avaflow, and HEC-RAS was employed to model process-5 (Fig. 2). In scenario 1 (SC-1), an avalanche of 25×10⁶ m³(high magnitude),similar to the Chamoli event (Shugar et al., 2021), was simulated from the headwall toward the lake. Similarly, in scenarios 2 (SC-2) and 3 (SC-3), avalanches with volumes of 12.5×10⁶ m³ (moderate magnitude) and 6.25×10⁶ m³ (low magnitude), were modeled. GLOF outflow hydrographs were extracted immediately downstream of the lake (Section 1; Fig. 11) for all GLOF modeled process-chain scenarios. To represent a conservative approach to map downstream GLOF inundation and flow hydraulics, we consider the high-magnitude (SC-1) as a worst case among all three scenarios considered in the study. The GLOF-process outflow hydrograph of the high-magnitude scenario (SC-1) was used as the upstream boundary condition for downstream GLOF routing using HEC-RAS (see https://www.hec.usace.army.mil/, last access: 10 July 2025). A model domain extending ∼ 85 km downstream from the lake outlet was delineated using a 1 km buffer on either side of the stream. We considered only water routing downstream of the lake. The same approach was followed by Sattar et al. (2025b) for reconstructing the South Lhonak GLOF event in Sikkim Himalaya (Sattar et al., 2025b). A computational mesh with a 12.5 m spacing was generated. Due to a lack of vegetation along the flow channel, a Manning roughness coefficient of 0.04 was used, consistent with previous Himalayan GLOF studies (Rinzin et al., 2026). All other hydraulic and computational parameters were kept default as given in Rinzin et al. (2026). The hydraulic behavior of the GLOF wave in the worst-case scenario (SC-1) was analyzed at different cross-sections (section 1, section 2, and section 3) along the flow channel using modeled output of flood discharge, flow depth, velocity, and inundation extent (Fig. 11). The first cross-section (Section-1) is located immediately downstream of Lake-B. Section 2 is located 49 km downstream of the lake, and section 3 is located immediately upstream of Qaderna village, ∼ 65 km downstream of the lake. Further, the modeled GLOF inundation was overlaid onto high-resolution Maxar imagery (© Google Earth) to manually map exposed downstream infrastructure including bridges, buildings, road segments, and other critical assets.

Figure 3 Evolution of the proglacial lake (Lake-A) and surrounding periglacial environment at GL-A from 1999 to 2024. Continuous sediment infill since 1999 led to bifurcation of the lake (Lake-A) by 2008, with a connecting channel visible in 2009; the main glacier shows consistent retreat throughout the period; the 2020 imagery captures the avalanche event and the fully drained lake following the GLOF. (Background images: Landsat imagery from 1999 to 2015 and © PlanetScope from 2016–2024).

4 Results

4.1 Contrasting evolution of glacial lakes Lake-A and Lake-B

Analysis of multi-temporal satellite imagery indicates a clear trend of glacier retreat within the glacier complex, leading to the formation and progressive expansion of glacial lakes (Figs. 3 and 5). In 1999, both glaciers, GL-A and GL-B, exhibited proglacial lakes, while meltwater from GL-D and GL-E contributed to GL-A's lake system. For GL-B, differential retreat was observed along its margins, with pronounced retreat along the left margin compared to the right (Fig. 4). This asymmetrical retreat appears to be driven by enhanced water-ice interactions, as meltwater inputs from GL-C and GL-A's lake converged into GL-B's lake from the left margin (Fig. 4).

Figure 4 The evolution of the proglacial lake (Lake-B) and associated processes in the periglacial environment of GL-B, from 1999 to 2024; the imagery shows the debris associated with the 2005 rock-ice avalanche, crevasses, glacier ice calving, and glacier retreat; ice floating on the glacial lake (Lake-B) due to calving is also visible in the imagery. (background imagery: Landsat for 1999–2015 and © PlanetScope for 2016–2024).

Over the last two decades, the lake associated with GL-B (Lake-B) has expanded from 0.07 km² in 1999 to 0.30 km² in 2024 (Figs. 4 and Table S3). Lake-A (associated with GL-A) expanded from 0.05 km² in 1999 to 0.06 km², before the 2020 GLOF event (Figs. 3, S1 in the Supplement and Table S3). The smaller area expansion of Lake-A is also due to the partial sediment influx from glacier-fed streams entering the lake on the northern bank. The sediment influx from the meltwater streams of GL-D and GL-E resulted in sediment deposits within the Lake-A basin (Figs. 3 and 5). Temporal area analysis of the sediment deposit revealed a progressive and significant expansion over the study period (1999–2020) (Fig. 5). In 1999, the sediment deposits covered an area of 0.0173 ± 0.008 km², which increased steadily to 0.0602 ± 0.004 km² by 2020 (expansion of ∼ 248 %) (Table S4). The most substantial relative expansion occurred between 1999–2004, with an increase of ∼ 71 % (reaching 0.03 ± 0.012 km² in 2004), likely driven by enhanced sediment supply from GL-D meltwater streams and active sedimentation processes (Fig. 5). Another notable expansion phase occurred between 2004–2008, indicating a ∼ 36 % increase (0.04 ± 0.013 km² in 2008), attributed to the continued glacier retreat and resulting sediment mobilisation. Between 2008 and 2012, the expansion rate of the aerial extent of the sediment deposit reduced to ∼ 10 % (0.45 ± 0.014 km² in 2012), indicating either decreased sediment influx or increasing stabilisation of the depositional environment. Minor yet consistent expansions were observed during 2012–2016 (∼ 17 %) and 2016–2020 (∼ 15 %).

Figure 5 Sedimentation in the Lake-A; spatial extents of sediment deposits from 1999, 2004, 2008, 2012, 2016, and 2020 are overlaid onto 2020 PlanetScope imagery (background imagery: © PlanetScope, 1 September 2020); meltwater streams from GL-D transported sediment into the lake, which initially existed as a single unit until 2008 (Fig. 3), later sediment deposition divided the lake into the fore-basin and main-basin.

This gradually growing sediment deposits led to splitting of Lake-A into two parts: the fore-basin, located between the moraine dam and the deposits (glacier-detached since 2005), and the main basin, located between sediment deposits and GL-A glacier (still proglacial in 2024). While the fore-basin is visible in satellite imagery as early as 1999, the main basin started to evolve from small supraglacial ponds in 2005–2006 (Fig. 3). The fore-basin experienced a gradual decrease in area between 2013 and 2019, driven largely by sediment influx from GL-D and GL-E streams draining a catchment of ∼ 8.0 km², most of which remains glacierized. After the 2020 avalanche and partial breach of the dam, the fore-basin has remained water-free. In contrast, the main basin experienced continued areal expansion driven by glacier retreat both before (2006–2020) and after (2020–2024) the avalanche event. Notably, post-avalanche imagery reveals the main basin was completely filled with ice shortly after the event.

The area of Lake-A can be calculated as a sum of the areas of the main basin, the fore-basin, and the area of the sediment deposits, giving a total area of ∼ 155 000 m². This lake area can be converted to lake basin volume using empirical equations (Sect. 3.3 and Table S5). The calculated mean lake basin depths ranged from 10.7 to 17.4 m and total basin volumes of 1 658 500 and 2 697 000 m³. The estimated volume of only the main basin is between 571 200 and 850 000 m³. These values suggest that the total volume of sediment infill in Lake-A between 1999–2024 is between 808 500 and 2 125 800 m³. Considering that sediment accumulated over 25 years, the mean annual sediment flux is in the order of 0.32–0.85×10⁵ m³ yr⁻¹ (or 0.04–0.11×10⁵ m³ yr⁻¹ km⁻²). If such sediment fluxes continued into the present-sized main basin, the lake may be fully infilled within the next few years to decades. This highlights the short geomorphic lifespan of glacier-detached high-mountain lakes in actively evolving periglacial environments (Carrivick and Heckmann, 2017; Emmer, 2024).

Figure 6 Geomorphic process in the proglacial environment of Lake-A and Lake-B. (a) pre-GLOF scenario where the lake is visible with sediment deposits dividing it into two, and (b) post-GLOF scenario where avalanche deposits are visible, and the lake is completely drained with no water present. Erosion immediately downstream of the lake and deposition of debris is also visible (Background imagery are © PlanetScope).

4.2 Past avalanches in the region

4.2.1 Observations

The September 2020 avalanche, observed in PlanetScope imagery dated 8 September 2020, triggered a GLOF that led to the complete drainage of Lake-A. The rock-ice avalanche initiated from the accumulation zone of GL-A at an altitude of ∼ 4890 m a.s.l. This ice avalanche had a total runout distance of 2355 m, with a mapped inundation area of ∼ 0.78 km² (Figs. 3 and S1). The avalanche impact completely displaced the lake water and sediment deposits within Lake-A, generating a downstream GLOF (Figs. 6 and S1). Following the event, satellite imagery revealed breaching of the Lake-A moraine, extensive GLOF outwash, and subsequent deposition in the downstream channel. However, this GLOF did not lead to a cascading GLOF from Lake B, located immediately downstream. In the months after the GLOF, melting of the ice and debris deposited within the drained lake and the meltwater inputs from the various glaciers that feed the lake, resulted in the reformation of a small proglacial lake in direct contact with the glacier terminus of GL-A (Figs. 3 and S1).

In September 2005, satellite imagery captured an avalanche originating at an altitude of ∼ 5150 m, from the headwall of GL-B, extending ∼ 2620 m (inundation area of ∼ 0.81 km²) into its ablation zone (Fig. 4). However, the runout terminated before impacting the proglacial lake (Fig. 4). During this period, the lake remained relatively small (∼ 0.084 ±  0.057 km² in 2005) but showed gradual lateral expansion due to frontal calving of GL-B. Concurrent with these changes, one tributary of GL-B began detaching from the main trunk. Between 2014–2015, imagery revealed lake water intrusion along the right margin of GL-B's terminus (Fig. 4). Also, extensive crevassing across the glacier terminus culminated in a significant calving event in 2016. The following year, crevasse propagation intensified, leading to another major calving episode in 2017 (Fig. 4, in 2016–2017), which further contributed to lake expansion. The floating ice debris observed in the lake in subsequent imagery confirms that these calving events occurred. Lake-B has expanded ∼ 3.6 times (from 0.09 ± 0.057 km² in 2005 to 0.306 ± 0.008 km² in 2024) since the September 2005 avalanche event. Owing to this growth in lake size, the hazard potential of the lake needs to be assessed for its present condition (see Sect. 4.3). It is also to be noted that another rock-ice avalanche occurred in March 2020 on GL-F (see Fig. 8), that initiated from the headwall at ∼ 5340 m with a runout distance of ∼ 2470 m, affecting ∼ 1.24 km² of the downstream region, in the upper ablation zone. Although this avalanche did not reach a proglacial lake, its scale and geomorphic impacts highlight the ongoing slope-instability processes active across the wider glacier complex.

Figure 7 Pre-GLOF PlanetScope imagery showing (a) the source zone of the 2020 avalanche from GL-A that triggered the GLOF event; (b) the mapped 2020 avalanche inundation; (c–f) reconstructed avalanche runouts (R1–R4) for different release volumes. Scenarios R1 and R2 exhibit the best fit with the actual mapped avalanche inundation; (all background imagery: © PlanetScope).

4.2.2 Avalanche reconstruction

The September 2020 avalanche reconstruction (on GL-A), based on various release depths, led to distinct variations in modeled runout behavior (Fig. 7 and Table S6). In scenario R1, with a release volume of 8.53×10⁵ m³, the modeled avalanche aligned well with the mapped runout. Similarly, in scenario R2, using a slightly higher release volume of 1.13×10⁶ m³, the simulated runout also demonstrated a good match with the observed extent (Table S6). However, as the assumed release volumes were increased by ∼ 50 % (R3) and ∼ 75 % (R4), the modeled runouts exceeded the mapped deposition limits, producing unrealistic flow overshoot and indicating that these larger volumes are inconsistent with actual avalanche runouts.

Figure 8 (a–c) Pre-event Google Earth imagery showing the lateral crack that developed in the upper reaches of GL-F, causing the 2020 avalanche (Imagery © CNES/Airbus, Map data © Google). (d) Pre-event and post-event PlanetScope imagery depicting the avalanche source zone, inundation area, and reconstructed runouts (R1–R4). Scenario R4 showed the best fit with the actual mapped avalanche inundation. (background imagery is © PlanetScope).

For the March 2020 avalanche that occurred on GL-F, reconstruction using varying release depths similarly produced differing runout patterns (Fig. 8 and Table S6). In scenario R1, with a release volume of 8.69×10⁵ m³, the simulated runout did not adequately match the observed avalanche extent. Likewise, scenario R2, with a release volume of 1.3×10⁶ m³, also failed to reproduce the mapped runout. However, when the release volumes were increased to 1.74×10⁶ m³ (R3) and 2.17×10⁶ m³(R4), the R4 scenario provided the best overall fit to the mapped avalanche extent (Fig. 8), suggesting that a larger initial failure volume is required to reproduce the observed deposition pattern for this event.

Figure 9 (a) Spatially distributed glacier bed elevation (in m a.s.l.) of GL-B showing the lake extent in 2000 (blue) and the glacier bed overdeepening extent (green); (b) Spatially distributed depth of the frontal overdeepening of lake GL-B; (c) Depth profile of the overdeepening (AA') along the central line of the lake GL-B and extrapolated bed profile based on the bed slope extending till the frontal moraine of the lake (BB').

4.3 Future avalanche-induced GLOF process chains from Lake-B and downstream exposure

Mapping of Lake-B showed exponential growth (∼ 308 %) between 1999–2024 (see Sect. 4.1; Table S3, and Fig. 4). Also, GL-B, which is associated with this lake, witnessed an avalanche in 2005 (Fig. 4). To address the future GLOF hazard of the lake, the maximum possible future lake extent is considered for GLOF modeling in this study. For this, the modeling of frontal overdeepening of GL-B (Fig. 9a) and frontal overdeepening of Lake-B between 2000 and 2022 shows a maximum depth of ∼ 60 m, with a mean depth of 28.5 m, where a deep glacier depression is observed (Fig. 9b) (Sect. 3.3). Based on this, the maximum lake volume was calculated as 7.7×10⁶ m³. The lake volume estimated using V–A scaling approaches was $\sim \mathrm{7.4}\times {\mathrm{10}}^{\mathrm{6}}$ m³, and comparison between the modeled volume and the ensemble mean of scaling results showed a difference of ∼ 4 %. However, for GLOF modeling, we adopted the ice-thickness-based volume estimate from Farinotti et al. (2019), as it provides robust and widely validated values and has been used in recent Himalayan lake studies (e.g. Sattar et al., 2023). Figure 9c also reveals that the frontal zone of Lake-B is relatively shallow (∼ 13 m depth), whereas the central part reaches depths of ∼ 60 m. The depth variability reflects the gently dipping lake-bed toward the frontal portion. Such glacial-lake beds are typical in nature, as also seen in bathymetric surveys of glacial lakes in the Himalaya, where the lakes are deeper near the glacier terminus and get shallower towards the front (Das and Ramsankaran, 2025; Haritashya et al., 2018; Watson et al., 2020).

In the worst-case (large-magnitude) GLOF process chain triggered by a 25×10⁶ m³ rock–ice avalanche (SC-1) entering Lake-B, peak discharges of 2214 m³ s⁻¹ (liquid phase: lake water) and 1316 m³ s⁻¹ (solid phase: avalanche debris) were generated at section-1, located immediately downstream of the lake (Sect. 3.3; Fig. 11a). 97.5 % of the water present in the lake is potentially drained out of Lake-B. In the moderate-magnitude GLOF process chain (SC-2), initiated by a 12.5×10⁶ m³avalanche, the peak discharge at Section-1 was 1644 (water) m³ s⁻¹, representing a 25 % reduction and draining 91 % of the lake volume. For the low-magnitude scenario (SC-3), the peak discharge was 1260 m³ s⁻¹, a 43 % reduction relative to SC-1, draining 85.71 % of the lake water. The avalanche volume deposited in the lake was 5×10⁶ m³ for SC-1, 3.1×10⁶ m³ for SC-2, and 1.7×10⁶ m³ for SC-3 (Table 1). For all three modeled scenarios, the avalanche mass flow continued beyond the lake and was subsequently deposited downstream. 18.4×10⁶ m³ was deposited downstream of the lake in SC-1, 8.7×10⁶ m³ in SC-2, and 4.2×10⁶ m³ for SC-3 (Fig. 10; Table 1).

GLOF routing in SC-1 showed floodwater propagating downstream and reaching Section-2 after ∼ 5 h and Section-3 after ∼ 6 h after the initiation of the GLOF event. The modeled maximum flow depths were 2.8 and 2 m at Section-2 and Section-3 respectively. The peak discharge attenuated substantially along the flow path due to channel resistance, reducing to 621 m³ s⁻¹ at Section-2 and further to 535 m³ s⁻¹ at Section-3 (Fig. 11b). By the time the flood reaches Section-3 near the populated settlements, the peak discharge is significantly reduced by 73 % compared to Section-1. This emphasizes the critical role of flow attenuation over distance and the influence of valley morphology, channel storage, and flow resistance along the GLOF pathway.

The modeled worst-case GLOF scenario (SC-1) indicates that a substantial number of infrastructures are exposed, particularly those situated along the banks of the Marusudar River near Yordu and Qaderna villages. The analysis reveals that 9 bridges, 40 buildings, and multiple road networks are exposed to GLOF-associated flooding (Fig. 11). GLOF flow depth depths and flow velocity are reaching up to 2–4 m and 4–6 m s⁻¹ at places where infrastructure is located (Fig. 12). GLOF velocities are maximumin the steep channel, immediately downstream of the lake (Fig. 12).

Table 1 Percentage of water released and the total deposition volume of avalanche in the lake for scenarios SC-1, SC-2, and SC-3.

Figure 10 Spatial distribution of avalanche deposition downstream of the Lake-B in scenarios (a) SC-1, (b) SC-2, and (c) SC-3; avalanche deposition volumes are shown in bar charts and are given in Mm³ (×10⁶ m³).

Figure 11 (a) Modeled GLOF process-chain discharge and flow height versus time plot of two phases: avalanche and lake water at Section 1 for SC-1; (b) modeled GLOF inundation from the lake to downstream settlements at Yordu and Qaderna, in the Warwan subbasin; GLOF hydrographs are extracted immediately downstream of the lake (Section 1) for all three modeled process chains (SC-1, SC-2, and SC-3); routed hydrographs for the worst-case GLOF (SC-1) extracted at Sections 2 and 3, at Anyar and before Qaderna; buildings and bridges are overlaid on the flood inundation (SC-1) to identify the exposed infrastructure; subsets B1 and B2 show flow depth and flow velocity overlaid on high-resolution imagery. (Imagery in the subsets © Maxar Technologies/CNES/Airbus, Map data © Google).

Figure 12 The elevation profile of the GLOF flow channel and the modeled average GLOF velocities and depth at 5 km intervals; potential zones of erosion (E) and deposition (D) are marked based on GLOF hydraulics and the slope of the valley; (b–o) high-resolution imagery showing previous mass movement events and their deposits at various locations along the potential GLOF path from Lake-B. (Background imagery in the panels b–o © Maxar Technologies/CNES/Airbus, Map data © Google).

5 Discussion

5.1 Unreported GLOFs in remote high mountains

The HKH is experiencing warming nearly twice the global average, providing favourable conditions for the formation and rapid expansion of glacial lakes (Wang et al., 2015). With rising temperatures, receding glaciers, expanding glacial lakes, and degrading permafrost, GLOFs have become an increasingly frequent phenomenon across the HKH (Pepin et al., 2025). The frequency of GLOFs is expected to rise further as a lagged response to the rapid expansion of glacial lakes, drawing parallels to the lagged response to post-Little Ice Age warming during the 1930s (Harrison et al., 2018). In the remote Himalayan region, GLOFs and mass movement events often remain unreported (Veh et al., 2019, 2022), due to their limited downstream impact and subtle geomorphic signatures. Zheng et al. (2021a) reported that at least 21 GLOF events occurred per year in the region from 1900 to 2016. Unreported GLOFs are events that occurred in the past but were not documented or reported due to their remote location or minimal downstream impact. Also, due to minimal to no impact on downstream population or infrastructure, such GLOF events are often not reported at the time they occur; they are identified retrospectively based on geomorphic processes and signatures. A regional inventory of unreported GLOFs revealed that 57 % (101 events) occurred in the Himalaya, highlighting the region as a global hotspot for unreported GLOFs (Zheng et al., 2021a).

In the case of Lake-A, the unreported 2020 outburst can be attributed to multiple interlinked factors (Fig. 6). First, the outburst released a low water volume, as sediment infilling and reduced meltwater inputs had substantially decreased the lake's storage capacity prior to failure (Emmer, 2024). Consequently, the resulting flood lacked the hydraulic force necessary to induce significant geomorphic disturbances such as extensive channel erosion, debris fans, or bank destabilization that typically serve as visible evidence of past GLOFs (Lützow et al., 2023b). Second, the absence of any notable damage to downstream infrastructure or settlements contributed to its invisibility in hazard records and local knowledge systems (Emmer et al., 2020; Emmer and Cochachin, 2013; Zheng et al., 2021a). It is also highlighted that small-volume outburst floods are often unreported, especially in sparsely populated regions with minimal monitoring networks (Emmer and Vilímek, 2014). Additionally, if the outburst flow is rapidly attenuated within short downstream distances due to channel storage, infiltration, or debris-dominated flow regimes, its cascading impacts may be negligible and undetectable through satellite-based geomorphic change analysis (Allen et al., 2022; Dubey et al., 2024).

The September 2020 GLOF from GL-A, initiated by an ice avalanche originating from its headwall, exemplifies the complex interplay between cryospheric processes and cascading hazards in high mountain regions (Sattar et al., 2022). Although this event breached a small proglacial lake (Lake-A) (Figs. 6 and S3), its downstream impact remained negligible. This reduced storage capacity resulted from significant sediment infilling and meltwater contributions from surrounding glaciers, which progressively lowered the effective reservoir volume prior to the outburst. Similar sediment infilling dynamics have been reported by Emmer (2024), who suggested that such processes can rapidly diminish lake storage capacity and, in some cases, lead to the disappearance of glacial lakes, thereby reducing their hazard potential (Emmer, 2024).

However, the reformation of the lake post-2020 (Figs. 3 and S1), coupled with ongoing sediment influx from adjacent glacier meltwater streams such as GL-E (Figs. 6 and S4), presents a renewed risk for future outbursts (Sattar et al., 2025b). While sedimentation may temporarily reduce hazard potential by lowering lake volume (Emmer and Cochachin, 2013; Worni et al., 2013), it can also create a false sense of stability if dynamic geomorphic and cryospheric processes are not regularly monitored (Emmer et al., 2025; Sattar et al., 2025b). Of particular concern is permafrost degradation, known to destabilize headwalls and slopes, facilitating rock and ice avalanches that can trigger lake breaches (Gruber, 2012; Haeberli et al., 2007; Hjort et al., 2022; Huggel et al., 2008).

The proximity of Lake-A to Lake-B, only a few hundred meters downstream (Figs. 11; S4), further elevates the potential of cascading hazard. Unlike Lake-A, Lake-B shows no visible sediment infill (Figs. 5 and 6) and has expanded significantly (∼ 308 %) over 25 years, retaining a substantial water volume (Table S3). This adjacency raises concern for cascading hazards: a slope failure or avalanche into Lake-A could produce sufficient flood momentum to impact Lake-B, potentially triggering a sequential outburst. Ice-calving and avalanche-induced displacement waves are known triggers of outburst floods (Westoby et al., 2014), while lakes with minimal sediment infill remain at higher hazard potential due to larger water storage (Worni et al., 2013).

Observations confirm that the GLOF and debris flow from Lake-A in 2020 did not affect Lake-B's geomorphology, as no evidence of flooding or erosion was detected downstream of Lake-B (Figs. 6 and S4). Nevertheless, the compounded-risk scenario remains significant where a single mass movement could (i) trigger a GLOF from Lake-A, (ii) transport substantial debris loads that enhance downstream flood intensity, and (iii) destabilize the moraine dam of Lake-B through either hydraulic or mechanical erosion. Such multi-lake cascading hazards can produce impacts far beyond those of single-lake outbursts (Harrison et al., 2018). Lake-B did experience a rock/ice avalanche in 2005; however, the lake was smaller at that time, and the avalanche stopped in the ablation zone before reaching the lake (Fig. 4). In contrast, by 2024, Lake-B has attained its full extent, and if an avalanche of similar magnitude to that of 2005 were to occur under current conditions, it would directly impact the lake, potentially triggering a hazard cascade downstream (Fig. S4).

Time series analyses show that Lake-B receives significant meltwater from GL-C and GL-A (Fig. 6), and PlanetScope imagery indicates large transverse crevasses and floating ice chunks within the lake, signifying active calving (Fig. 4). A large calving event with higher magnitude to those in 2016 and 2017 (Fig. 4) could rapidly displace water, generating displacement waves capable of overtopping the moraine dam and triggering an outburst (Westoby et al., 2014).

The Lake-B poses a potential GLOF risk to downstream villages, including Metwan village, Quadarna, Anyar, and Yordu, located approximately 50–55 km downstream of the lake. Two major hydropower project sites (Bursar; 800 MW and Pakal Dul; 1000 MW) are around 90 km away from the lake. Though potential GLOF flow would reach a steady flow at these locations, transported debris and water volume can adversely impact these sites. High-resolution imagery shows that several settlements developed along the river in the span of two decades (Fig. 11). Additionally, new hydropower projects could lead to further population growth in this valley to meet the region's energy demand, potentially exposing these settlements to considerable GLOF risks in the future. This spatial pattern of exposure underscores the urgent need for community-level early warning systems, strategic GLOF mitigation for highly exposed settlements, and resilient infrastructure design to mitigate the potentially devastating impacts of future GLOF events in the region.

5.2 Permafrost degradation and avalanche formation in the Warwan sub-basin

The progressive increase in regional temperature directly influences permafrost stability in high-altitude glacierized catchments (Gruber et al., 2017), including the glacier complex studied here. This warming leads to the progressive degradation and thawing of ice-rich permafrost interface along steep headwalls and valley flanks (Gruber and Haeberli, 2007; Nicolsky and Romanovsky, 2018). Also, thawing permafrost results in a loss of cohesive strength within rock-ice mixtures, enhancing slope instability, exacerbating surface deformation, and promoting rockfall and rock-ice avalanches (Gruber and Haeberli, 2007).

The avalanche events reported in this study might also be linked to the ongoing permafrost degradation. Figure S2 shows permafrost probability and active layer temperatures (Ran et al., 2022) in the catchment, both of which indicate increasingly unstable ground conditions. This information is particularly relevant for the September 2005 and September 2020 events, which occurred on GL-B and GL-A, respectively. Both glaciers are fed by the steep headwalls underlain by permafrost, where warming-induced reduction in intergranular ice cohesion could have triggered detachment of ice masses and rock blocks (Fischer et al., 2012; Haeberli et al., 2010). Such permafrost-related slope failures can further modify glacier dynamics by depositing debris onto glacier surfaces, increasing surface loading, and modifying ice flow patterns, as observed during the Blatten rock-ice avalanche in the Swiss Alps in 2025 (Islam et al., 2025; Yang et al., 2025a).

We reconstructed the avalanche events, estimating release volumes of 1.13×10⁶ m³ (September 2020) and 2.17×10⁶ m³ (March 2020) (see Sect. 4.2.2). These magnitudes are comparable to GLOF-triggering avalanches in other mountain ranges such as Cordillera Blanca (0.2–0.4×10⁶ m³) (Worni et al., 2014) and the Tibetan plateau (3.8×10⁶ m³) (Yang et al., 2025b). In IHR, studies indicate that avalanches with release volumes exceeding 0.5×10⁶ m³are typically required to trigger a GLOF (Abhinav and Sattar, 2025; Rounce et al., 2017). In this study, avalanche reconstruction relied primarily on mapped spatial extents, due to the absence of field-measured snow/ice thickness, deposit depth, or density data. As a result, the reconstructed release depths and velocities carry uncertainty, and matching modeled versus observed deposit characteristics remains a challenge. Future studies integrating drone surveys, LiDAR, or field-based measurements of avalanche debris depth and flow indicators would greatly improve parameter estimation (Abhinav and Sattar, 2025).

Permafrost degradation is expected to continue under current warming trends, further weakening slope stability and increasing the likelihood of repeated avalanche activity in the catchment. The combination of repeated ice avalanches (volumes $>\mathrm{0.5}\times {\mathrm{10}}^{\mathrm{6}}$ m³), steep permafrost-degrading headwalls, and expanding glacial lakes in the sub-basin poses a serious and evolving GLOF hazard. As Lake-A has re-formed following the 2020 event, the renewed presence of a lake increases the likelihood of future avalanche-triggered GLOFs, raising concern over potential cascading hazards in the region.

5.3 Potential for future debris flow in the valley

The catchment shows clear geomorphic evidence of past mass-movement processes, including debris flows and rock slides, particularly along the river course, which originates from the proglacial lake (Fig. 12) (Rai et al., 2025). Historical slope failures, including large rock slides (Fig. 12d, g, h), appear to have periodically disrupted river flow, suggesting a long history of landscape instability. More recent debris flows, likely triggered by intense monsoonal precipitation (Fig. 12j–n), further demonstrate the catchment's sensitivity to hydrometeorological forcing. These processes have contributed to a significant accumulation of unconsolidated sediments along the valley, increasing the potential for future sediment remobilization. In this context, a future GLOF could act as a high-magnitude trigger (multi-fold hazard cascade), mobilizing existing debris and initiating extensive downstream debris flows (Emmer et al., 2025; Sattar et al., 2023; Yanites et al., 2025). The combination of steep topography, readily available loose material, and hydrological triggers indicates that the catchment is already predisposed to cascading hazards. Modeled GLOFs from Lake-B shows that flows would encounter alternate steep and gentle slopes (Fig. 12a). The major portion of the valley along the GLOF flow path is steep, which will allow high-energy flows and thus more erosion. Thus, any future GLOF event in the valley is unlikely to occur in isolation; instead, it may initiate a cascade of hazards, including debris flows, temporary channel blockages, and subsequent downstream flooding. However, in the current study, we have modeled only the downstream water routing, as modeling the complex interplay between sediment erosion and deposition in the case of a debris flow from the lake is challenging and requires further investigations. This underscores the importance of assessing compound and cascading hazard interactions rather than focusing solely on isolated processes when evaluating risk in glacierized catchments.

5.4 Socio-Economic factors of GLOF exposure and vulnerability

The GLOF inundation map (Fig. 11), prepared based on a potential worst-case GLOF scenario, demonstrates that most infrastructures located along the flow channels of the Fariabaad Naala and Marusudar River fall within the GLOF inundation zone. Notably, road networks and residential buildings along the riverbanks of Fariabaad Naala, at Anyar and Rinaie villages, fall within these GLOF zones due to their immediate proximity to the river course.

The exposure and vulnerability of communities to GLOFs are deeply rooted in the socio-economic characteristics of affected regions. Numerous studies (Carrivick and Tweed, 2016; Khanal et al., 2015) highlight that settlements in high-mountain environments are inherently exposed due to their proximity to glacial lakes and dependence on valley-bottom resources (Muhammad et al., 2021). The concentration of infrastructure, such as hydropower facilities, bridges, and road networks, along river corridors significantly amplifies economic exposure, as seen in previous GLOF disasters in Bhutan (Komori et al., 2012), Nepal (Cook et al., 2018; Sattar et al., 2025a), and more recently in Uttarakhand (Shugar et al., 2021), Karakoram (Muhammad et al., 2021), and Sikkim Himalaya (Sattar et al., 2025b).

Vulnerability, however, is shaped not only by physical exposure but also by socio-economic capacities such as income level, education, access to resources, social networks, and institutional preparedness (Huggel et al., 2015). Marginalised and economically weaker communities are disproportionately vulnerable, as their limited adaptive capacity restricts proactive risk reduction measures (Carey et al., 2012). For example, dependence of the local economy on agriculture along floodplains and limited access to early warning systems increase both direct and indirect risks (Shrestha et al., 2010). In the Himalayan region, increasing agricultural expansion and newly built infrastructure along the floodplain exacerbate exposure, placing more households and assets directly in harm's way.

For hazard mapping at villages such as Qaderna, Yurod, Anyar, and Rinaie, which are located along the shoreline of the flow channel, the baseflow contribution was assumed negligible relative to the hydraulic force generated by the outburst flood. However, it is important to acknowledge that in a modeled high-magnitude scenario, the addition of baseflow to the GLOF peak discharge would cumulatively increase the total flow volume, potentially resulting in elevated water levels and enhanced flood stage (Carrivick and Tweed, 2016; Veh et al., 2020). This could marginally expand the inundation extent, particularly in confined channel reaches or low-lying areas adjacent to the river, thereby intensifying the hazard levels faced by these settlements and increasing the vulnerability of infrastructure situated close to the riverbanks.

GLOF modeling of Lake-B shows the arrival time at Anyar and Qaderna is 5 and 6 h, respectively, in SC-1 (worst-case scenario). We understand that the GLOF hydraulics and arrival times are sensitive to the varied process chains and the modeling approaches, which lead to different arrival times apart from those modeled in this study. However, an early warning detection system that can include discharge sensors, and/or water level sensors, and/or seismic stations, closer to the source (Lake-B) will potentially be able to provide a good lead time in case the of GLOF warnings. However, effective GLOF adaptation requires integrated approaches that combine glacial lake hazard assessments with community-based vulnerability evaluations and social protection policies. Such approaches can reduce exposure while enhancing resilience in vulnerable high-mountain communities. All downstream infrastructure, including settlements, roads, bridges, and hydropower projects, is located along the potential GLOF path from Lake-B.

6 Conclusions

This study provides a comprehensive understanding of the evolving cryospheric and geomorphic processes in the Warwan sub-basin over the past 25 years. Multi-temporal satellite analyses revealed significant glacier retreat, particularly for GL-A and GL-B, contributing to the expansion of associated proglacial lakes. Lake-B exhibited a remarkable growth of ∼ 308 % between 1999 and 2024, primarily due to glacier retreat, meltwater inflow, and successive calving events. In contrast, Lake-A showed a modest expansion (∼ 26 %) until it was completely drained in the 2020 GLOF event, triggered by an ice avalanche entering the lake. Permafrost distribution in the region suggests that permafrost warming can play a significant role in triggering mass-movement processes and promoting the destabilization and degradation of steep slopes. Future GLOF modeling suggests Lake-B poses a substantial downstream GLOF risk in the valley. Despite attenuation along the flow path, vulnerable downstream settlements, including Yordu and Qaderna, face significant flood exposure. GLOF exposure assessment reveals that numerous bridges, residential structures, and road segments lie within modeled inundation zones. The possibility of GLOF to debris flow cascades is very high in the valley owing to the available sediment for remobilization, which in turn can pose a higher risk to downstream infrastructure and hydropower.

The findings emphasize the dynamic interplay between glacier retreat, sedimentation, lake expansion, and avalanche activity in shaping GLOF hazards in the Warwan sub-basin. The results underscore the critical need for integrated monitoring of cryospheric processes, improved early warning systems, community-based preparedness, and robust mitigation strategies to safeguard life, livelihoods, and infrastructure in the region.