Multi-Temporal Assessment of Glacier Change in the Arun and Tamor Sub-Basins of the Kosi Basin, Eastern Himalaya (2010–2025)

High-mountain glaciers are vulnerable to changing climate. Their size depends directly on temperature, snowfall, radiation balance, and surface processes like debris cover and ice movement. Because of this sensitivity, glacier extent changes are often used as reliable indicators of climate trends and environmental shifts in mountain areas [1-4]. Over the past few decades, many glaciers worldwide have noticeably shrunk in size and thickness [3, 5, 6]. This shrinkage has been particularly clear in the Himalayan region, where complex geography and monsoon-driven weather led to varied glacier responses [2, 7]. A decline in glacier size changes the seasonal meltwater supply, affects how sediments move, and alters the water systems that support ecosystems and human communities [8, 9].

The Himalayan Mountain belt is home to thousands of glaciers in various climates and elevations, constituting one of the most extensive ice reserves beyond the polar areas [5, 7]. Known as a key freshwater source outside the polar areas, this glacial landscape is vital for keeping rivers flowing year-round in major cross-border river systems [9]. Thus, understanding how glacier boundaries change over time is important for assessing long-term water supply, landform changes, and potential future environmental dangers.

Recent advances in Earth observation technology have made it possible to monitor glacier changes across vast and hard-to-reach areas [10, 11]. Satellite images combined with Geographic Information System (GIS)-based analysis enable accurate mapping of glacier edges and allow for comparisons of glacier size over different periods [6, 12]. Analyzing these changes over time helps reveal patterns of glacier retreat or stability and measures the extent of change at both basin and regional levels [3]. In this context, this study examines changes in glacier extent and associated spatial characteristics from 2010 to 2025 within the chosen basin. By outlining glacier boundaries using satellite data from various times, the study assesses recent trends in glacier size and enhances understanding of ongoing changes in the Himalayan cryosphere.

The current study focuses on the Arun and Tamor sub-basins of the Kosi River basin, a major Himalayan river system spanning parts of China (Tibet), Nepal, and northern India [7, 8]. The river starts in high-altitude glacial zones of the Tibetan Plateau and the eastern Himalaya, flowing southward before joining the Ganga River [8]. The basin has a complex hydrology, marked by large sediment flows and changing river shapes [7]. The basin consists of seven main tributaries: Arun, Tamor, Dudhkosi, Tamakosi, Sunkosi, Indrawati, and Likhu. Together, they form the Sapta Kosi system [8]. Among these, the Arun and Tamor sub-basins have considerable glacier cover at higher elevations and significantly contribute to runoff, especially during the melt season [8, 9].

Geographically, the basin showcases steep altitude differences, ranging from snow-covered mountain peaks to deep valleys and mid-mountain areas. Steep slopes, active tectonic activity, and widespread glacial features dominate the upper catchments, leading to high geomorphic activity and rapid sediment transport. The local climate is mainly impacted by the South Asian summer monsoon, which brings yearly rainfall [7]. High-elevation areas see long-lasting snow accumulation and cooler temperatures that favor glacier growth, while lower elevations experience warmer conditions that increase seasonal melting [9]. This mix of strong monsoon effects, steep terrain, and glacier-covered headwaters makes the Arun and Tamor basins ideal for studying changes in glacier size over time. Figure 1 shows the location of the study area in the Kosi Basin.

Multi-temporal Landsat images were used to analyze changes in glacier size within the study basin. Landsat Thematic Mapper imagery acquired during 2010–2011 and Landsat Operational Land Imager/Operational Land Imager-2 imagery acquired in 2025 were obtained from the U.S. Geological Survey EarthExplorer portal. Table 1 summarizes details of satellite datasets used. The basin boundary shapefile was obtained from the HydroSHEDS database. Google Earth Pro imagery was additionally used as reference data to support the visual interpretation and validation of the delineated glacier boundaries.

