1. Introduction
The Himalaya is a critical global water tower, storing approximately 0.77% of the world’s freshwater in the form of glaciers and snow cover (Azam, Reference Azam2021). These high-altitude mountain glaciers and ice caps, separated from the down valley settlement, provide essential freshwater resources for over a billion downstream inhabitants (Jansson and others, Reference Jansson, Hock and Schneider2003). Several studies, primarily based on remote sensing techniques, have documented the continuous depletion of the Himalayan glaciers due to the ongoing warming, particularly over the past four decades (Bolch, Reference Bolch2012; Brun and others, Reference Brun, Berthier, Wagnon, Kääb and Treichler2017; Maurer and others, Reference Maurer, Schaefer, Rupper and Corley2019; Shean and others, Reference Shean2020; Hugonnet, Reference Hugonnet2021; Garg and others, Reference Garg2022). However, field-based validations of remote sensing estimates are still limited in the Himalayan glaciers, with in-situ glacier observations (mass balance) restricted to only 35 glaciers (Azam and others, Reference Azam2018; Azam, Reference Azam2024). Glaciers across the Himalaya exhibit a consistent negative mass balance trend, although with spatiotemporal heterogeneity induced by regional climate and debris cover. For instance, Dokriani and Chorabari glaciers in the Garhwal Himalaya show an average loss of −0.72 and −0.32 m w.e. a−1, while the Chhota Shigri glacier in the western Himalaya records −0.46 m w.e. a−1. Mera and Pokalde Glaciers in Nepal experience negative mass balances of −0.23 and −0.79 m w.e. a−1, respectively (Mandal, Reference Mandal2020; Vishwakarma, Reference Vishwakarma2022). At regional scale, the average geodetic mass balance across the Himalaya reveals significant ice loss over the past 40 years (1975–2016), with ice loss doubling post-2000 (Maurer and others, Reference Maurer, Schaefer, Rupper and Corley2019). Negative mass balance and continuous glacier retreat generate various surface morphological features such as supraglacial ponds, ice cliffs, and debris-covered areas over the ablation zone of glaciers (Zhang, Reference Zhang2024). Particularly, glacier retreat often leads to the formation of proglacial lakes as meltwater fills depressions left behind by receding ice (Basnett and others, Reference Basnett, Kulkarni and Bolch2013; Patel and others, Reference Patel2017; Sattar, Reference Sattar2023). These proglacial lakes typically begin as small ponds near the glacier’s terminus, dammed by end moraine ridges (Carrivick and Tweed, Reference Carrivick and Tweed2013). Over time, with the glacier’s backwasting and favorable topographic condition, these ponds often expand to form larger water bodies in direct contact with the glacier stream and frontal ice cliffs (Richardson and Reynolds, Reference Richardson and Reynolds2000; Sato and others, Reference Sato, Fujita, Inoue, Sakai and Karma2022). Their formation is generally associated with glaciers lacking efficient drainage pathways for meltwater. The rate of proglacial lake growth is directly linked to the rate of glacier retreat and mass loss, with increasing lake area typically observed alongside a negative mass balance(Ives and others, Reference Ives, Shrestha and Mool2010; Pandit and Ramsankaran, Reference Pandit and Ramsankaran2020; Gantayat and Ramsankaran, Reference Gantayat and Ramsankaran2023). In some cases, the proglacial volume can even reach or surpass the volume of the remaining glacier ice (Pandit and Ramsankaran, Reference Pandit and Ramsankaran2020). Several studies have highlighted the potential for ongoing continuous mass loss of Himalayan glaciers and proglacial lake expansion to create new glacial lake outburst flood (GLOF) hazard zones (King and others, Reference King, Bhattacharya, Bhambri and Bolch2019; Wang, Reference Wang2020; Sattar, Reference Sattar2023).
Glaciers that feed the proglacial lakes undergo a rapid change in surface morphology, such as the formation of a calving front, enhanced retreat rate and higher rate of ice mass loss (Shukla and others, Reference Shukla, Garg and Srivastava2018; Mohanty and Maiti, Reference Mohanty and Maiti2022; Peacey and others, Reference Peacey, Reynolds, Holt and Glasser2023). Notably, calving processes can contribute significantly to annual ice ablation and have often been overlooked as a component of mass balance assessments in the Himalaya (Zhang, Reference Zhang2023). This phenomenon is highly site-specific, and few studies have specifically focused on its impact (Watson and others, Reference Watson2020). A few regional-scale studies have successfully captured the general trend of mass loss associated with lake-terminating glaciers (Zhang and others, Reference Zhang2019; Zhang, Reference Zhang2023). This knowledge can be used to refine future studies and provide more accurate estimates of the current and future state of glaciers, as well as the response of proglacial lakes to climate change (Kaushik and others, Reference Kaushik, Rafiq, Joshi and Singh2020; Watson and others, Reference Watson2020). Additionally, research has explored the role of lake width (lake–ice interface) in influencing mass balance dynamics of lake-terminating glaciers, investigating whether wider lakes promote more negative mass balance (Kirkbride and Warren, Reference Kirkbride and Warren1997).
Motivated by these observations, the present study focuses on the lake-terminating Gepang Gath glacier in the western Himalaya. Our primary aim is to assess the glacier mass balance trends, quantify retreat rates, and understand how proglacial lake expansion drives ice loss at the glacier front. The specific objectives of this study are: (a) to measure the in-situ (glaciological) mass balance of the Gepang Gath glacier during the last decade, (b) to compare the glaciological mass balance estimates with geodetic data for the corresponding period, (c) to quantify long-term changes in proglacial lake between 1962 and 2023 and (d) to assess the impact of glacier–lake interactions on the overall mass balance and state of the glacier. To achieve these, we employed an integrated approach of in-situ measurements with the multisensor and multitemporal remote sensing-based observations acquired between 2014 and 2023.
2. Glacier and Climate Setting
The Gepang Gath glacier is situated within the Chandra Basin of the western Himalaya (Himachal Pradesh, India) (Fig. 1). The total catchment area of the basin is approximately 2446 km2. Notably, the Chandra Basin hosts a significant glacial coverage, with a total of 211 glaciers occupying an area of ∼631 km2 (Oulkar and others, Reference Oulkar, Sharma, Laha, Pratap and Thamban2024). The meltwater generated from these glaciers contributes to the Chandra River. At an elevation of 2850 m a.s.l. (near Tandi), the Chandra River merges with the Bhaga River to form the Chandra-Bhaga River. The region falls in the higher Himalaya under the monsoon-arid zone (Bookhagen and Burbank, Reference Bookhagen and Burbank2006).
Location map of the Gepang Gath catchment in the Chandra Basin, Himachal Pradesh, western Himalaya, India. The location of the main map (i.e. Gepang Gath glacier) is shown in the inset. A dense Differential Global Positioning System (DGPS) survey (yellow line) was done in 2021 and 2022, mainly over the lower ablation area and the peripheral of exiting proglacial lakes. Debris thickness was measured at 50 locations (blue dots), and the point surface mass balances at each site (red bar-ablation; yellow star-accumulation). For the presentation, only the locations of the 2018–2019 stakes are plotted. The background image is Sentinel 2 imagery of October 2018 at 10 m resolution. The glacier outline (green polygon) corresponds to the same image.
2.1. Gepang Gath glacier and proglacial lake
The Gepang Gath glacier is a compound-type valley glacier consisting of two tributary glaciers that merge with a main trunk glacier at an elevation of 4200–4300 m (Fig. 1). The main trunk glacier exhibits debris-free ice down to an elevation of 4200 m, while the debris-covered tributaries extend to a higher elevation of 4400 m. Details of the glacier characteristics are provided in Table (S1). Years of glacier recessions have created a topographic depression, now occupied by a proglacial lake having a total length of ∼2.76 km (as of 2023, Fig. 1). The lake’s spillway, located at an elevation of 4080 m, has maintained its position since at least 1962 (Rathore and others, Reference Rathore2015). Observations suggest that the right side of the terminal moraine is more prone to erosion by overtopping flow events (Worni and others, Reference Worni, Huggel and Stoffel2013; Sattar, Reference Sattar2023). The meltwater from Gepang Gath glacier directly feeds the lake (Fig. 1). Additionally, six small streams, mainly on the northern side of the glacier, contribute to the lake’s water volume (Fig. 1). Currently, the glacier terminus has transformed into a concave-shaped calving front with a total width of ∼900 m (Fig. 1). Because of calving processes, the lower ablation area of ∼1 km length is characterized by numerous transvers crevasses. The warm lake water deep into the ice, which may also afloat the glacier ice is yet to be confirmed.
3. Data and Methodology
Building on the previous studies (Patel and others, Reference Patel2017, Reference Patel2021), we digitized the area of Gepang Gath proglacial lake from 1962 (based on Survey of India toposheets) to 2017 (using Sentinel-2 Multispectral Instrument (MSI) satellite imagery). Subsequently, we obtained proglacial lake area data for 2020 from (Sattar, Reference Sattar2023) and have conducted a ground-based DGPS survey of the lake periphery in 2022. We have been monitoring the glacier-wide mass balance measurement since 2014, and this study presents glacier-wide mass balance data up to 2023 (Table S2). Our analysis focused on changes in the calving front, including cumulative and yearly retreat rates, as well as the change in the frontal glacierized area. The geodetic mass balance and the total volume loss at the calving front were derived by comparing the September 2013 and October 2023 digital elevation models (DEM) at 30 m resolution.
3.1. Kinematic DGPS survey
We made differential kinematic surveys using Trimble R9 carrier phase dual-frequency GPS receivers at a 1 s logging interval over the terminus area, the periphery of proglacial lakes and the lower ablation area of the Gepang Gath glacier during the summer field seasons of 2021 and 2022 (Fig. 1). Three rover receivers were installed on each of the wooden boxes that were carried as a backpack. Data were collected at 1 s intervals, giving the nominal data spacing of about 1 m, the data density was populated while the rover was stationary. The coordinate position was later processed in kinematic mode using a Canadian precise point-processing service (http://webapp.geod.nrcan.gc.ca/geod/tools-outils/ppp.php). The uncertainty in the kinematic survey over the glacier catchment, after post-processing, is ±10 cm in horizontal positioning and ±20 cm in vertical measurement.
3.2. Debris thickness measurement
Debris thickness across the ablation zone of the Gepang Gath glacier was measured manually at 50 points through debris pit measurements. Debris-free ice was assigned to zero debris thickness. Measurements were only conducted at locations where ablation stakes were installed every year. A handheld Garmin GPS unit is used to determine the location with a horizontal accuracy of ±3 m. The error in debris thickness is associated mainly with the height measurement from the ice surface to the debris surface. To minimize the error, we estimated depth on all sides and took the mean value as the final debris thickness. The standard deviation (std) of these values was assumed to be the error in debris thickness.
3.3. Frontal length and area
The annual retreat rate at the front was estimated using Sentinel-2 MSI imageries from 2014 to 2023. The images were obtained from the Copernicus Open Access Hub (https://scihub.copernicus.eu). The images are selected with a minimum cloud cover of <10%, mainly during September and October. The calving front of the glacier was identified and established as the interface between the glacier and the proglacial lake. The frontal outline was initially delineated for the year 2014 and then adjusted according to the temporal satellite images for the years 2015 to 2023.
3.4. Glaciological mass balance
For the mountain glacier, the glaciological mass balance method is well described (Østrem and Brugman, Reference Østrem and Brugman1966; Cogley, Reference Cogley2010). We determined two components for Gepang Gath glacier-wide mass balance (i.e. surface ablation and surface accumulation). It represents the net difference between accumulation and ablation processes across the entire glacier and is expressed in meters of water equivalent (m w.e.). For surface ablation estimation, 6–14 bamboo stakes were drilled into the ice at several points over the glacier surface along with different altitudinal zones for every balance year (i.e. 1 October to 30 September) (Fig. 1 and Table S2). We followed the fixed date system, where 30 September was taken as the end of the mass balance year. Therefore, stake networking was done mostly in September. Surface ablation at a point location was calculated by the sum of the exposed length of stake measured over a balance year. To estimate ice ablation in meter water equivalent from the stake height change, a mean density of ice 870 ± 25 kg m−3 was considered (Pratap, Reference Pratap2019).
Surface accumulation was measured every year at the end of the ablation period (i.e. 15–30 September) via measuring snow density and thickness using a snow corer (Table S2). The thickness of snow accumulation was adjusted with density for total site-specific surface accumulation. The snow/firn density varied from 470 kg m−3 to 570 kg m−3. Due to sparse measurements and steep gradients in the accumulation area, we applied a linear increase between measured points and extrapolated up to 5000 m a.s.l. Beyond 5000 m a.s.l., particularly in steep slopes (covering 19% of the area), a constant extrapolation was applied to estimate total accumulation.
The mean annual altitudinal mass balances were estimated for each 50 m elevation band by averaging all point mass balance data within that band. Finally, the glacier-wide mass balance (Ba) is calculated by integrating these mean annual altitudinal mass balance values over the entire glacier surface area using the following equation:
\begin{equation}B{\text{a}} = \frac{1}{S}\mathop \sum {b_{\text{z}}}{s_z}\end{equation} where Ba is the glacier-wide mass balance (m w.e. a−1) and bz is the mean altitudinal mass balance of the altitudinal area of Sz (m2). 
$S$ is the total glacier area (m2).
The uncertainty in glaciological mass balance depends on uncertainties of in-situ measurements of point surface ablation and accumulation, and the extrapolation of point mass balance to the entire glacier area. Following the refs (Kenzhebaev and others, Reference Kenzhebaev2017; Stumm and others, Reference Stumm, Joshi, Gurung and Silwal2021; Pratap and others, Reference Pratap, Sharma, Patel, Singh, Oulkar and Thamban2023; Azam, Reference Azam2025), we have estimated random and systematic errors for the glacier-wide annual mass balance.
The error in annual surface ablation (Δɑ a) of each point measurement was derived from error in stake height change measurement (Δh), ice density used to convert height change into water equivalent (Δd), and error in debris thickness from the ice surface to debris surface (Δdt) owing to varying debris thickness (Dt). The final error associated with the surface ablation measurements of individual stakes during the mass balance year was calculated as follows (Pratap and others, Reference Pratap, Sharma, Patel, Singh, Oulkar and Thamban2023):
\begin{equation}\frac{\Delta a_a}{a_a}=\text{ }\sqrt{\left(\frac{\Delta_h}H\right)^{2\text{ }}+\left(\frac{\Delta_d}D\right)^{2\text{ }}+\left(\frac{\Delta_{dt}}{Dt}\right)^{2\text{ }}}\end{equation}This results in a mean uncertainty in net surface ablation of ±0.19 m w.e. a−1 for the period from 2015 to 2023.
The error in surface accumulation measurements arises from snow core density variations, measurement errors and surface roughness effects (Thibert and others, Reference Thibert, Blanc, Vincent and Eckert2008; Azam, Reference Azam2025). While point accumulation measurements provide localized estimates, surface roughness and topographic variations introduce additional uncertainties. Furthermore, because total snow accumulation relies on sparse data limited to a specific elevation, this introduces uncertainties in the derived results, particularly the spatial distribution of snow accumulation. The overall uncertainly in net surface accumulation is estimated by
\begin{equation}\Delta ac = \sqrt {{\Delta _{rgh}}{d^2} + \Delta {d^2}{\text{ }}{l^2}} \end{equation} where 
${\Delta _{rgh}}$ is surface roughness 0.35 m w.e. taken from (Thibert and others, Reference Thibert, Blanc, Vincent and Eckert2008), 
$d$ is the mean density of snow/firn accumulation, 
$\Delta d$ is uncertainty in snow/firn densities along the length (l) of the core. The uncertainty in net surface accumulation was on average of ±0.22 m w.e. during 2015–2023.
Considering the uncertainties in net surface accumulation and net surface ablation, we estimated the glacier-wide mass balance uncertainty to range from ±0.20 to ±0.39 m w.e. a−1 with a mean annual uncertainty of ±0.30 m w.e. a−1 during 2014/15–2022/23. Uncertainty in mass balance measurements for individual Himalayan glaciers using the glaciological method varies between ±0.27 and ±0.53 m w.e. (Soheb and others, Reference Soheb2020; Angchuk, Reference Angchuk2021; Wagnon, Reference Wagnon2021). This range is consistent with the uncertainty estimated for the Gepang Gath glacier.
3.5. Geodetic mass balance
The geodetic elevation change of the Gepang Gath glacier was computed by comparing the Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER) DEMs of 17 September 2013 and 30 October 2023. The ASTER DEMs were acquired from NASA’s Earth Observing System Data and Information System (https://www.earthdata.nasa.gov). Before comparison, both the DEMs underwent pre-processing steps to ensure accuracy. First, any erroneous peaks and sinks identified in the DEMs were removed using the thin plate spline method. Then, the 2023 ASTER DEM was taken as the reference DEM and the 2013 ASTER DEM was co-registered to it. The horizontal congruence between both the DEMs was achieved by minimizing the standard deviation of the elevation differences over the stable terrain (Nuth and Kääb, Reference Nuth and Kääb2011; Berthier and others, Reference Berthier, Cabot, Vincent and Six2016). The stable terrain was achieved by masking the glaciers (using RGI7.0 glacier outlines), steep slopes (<4° and >45°) and elevation difference outliers (>±100 m) (Berthier and others, Reference Berthier, Cabot, Vincent and Six2016; Garg and others, Reference Garg, Shukla and Jasrotia2019). Subsequently, the DEMs were checked for vertical biases such as along and cross-track biases originating from satellite attributes. To correct them, the azimuth of ASTER’s ground track was used to rotate the coordinate system and then the elevation differences over stable ground were used to correct the biases (Gardelle and others, Reference Gardelle, Berthier, Arnaud and Kääb2013). Then, considering the fact that the elevation changes are not homogeneous across the glacier surface and varies as a function of altitude (Pieczonka and Bolch, Reference Pieczonka and Bolch2015), a method from (Zhou and others, Reference Zhou, Li, Li, Zhao and Ding2018) was applied to identify outliers at each particular altitude:
\begin{equation}{E_n} = {\text{ }}\frac{{{E_{mx}} - {E_{gl}}}}{{{E_{mx}} - {E_{mn}}}}\end{equation}
\begin{equation}{{\Delta }}{h_{mx}} = A{\text{*}}{E_n}\end{equation}where En is normalized elevation, Emx and Emn are the maximum and minimum glacier elevation, Egl is the elevation of individual glacier pixels, Δhmx is the maximum allowable thickness change at a particular elevation, and A is the empirical coefficient. The empirical coefficient A represents the maximum elevation change at the glacier front which was set to 65 m in order to preserve the true elevation change values while identifying the outliers (Zhou and others, Reference Zhou, Li, Li, Zhao and Ding2018; Garg and others, Reference Garg2022). The eliminated values using Δhmx were replaced by the mean of corresponding 25 m altitude bins.
The volume change (i.e. area-weighted mean) for each 25 m altitude bin was computed by multiplying the area and elevation difference of respective bins. Finally, the geodetic mass balance (B geo) was calculated using the ice density value of 870 ± 25 kg m−3 for the ablation zone and snow density of 570 kg m−3 for the accumulation zone as suggested for the study region (Pratap, Reference Pratap2019).
To quantify geodetic mass balance uncertainties, first, we used the standard deviation of elevation differences over non-glaciated stable terrain (
${\delta _{ng}}$) to approximate the error in elevation change computations (δhng) by dividing it by the square root of a total number of effective pixels (Neff) (Gardelle and others, Reference Gardelle, Berthier, Arnaud and Kääb2013; Zhou and others, Reference Zhou, Li, Li, Zhao and Ding2018). The Neff was obtained using a decorrelation length of 500 m (Brun and others, Reference Brun, Berthier, Wagnon, Kääb and Treichler2017; Shean and others, Reference Shean2020). Finally, the overall uncertainty in mass balance computation 
$\left( {{\delta _{Bgeo{\text{ }}}}} \right)$ was quantified by incorporating the error in density assumption 
$({\rho _{si}})$ as per the following:
\begin{equation}{\delta _{Bgeo}} = {\text{ }}\sqrt {{{({\delta _{hng}})}^{2{\text{ }}}} + {{({\rho _{si}})}^{2{\text{ }}}}} \end{equation}4. Results
4.1. Evolution of the proglacial lake and the surrounding area
The proglacial lake area of Gepang Gath glacier has changed significantly since the 1960s (Fig. 2). As the glacier has retreated over time, the proglacial lake has also expanded dramatically. We combined all published data with our kinematic DGPS survey around the lake perimeter in 2022 (Figs. 1 and 2). Over the past four decades (1979–2017), the lake area has expanded from 0.17 ± 0.06 km2 to 0.84 ± 0.04 km2. Incorporating the 2017 lake area of 0.84 ± 0.04 km2, the DGPS survey conducted in October 2022, and Sentinel-2 imagery from October 2023, the total area is recorded to have expanded to 1.21 ± 0.05 km2. A significant expansion rate in recent years was recorded between 2017-2023 when the area of the lake expanded by ∼0.06 km2 a−1. The trend in lake expansion shows a steady increase in lake area from 1962 to 2000, a moderate rise from 2000 to 2010, and a steep increase from 2010 to 2023, reflecting the accelerated expansion of the proglacial lake in recent times (Fig. 2). The higher expansion rate and steep slope post-2010 is linked to the larger width of the glacier terminus loss and the interaction of lake water with the glacier ice at the bottom.
Time series progression of the proglacial lake area from 1962 to 2023 based on the published and the year 2023 lake area data (area in m2 provided in the literature was converted to km2). Dots show the lake area observed in a specific year. The solid color line shows the positive correlation of lake area expansion between the time intervals of (i) 1962–2000, (ii) 2000–2010, and (iii) 2010–2023. Standard errors (uncertainties) are included for each rate of change in the lake area.
4.2. Glacier surface morphology, frontal length and area change
The combined analysis of the satellite images and debris thickness estimation revealed a complex surface topography within the ablation zone of the Gepang Gath glacier. While the overall surface slope in this zone is notably gentle (less than 20°), there are small-scale undulations with steeper slopes (40°–60°) reaching up to 100–200 m in length. These undulations have an elevation difference of 10–15 m, primarily attributed to the formation of ice cliffs and supraglacial channels (Figs. 3a,b). Several surface features, such as the transverse crevasses, ice cliffs, and debris tills, are also evident over the ablation area (Figs. 3b,d). Two tributaries and a main centreline further characterize the glacier surface topography alongside the debris-covered region, where debris thickness reaches up to 35 cm, a portion of debris-free ice extending along the central line up to an elevation of 4250 m (Fig. 3a). The right and left lateral moraines are elevated by ∼100 m to the present glacier surface. The moraines are highly eroded and fragile. The supraglacial lakes are not prominent across this glacier’s ablation area. Furthermore, there has been no notable increase in supraglacial lakes over time, even in the context of an increase in debris coverage.
Characteristics of Gepang Gath glacier and its proglacial lake. (a) In-situ debris thickness measurement sites with the size of the orange color dots representing the debris thickness over the ablation area. Blue dots represent the debris-free sites based on the stakes installed during 2014–2023. Debris thickness measurements were mostly made in the year 2019, hence plotted over the glacier outline derived in the year 2019. (b) The frontal recession of the Gepang Gath glacier is shown over the Sentinel-2 MSI image of 2023. The frontal length change was observed on a yearly basis from 2014 to 2023. (c–d) field photos of the calving front and proglacial lake and the ablation and accumulation area of the glacier.
The past terminal retreat of this glacier was earlier reported by (Kumar and others, Reference Kumar, Mehta and Shukla2021) for four distinct periods vis-a-vis 1989–2000 (14 m a−1), 2000-2009 (38 m a−1), 2009–2013 (46 m a−1) and 2013–2015 (50 m a−1). We further examined the frontal length and area change from 2014 to 2023. Our study integrated with (Kumar and others, Reference Kumar, Mehta and Shukla2021) reveals a frontal retreat of ∼53 m a−1 (Fig. 3b), indicating an accelerating retreat rate through time. Factors like increasing glacier widths (lake–ice-moraine interface) promote the loss of larger, thick-iced glacierized areas. The transverse profile (width) of the glacier’s front in 2005 was ∼760 m, while in 2014, the width increased to ∼840 m, and in 2023, it reached ∼940 m. This reflects the increased ice mass loss and the expansion of the lake with time (Fig. 3). The yearly length changes of the glacier correlate with the size of the calving ice and the width of the glacier’s frontal portion. Over the period from 2014 to 2023, the glacier has retreated about ∼482 m, and the total glacier area decreased by 0.40 km2, which is accompanied by the expansion of the lake. As the glacier terminates into the lake water, it has a significant impact on the dynamics of the glacier ice, as several transverse crevasses are visible in that terminal area (Figs. 3c,d). This glacier is expected to shrink in the future, and conversely, the lake is expected to expand further, either by calving at the front or by the higher negative mass balance of the glacier.
4.3. Glacier-wide mass balance, gradient, ELA and AAR
The annual glacier-wide surface mass balance (Ba) of the Gepang Gath glacier has predominantly exhibited a negative trend, averaging at −0.90 ± 0.30 m w.e. a−1, except for positive Ba of 0.10 ± 0.28 m w.e. recorded during the mass balance year 2018/19. These values are detailed in Table (1) and shown in Fig. 4. The highest negative Ba was observed during 2017/18, reaching −1.34 ± 0.22 m w.e., closely followed by the Ba − 1.30 ± 0.39 m w.e in 2021/22. The cumulative mass change from 2014/15 to 2022/23 was −7.66 ± 0.84 m w.e. In the lowest ablation area, which is mostly debris-covered and heavily crevassed, the point mass balance during 2014–2023 varies significantly from −5.00 to −1.00 to m w.e. a−1 (Fig. 4). These differential point mass balance rates within similar elevation bands relate to the presence of spatially varying debris thickness. The high interannual variability in ablation rates results from stake replacement under varying debris thickness and differences in seasonal meteorological conditions such as seasonal snow cover and temperature. Also, the spatial distribution of debris cover/thickness over the glacier likely varies over time under the influence of widening crevasses, enlarging ice cliffs, and growing supraglacial ponds. The stakes along the central line (debris-free ice) area showed higher ablation (>−5.0 m w.e. a−1). In comparison, stakes near the lateral moraines with an area under debris thickness of >20 cm had a significant reduction in ablation rate (<−2.0 m w.e. a−1). To show the mass balance gradient, we have presented the annual point ablation and accumulation rates as a function of elevation for every studied year (Fig. 4). The mass balance gradient (i.e. db/dz) shows a linear trend with the mean value of 0.83 m w.e. (100 m−1) (Table 1). The equilibrium line altitude (ELA) of the Gepang Gath glacier ranged from 4435 to 4640 m a.s.l. (average 4583 m a.s.l.). The accumulation area ratio (AAR) ranged from 42% to 62% (average 49%). Furthermore, our analysis revealed a balance budget ELA0 of 4445 m a.s.l. and an associated AAR0 of 62%. These values are consistent with the typical range of AAR0 (55–65%) for glaciers in the Himalaya.
(a) Hypsography, (b–i) point mass balance (ablation: orange dots, accumulation: yellow dots) and the mean annual altitudinal mass balances (black dots) at every 50 m elevation interval of Gepang Gath Glacier for the balance years between 2014/15 and 2022/23, with a gap year of 2019/20. The straight red represents the mass balance gradient (i.e. db/dz) as a linear fit with 95% confidence bond.
Annual net surface ablation (Aa), net surface accumulation (Ca) and glacier-wide mass balances (Ba), ELA and AAR were calculated using the glaciological mass balance method. The ablation and accumulation area, and mass-balance gradients are also presented for the Gepang Gath glacier from 2014/15 to 2022/23. There is no observation for the year 2019/20
4.4. Geodetic mas balance
Geodetic mass balance at 30 m resolution was calculated for the Gepang Gath glacier (Fig. 5a). The surface elevation change analysis from 2013 to 2023 indicates significant thinning of the Gepang Gath Glacier, particularly in the lower ablation zone (Fig. 5a). Notably, the elevation difference on the glacier was highly heterogeneous on the spatial scale. Surface thinning (negative values of surface elevation changes are referred to as surface lowering) is observed on ∼50% of the glacier surface, with 100% thickness loss occurring primarily at elevations below 4150 m. This excessive surface lowering at the terminus reflects the influence of the proglacial lakes which causes the terminus portion to melt or break off through calving (Fig. 5a). The maximum volume loss (area-weighted surface lowering) is observed at 4100–4200 m altitude band while the maximum volume gain is observed at 4700–4775 m a.s.l. Fig. 5b illustrates the variability in geodetic mass balance across different elevations, including the average geodetic mass balance computed for every 25-m elevation. The geodetic mass balance gradient shows significant negative values at lower elevations, indicating substantial mass loss due to melting (Fig. 5b). As elevation increases, the mass loss decreases, reaching values near zero or slightly positive at higher elevations. The total geodetic mass balance for the 10 year monitoring period was calculated as −6.08 ± 1.2 m w.e., which translates to an annual mass balance of −0.61 ± 0.1 m w.e. a−1.
(a) Surface elevation change profile on the Gepang Gath glacier between 2013 and 2023 and (b) average geodetic mass balance of every 25 m elevation band of Gepang Gath glacier.
5. Discussion
5.1 Comparison of in-situ and geodetic mass balances
The geodetic mass balance estimates for the Chandra basin are available in different time frames reflecting the mass loss of −0.68 ± 0.15 m w.e. a−1 from 1999 to 2011 (Gardelle and others, Reference Gardelle, Berthier, Arnaud and Kääb2013) and −0.65 ± 0.04 m w.e. a−1 from 2000 to 2012 (Vijay and Braun, Reference Vijay and Braun2016). While there is no specific geodetic mass balance was observed for the Gepang Gath glacier. Using ASTER data, we have derived surface ice thickness change for the period from September 2013 to October 2023. Further, a spatial surface mass balance was generated by interpolating all yearly point-wise surface mass balance data from 2014 to 2023 to measure spatially variable surface mass balance with respect to surface elevation. The comparison between surface mass balance and ice thickness loss (geodetic mass balance) reveals key differences in ice mass loss (Fig. 6a and b). The surface mass balance ranges from −7.5 to 2 m w.e., showing a gradual transition across the glacier, while the geodetic mass balance varies from −7 to 2 m w.e, exhibiting localized thinning patterns with more spatial variability. The lower glacier tongue shows significant negative values in both cases, but geodetic mass balance captures more heterogeneous thinning, especially around the lower ablation zone with significant surface lowering of up to 90 m (Fig. 6b). The upper glacier regions show lower mass loss in both datasets, yet surface mass balance appears more evenly distributed, whereas geodetic mass balance reflects topographic influences.
Surface mass balance of Gepang Gath glacier: (a) Annual mass balance (points) extrapolated over a digital elevation map to illustrate the spatial surface mass balance pattern from 2014 to 2023. (b) Geodetic mass balance for the period from Sept 2013 to Oct 2023, a zoomed view of the rectangle is shown in (c–e) with three years of cumulative ice volume loss estimated based on the annual retreat rate. (f) Total geodetic ice volume loss from the glacier front during 2014–2023.
As noted, the average annual glacier-wide mass balance (−0.90 ± 0.30 m w.e. a−1) of the Gepang Gath glacier is substantially high compared to the other valley glaciers in the Himalaya (Pratap and others, Reference Pratap, Dobhal, Bhambri, Mehta and Tewari2016; Mandal, Reference Mandal2020; Stumm and others, Reference Stumm, Joshi, Gurung and Silwal2021; Romshoo and others, Reference Romshoo, Abdullah, Murtaza and Bhat2023). We explore this further by focussing on the impact of calving due to the proglacial lake on the annual glacier-wide mass balance. The frontal ice loss by calving, which refers to the breaking-off or detachment of ice chunks from the glacier’s terminus, plays a significant role in determining the overall mass balance of a glacier (Watson and others, Reference Watson2020). When ice calves from the glacier to the lake, it contributes directly to the total ice loss and provides more surface area for the lake’s expansion. Recently, Zhang (Reference Zhang2023) showed that the previous estimate of the total mass loss of lake-terminating glaciers in the Himalayas is underestimated by 6.5 ± 2.1%, emphasizing the need for proper accounting of proglacial lake–ice interaction. Figs. 6c–f shows the trend of ice loss and total volume loss at the frontal part of the glacier between 2014 and 2023. The total volumetric frontal ice loss at the glacier is −21.77 × 106 m3 w.e., which includes a surface area loss of 0.42 km2 during the period from 2014 to 2023. Initially, the frontal ice loss from 2014 to 2017 (3 years) was −6.95 × 106 m w.e. a−1, followed by −6.60 × 106 m w.e. a−1 during 2017-2020. Subsequently, there was an enhanced ice loss rate of −7.78 × 106 m w.e. a−1 during 2020–2023 (Fig. 6e). Frontal ice calving contributed to increased surface ice loss, affecting the estimation of glacier-wide surface mass balance. Overall, the response of the calving front to the annual glacier-wide mass balance of the Gepang Gath Glacier indicates an increase in ice loss by approximately 17–22% during the period from 2014 to 2023.
5.2 Significance of frontal calving on glacier ice loss and interdependency of glacier retreat and lake stability
An increase in proglacial lake volume confirms the retention of glacial meltwater in the glaciated basin. The rise in glacial lake volumes, in turn, significantly contributes to the acceleration of glacier melting through frontal calving and lake-glacier interface processes. Most of these lake-terminating glaciers in the Himalaya have shown an enhanced retreat rate in the last two decades (Patel and others, Reference Patel2017; King and others, Reference King, Dehecq, Quincey and Carrivick2018; Sahu and Gupta, Reference Sahu and Gupta2020). When glaciers terminate into the lake water, they lose glacierized area to the lake, leading to a more negative mass balance. As these glaciers experience accelerated melting and retreat due to climate warming, the ice mass loss from the total glacier volume causes the rapid formation of large proglacial lakes (Zhang and others, Reference Zhang, Xie, Kang, Yi and Ackley2011; Furian and others, Reference Furian, Maussion and Schneider2022) a phenomenon that is vividly reflected in the Gepang Gath glacier catchment. (Prakash and Nagarajan, Reference Prakash and Nagarajan2018) inventoried 46 water bodies in the Chandra basin categorized as ice-dammed, bedrock-dammed, and moraine-dammed lakes. They identified seven lakes to be potentially dangerous sites in the Chandra Basin, with Gepang Gath and Samudra Tapu lakes classified as having high-outburst probability (Prakash and Nagarajan, Reference Prakash and Nagarajan2017, Reference Prakash and Nagarajan2018).
Our analysis of surface thickness change between 2013 and 2023 for the Gepang Gath glacier reveals significant changes, particularly in the lower part of the glacier. The surface changes of up to 90 m in 10 years indicate a notable total ice mass loss in the lower part of the glacier, coinciding with the expansion of a proglacial area. This clearly demonstrates a dynamic shift in the glacier’s morphology and volume change. Currently, the Gepang Gath lake encompasses about 1.21 km2 of area from the origin of the proglacial lake phenomenon, traced back to 1962, when the lake area was approximately 0.2 km2. This study clearly suggests that proglacial lakes significantly impact the overall ice mass loss. The presence of proglacial lakes at the glacier terminus also influences the glacier movement rates and its morphology, thus playing a significant role in overall glacier evolution (Benn, Reference Benn2012; Pronk and others, Reference Pronk, Bolch, King, Wouters and Benn2021). Considering these dynamics and the glaciological aspects of losing ice, it is evident that the glacier ice mass loss over the Gepang Gath glacier will intensify in coming years, and the proglacial lake width may amplify in the coming future until the proglacial lake reaches the full longitudinal extent along the glacier terminus. This, along with the fact that it is categorized as a potentially hazardous lake, necessitates the establishment of an early warning system for Gepang Gath Lake.
6. Conclusions
Our comprehensive study on the lake terminating Gepang Gath glacier, incorporating field-based measurements of mass balance and debris thickness along with DGPS survey and geodetic estimates, provides a better understanding of the glacier response towards climate warming and its interactions with its proglacial lake. The analysis revealed that the proglacial lake associated with the glacier has expanded from 0.20 to 1.21 km2 between 1962 and 2023, showing a nearly sixfold growth during this interval. This expansion is linearly correlated with time, showing accelerated expansion after 2010. The expanded lake area relates to glacier retreat, larger calving ice, and the frontal width of the glacier. The glacier exhibits predominantly negative glacier-wide mass balance values averaging at −0.90 ± 0.30 m w.e. a−1 between 2014/15 and 2022/23. The ice loss at the frontal part due to calving adds approximately 17–22% increase in more negative mass balance. Moreover, the glacial lake expansion coincided with a total glacier retreat of 480 m, area loss of 0.42 km2 and volume loss of −21.77 × 10 6 m3 w.e in the terminus region of the glacier.
Our findings emphasize the significant contribution of the proglacial lake in intensifying the process of ice mass loss from the Gepang Gath glacier. This influence involves both surface mass balance alterations and the frontal ice loss attributed to calving. Notably, a visible disequilibrium between the annual balance and frontal ice loss becomes evident. At that same time, the retreat of the glacier also enabled the systematic expansion of its proglacial lake. All our findings consistently point towards a mutually dependent and compounding effect of the proglacial lake in exacerbating the overall ice mass loss faced by the Gepang Gath glacier. Considering these, it is expected that the Gepang Gath glacier will increasingly lose its mass in the coming years, making its proglacial lake highly vulnerable. Therefore, both the glacier and its proglacial lake need continuous and systematic monitoring at a high spatiotemporal scale and necessitates the establishment of an early warning system for potential cryospheric hazards like GLOF.
Supplementary material
The supplementary material for this article can be found at https://doi.org/10.1017/jog.2025.31.
Data availability statement
Kinematic DGPS surface elevation data and debris thickness data are available at https://data.ncpor.res.in/static/datasets/GNSS_gepang_gath.zip. GIS shapefiles for Glacier and lake outlines, as well as surface ice thickness change raster files, will be available from the author on request.
Acknowledgements
We are thankful to Dr. Lavkush Kumar Patel, Dr. Ajit T. Singh, and our field support team, who have been involved in our field trips and field data collection. P.K. Garg thanks the Director, G.B. Pant National Institute of Himalayan Environment for facilitating the study. We thank the scientific editors, Dr. Fanny Brun and Dr. Nicolas Eckert, and acknowledge helpful suggestions by two anonymous reviewers, which have improved the manuscript. This is NCPOR contribution no. J-2/2025-26.
Author contributions
MT, PS, and BP defined the objectives and designed the study. BP, SNO, and PS collected the in-situ mass balance data. BP led the analysis and interpretation of mass balance and DGPS data, with contributions from PKG, SNO, and MT. PKG analyzed geodetic mass balance data. All authors participated in the preparation of the manuscript.
Funding statement
We thank the Ministry of Earth Sciences for financial support through the project ‘PACER—Cryosphere and Climate.’
Competing interests
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
Open access