All selected Landsat scenes were obtained under cloud-free or almost cloud-free conditions ($<$10% cloud cover) in the post-monsoon and dry winter window (October–January). The amount of seasonal snow cover observed post monsoon is minimum due to the late ablation season, and the December 2010 and January 2025 scenes were systematically checked using band ratios and high-resolution reference images to make sure that temporary winter snowfalls did not hide or distort the glacier margins. Standard U.S. Geological Survey Level-1 Precision and Terrain Corrected products, which are geometrically and radiometrically calibrated, were used to keep co registration errors under 0.5 pixels.

To handle radiometric differences between the Landsat 5 Thematic Mapper (2010) and Landsat 8/9 Operational Land Imager (2025) sensors, raw digital numbers were converted into top-of-atmosphere reflectance before creating false color composites. Since boundary delineation was performed by manual digitization of geomorphological features and visual interpretation, rather than automated processing of single-band spectral thresholds, small differences in radiometric calibration across sensor generations did not lead to systematic bias in mapping [13].

Multi-temporal Landsat images were pre-processed with GIS techniques to map changes in glacier size between 2010 and 2025 which included band mosaicking, clipping and generation of false color composites to improve the visibility of glacier ice, debris-covered zones, and the surrounding land. Glacier boundaries were mapped through manual digitization, visual interpretation, and high-resolution imagery for reference.

In order to avoid overestimating the area of the glaciers, the seasonal snow cover was separated from the perennial glacier ice by means of the ratios of the shortwave infrared and near-infrared bands, since glacier ice has a lower albedo than fresh snow. Temporary snow patches that had broken off from the main glacier masses were excluded following a check against historical images at the end of the ablation season in Google Earth Pro. In the case of heavily debris-covered glacier tongues where the optical contrast with the adjacent moraines is small, the boundaries were determined by referring to topographic and geomorphological features such as the crests of the lateral moraine ridges, slope breaks, supraglacial ponds, and the meltwater outlets at the terminus [14].

The volume of the glacier was calculated by means of the empirical area–volume scaling relationship, which is given in Equation (1) [15]. In the present study, the scaling constant $c = 0.034$ and the scaling exponent $\gamma = 1.375$ were used. These parameters were derived theoretically from mass and momentum conservation principles applicable to valley glaciers. The scaling constant $c$ acts as a proportionality factor and is assumed to be identical for all individual glaciers within the study area [15]. To obtain an accurate estimate of total glacier volume, Equation (1) was applied to each glacier polygon individually. The resulting volumes were summed to yield the total glacier volume of the study area.

where, $V$ is the estimated glacier volume (km$^3$), $A$ is the glacier surface area (km$^2$), $c$ is the scaling constant ($c=0.034$), and $\gamma$ is the dimensionless scaling exponent ($\gamma=1.375$).

The resulting glacier polygons were then used to calculate glacier area statistics within a GIS environment. The volume estimations and some spatial calculations were also aided by Python-based geospatial tools. A multi-temporal comparison of glacier sizes was conducted to measure spatial changes and assess glacier retreat trends in the basin.

The glacier boundaries, which were determined by visual interpretation and manual digitization, are affected by uncertainties due to the spatial resolution of the sensors, terrain shading, debris cover, and the presence of seasonal snow cover. In order to assess the uncertainty involved in the mapping, a standard buffer technique using a half-pixel criterion (±15 m) was employed along the glacier perimeters [16]. The total length of the glacier perimeter was 3,626.07 km in 2010 and 3,527.61 km in 2025. The estimated area uncertainty was therefore calculated as ±54.39 km$^2$ (±5.65%) for 2010 and ±52.91 km$^2$ (±5.73%) for 2025. The delineated boundaries were cross-checked and confirmed by means of high-resolution images in Google Earth Pro to identify glacier termini and distinguish stable bedrock from debris-covered glacier tongues [16]. Figure 2 shows the methodological workflow for glacier mapping and change detection.

An analysis of glacier polygons using Landsat images shows obvious changes in the extent and volume of glaciers in the area between 2010 and 2025. The total glacier area decreased from 962.15 km$^2$ in 2010 to 923.89 km$^2$ in 2025, which is a total reduction of 38.26 km$^2$ (3.98%). The estimated glacier volume also showed a steady decline from 85.584 km$^3$ in 2010 to 81.885 km$^3$ in 2025. During this time, the total number of glaciers recorded increased slightly from 295 in 2010 to 297 in 2025 ( Table 2). This small increase in the number of glacier units is mainly due to local glacier fragmentation, caused by structural thinning and retreat of the terminus, which cut off narrow ice bridges between the tributary arms and the main glacier trunks and thus produced isolated smaller ice masses rather than indicating any real expansion of the glacier extent.

To illustrate the glacier-change trend over time, Figure 3 presents a graph of total glacier area and estimated glacier volume for the two observation years.

There is a general decline in the total glacier area and the estimated ice volume between 2010 and 2025, as shown in the graph. Figure 4 shows the spatial distribution of glacier extent for both years. Glacier cover is mainly found in high-elevation areas, especially in the northeastern and northern sections. The analysis shows that many glacier margins have pulled back from 2010 to 2025. This has led to the shrinking of valley glaciers and a reduction in the connection of previously joined ice bodies. The increase in the number of glaciers, along with the overall decrease in glacier area and volume, suggests that larger valley glaciers are breaking into smaller, separate ice bodies.

Figure 5 highlights retreating margins and changes in the shape of selected glaciers. In several places, glacier tongues are moving landward and narrowing, which indicates ongoing surface melting and thinning. These physical changes signal continued glacier thinning, negative mass balance, and ongoing changes in glacier shape. These trends are consistent with the long-term glacier retreat observed in the Himalayan region.

The observed reduction in glacier area (3.98%) and volume loss throughout the basin are in agreement with the documented trends of regional warming and changes in the precipitation patterns in the Eastern Himalayas [17], even in the absence of direct in-situ meteorological station records, which remain scarce across these high-altitude headwater zones. The rate of glacier area reduction observed in this study (about 0.26% per year) shows strong overall consistency with the baseline regional inventories, such as the Randolph Glacier Inventory 6.0 [11] and the multi-temporal datasets of the International Center for Integrated Mountain Development [8, 14], as well as with the documented patterns of glacier retreat in adjacent sub-basins. Comparable rates have been reported in the Dudh Koshi basin of Nepal ($\sim$0.2%–0.4% per year) [16] and the Sikkim Himalaya ($\sim$0.2%–0.3% per year) [14], while extensive terminus retreat and mass loss have also been documented in parts of the Bhutan Himalaya [14, 17].

Recent regional climate assessments point to an acceleration of warming at high altitudes and a reduction in solid snowfall during the monsoon transition period, which in turn increases summer ablation and reduces winter accumulation [17]. Although the response of individual glaciers is still influenced by local microtopography and debris thickness [4, 14], regional atmospheric conditions remain the main cause of the decline in the cryosphere across the basin [17]. Further ice loss could modify the seasonal discharge patterns and increase the vulnerability to downstream hazards such as glacial lake outburst floods before leading to a long-term decrease in baseflow availability [8, 17]. These findings highlight the need for ongoing glacier monitoring using remote sensing and GIS methods to understand how the cryosphere is responding to environmental changes.

The present study examined changes in glacier size in the Kosi River Basin using Landsat satellite images and GIS analysis. The results show a decline in the total glacier area from 962.15 km$^2$ in 2010 to 923.89 km$^2$ in 2025, which represents a reduction of about 3.98%. The estimated glacier volume also decreases, indicating ongoing ice loss in the basin. Despite the overall reduction in glacier area, a small increase in the number of mapped glaciers from 295 to 297 suggests that larger glaciers are breaking up into smaller, more isolated pieces. Spatial analysis shows significant retreat of glacier margins and changes in their shape, especially in high-elevation areas. These changes highlight how sensitive Himalayan glaciers are to local climate changes, such as rising temperatures and shifts in precipitation.

The ongoing reduction in glacier size and volume may affect the long-term availability of meltwater and hydrological processes downstream and also increases the risk of potential hazard in the Kosi River system. Therefore, continuous monitoring of glacier changes through remote sensing and geospatial methods is important for improving our understanding of cryospheric responses to climate change and supporting future water resource assessments in Himalayan basins.