Himalayan glaciers are normally difficult to monitor through field observations because of highly rugged and extremely inaccessible mountainous terrain. Thus, using Landsat data (MSS, ETM and TM), changes in glacier area, length, and debris cover have been delineated in the Baspa basin, which is a highly glacierized sub-basin of the Satluj River in the western Himalaya. Out of the total 109 glaciers inventoried through Landsat TM imagery (2011), 36 glaciers were found to be heavily debris covered (32.5 ± 2.0%). A shrinkage in glacier area of 41.2 ± 10.5 km2 (i.e., 18.1 ± 4.1%) at a rate of 1.18 ± 0.3 km2 a-1 from 1976 (227.4 ± 9.4 km2) to 2011 (186.2 ± 3.7 km2) has been recorded. The overall glacier retreat studied for 33 glaciers varied from 3.3 ± 0.03%, that is, 0.87 ± 0.06 km at a rate of 17.2 ± 1 m a-1 to 30 ± 6.6%, that is, 0.60 ± 0.04 km at a rate of 24.8 ± 0.2 m a-1. Consequently, the debris cover has increased by 23.5 ± 1.4 km2 (16.3 ± 3.8%) from 1976 to 2011. Overall, the clean, small sized, low-altitude glaciers with south to southwest aspect and relatively steep slope have lost maximum area, which indicated a major control of these factors on the glacier changes. Simultaneously, a trend estimation of observed climatic data (1976/1985–2008) of three meteorological stations (Sangla, Rakcham, and Chitkul) using Mann Kendall test, Sen's Slope estimator and linear regression test revealed an increase in temperature and rainfall while a decline in snowfall. Importantly, the Tmin has increased significantly at 95% confidence level during all the studied periods. The mean annual Tmin and Tmax indicated a rising trend at a rate of 0.076 and 0.071 °C a-1. Thus, the changes in temperature and precipitation may be the major causes of accelerating the glacier ablation. The higher area changes (53.0 ± 0.4%), of small glaciers <0.5 km2 mark their sensitivity to climatic changes especially rising temperature. Under the warming climate, formation and progressive expansion of glacial lakes is expected because of the glacier recession in the basin. For instance, the Baspa Bamak Proglacial Lake at the snout of Baspa Bamak glacier has expanded continuously from 2000 onward.
In the Hindu-Kush-Himalayan (HKH) region, an area approximately 40,800 km2 is covered by the glaciers (Bolch et al., 2012) with a total of 9675 glaciers in the Indian Himalaya (Raina and Srivastava, 2008). Many glaciers are classified as compound valley type glaciers, which are formed by the convergence of two or more tributary glaciers (Müller et al., 1977; Raina and Srivastava, 2008). The region is the birthplace of some of the largest rivers in Asia, such as Ganga, Brahmaputra, and Indus (Immerzeel et al., 2010; Thayyen and Gergan, 2010) and provides a huge freshwater resource. The glaciers influence the overall runoff in lowland rivers; recharge river-fed aquifers; provide water for hydropower, agriculture, and ecosystems; and eventually contribute to sea-level rise (Dyurgerov and Meier, 2005). However, a majority of the earlier studies showed that the Himalayan glaciers are diminishing and have been in a state of recession since the 1850s (Bhambri and Bolch, 2009; Bahuguna et al., 2014; Mir et al., 2014a).
In this region, it is reported that the glacier retreat begun with more gradual climate warming since the Little Ice Age (LIA) (Mayewski and Jeschke, 1979; Bhambri and Bolch, 2009). During recent decades, accelerated rates of glacier changes have been reported in the Himalayan region (Zemp et al., 2009; Chand and Sharma, 2015). Several studies on the Himalayan glaciers have indicated a receding-trend over the past few decades (Kulkarni et al., 2007; Bolch et al., 2008; Chand and Sharma, 2015). Overall, the shrinkage of Himalayan glaciers since the 1970s (Yao et al., 2012) has been attributed to rising temperatures and decreasing precipitation (Yao et al., 2012; Yang et al., 2014; Mir et al., 2014a). A number of studies suggest the primary source for melt energy and fluctuations in glaciers through the years is mainly solar radiation, temperature, and precipitation (Greuell and Smeets, 2001; Oerlemans, 2005; Mir et al., 2014a, 2014b). Glaciers being physically complex and dynamic systems are sensitive to climate change (Benn and Evans, 2010) and are therefore considered as key indicators for assessing climate change in the Himalayan region, especially in areas without climatic observation instrumentation facilities (Barry, 2006; Bolch et al., 2012; Mir et al., 2014a).
In contrast to the global trend of declining glaciers, the Karakoram glaciers are reported to be in a state of expansion (“Karakoram Anomaly”) with many glaciers having stable fronts since 2000 (Hewitt, 2005; Hewitt,2011; Kääb et al., 2012; Bolch et al., 2012; Bhambri et al., 2013; Gardelle et al., 2013; Bahuguna et al., 2014; Mir and Majeed, 2016). This irregular behavior of Karakoram glaciers in general is attributed to local/regional topography (Haeberli,1990), the local/regional climatic system (Kargel et al., 2005), glacier hypsometry (Furbish and Andrews, 1984), the characteristics and thickness of supraglacial debris cover on the glacier surface (Bolch et al., 2008; Scherler et al., 2011), the glacier size and accumulation area ratio (AAR) (Kulkarni et al., 2007), contributions from tributary glaciers (Nainwal et al., 2008a), and their geometrical/morphological properties (Mehta et al., 2014). Similarly, the phenomenon of the stable or advancing state of Karakoram glaciers (e.g., Kääb et al., 2012; Yao et al., 2012; Gardneret al., 2013; Neckel et al., 2014) is considered to result from increased snowfall and cooling at high altitudes (Hewitt, 2005; Yang et al., 2014; Wiltshire, 2014).
Generally, the Himalayan glaciers are losing an average of 0.4% of area per year (Bolch et al., 2010, 2012; Bhambri et al., 2011; Kamp et al., 2011; Kulkarni et al., 2011, Chand and Sharma, 2015). However, in Xiongcaigangri and neighboring areas, glaciers have retreated less (Brahmbhatt et al., 2015; Li et al., 2015). But, in the Tirungkhad basin located in western Himalaya, a higher shrinkage of 29.1 km2 (i.e., 26.1% reduction) has been reported from 1966 to 2011 (Mir et al., 2014a, 2014b). Similarly, in the Zanskar valley of western Himalaya, a shrinkage of 212 glaciers resulted in a loss of 57 km2 (i.e., 8% reduction); however, this was partially offset by a 42 km2 (or 6%) increase in other glaciers. Overall, a glacial area decrease of 15 km2 (or 2%) between 1962 and 2001 has been reported (Ghosh et al., 2014), of which 6.5 km2 (or 5.6%) occurred from 1980 to 2010 in the Shyok Basin (Bajracharya et al., 2015).
Like other regions, most of the glaciers in the Himalayan region are retreating because of accelerated global warming during the past century. Over the past century, an increase of about 1.7 °C in temperature and decreasing trend in precipitation (Bhutiyani et al., 2009; Mir et al., 2015a, 2015b) are reported in the western Himalayan region. However, the response of a glacier to climate change depends on its geometry and on its climatic setting (Oerlemans, 2005). The debris cover also affects glacier response to climate change by altering surface ablation rates and spatial patterns of mass loss (Nakawo et al., 1999; Benn and Lehmkuhl, 2000; Kirkbride, 2000; Benn et al., 2003; Scherler et al., 2011). In the Himalaya, debris-covered glaciers are widespread, influencing regional-scale patterns of glacier changes and ice dynamics (Bolch et al., 2008; Quincey et al., 2009; Scherler et al., 2011). The changes in Himalayan glaciers are expected to have a significant impact at the catchment scale, including long-term reduction of water resources. The continuous recession of these glaciers can lead to increased probability of glacier lake outburst flood (GLOF) events because of the development of lakes behind weak moraine dams (Yamada, 1998; Richardson and Reynolds, 2000). During recent years, hazardous situations in the downstream areas as the result of sudden collapse of ice dams have also been given ample attention (e.g., Hewitt, 1982; Raina and Srivastava, 2008).
With the reported impact of climate change on the glaciers in recent decades (IPCC, 2013), the Himalayan cryosphere system has received greater attention for monitoring and observations (Kääb et al., 2003; Wang et al., 2014). The changes occurring in the Himalayan glaciers are being discussed extensively, as it can affect the supply of water to a large number of people dependent on waters of the major rivers in the region. Similarly, during the recent past, a number of studies have been carried out for the Baspa basin in the western Himalaya. For instance, Kulkarni and Alex (2003) have mapped 19 glaciers based on Survey of India (SOI) topographic maps and Indian Remote Sensing, Linear Imaging Self Scanning Sensor (IRS-LISS-III) satellite data from 1962 to 2001. Retreat of 8 glaciers of the basin were also estimated from 1962 to 1998 using the SOI topographic maps and high resolution satellite stereo data from Indian Remote Sensing Satellite-1C (Kulkarni and Bahuguna, 2002). In another study, 19 glaciers have been mapped and the mass balance studies have been carried out using WiFS (Wide Field Sensor) images of IRS for 2001 and 2002 (Kulkarni et al., 2004). Similarly, Gaddam et al. (2016) carried out studies on the retreat and mass loss of 19 glaciers based on the SOI maps and Landsat datasets. Deota et al. (2011) has generated geomorphic information on the retreat of 22 glaciers based on IRS 1D LISS III + PAN-merged data of 2001, and Resourcesat, LISS III data of August 2005. In addition, the Geological Survey of India (GSI) has published an inventory of 89 glaciers for the Baspa basin (Raina and Srivastava, 2008) based on the old SOI maps (1962) and satellite imageries. For the Baspa basin, the glacier outlines are also provided in the Randolph Glacier Inventory (RGI) v 5.0 (Pfeffer et al., 2014). The Baspa basin is a part of northwestern India, for which the RGI v 5.0 has utilized the glacier inventory data compiled by the GlobGlacier project of the European Space Agency (ESA) (Paul and Andreassen, 2009) from Landsat ETM+ of 2001 and ALOS PALSAR data (Frey et al., 2012).
However, for the entire Baspa basin, there is no study to date that reports recent glacier changes and their correlation with meteorological variables (e.g., temperature and precipitation) and other related factors such as size, topography, and debris cover. Thus, the present study is conducted for (1) preparation of a current glacier inventory and study of glacier changes, (2) analysis of various meteorological and related factors, and (3) the relationship of these factors with the glacier changes.
The Study Area
The Baspa river is a major tributary of the Satluj River basin. Generally, the Satluj basin has a diverse climate as a result of the wide range of altitudes and precipitation patterns as it covers outer, middle, and greater Himalayan ranges. Around the lower part of the basin, tropical and warm temperate climate prevails, whereas the middle part has a cold temperate climate, and the upper part is very cold. In the upper part of the basin, most of the precipitation is produced by the westerly weather disturbances (Dimri and Mohanty, 2007; Mir et al., 2015a) that bring snow during winter months.
The Baspa river basin is located at the middle of the Satluj basin in Kinnaur District, Himachal Pradesh, India, and spreads from 31°05′ to 31°30′N latitudes and 78°00′ to 78°50′ E longitudes (Fig. 1). The Baspa river originates from Arsomang and Baspa Bamak glaciers and travels ∼72 km in west-northwest direction through an elongated ∼75-km-long and ∼18-km-wide valley before joining with the Satluj River at Karcham at an elevation of 1793 m a.s.l. The valley is narrow in the lower stretch but in upper reaches it is wide and characterized by deglaciated wide open terraces (Raina and Srivastava, 2008). The coarser clastics derived from the valley slopes occur along the valley bottom. This basin occupies ∼1100 km2 of area in the highly glaciated (17%) basin. Most of the glaciers are located on the northern slope of the Pir Panjal mountain range of the Greater Himalaya (Raina and Srivastava, 2008). The elevation of the Baspa basin ranges from 1793 m to 6442 m a.s.l. The main glaciers present are Baspa Bamak, Shaune Garang, Jorya Garang, and Karu. The glaciers are of varied sizes and vary in form from cirque to compound valley glaciers. The stream flow is mostly generated from snow and glacier melt runoff with a significant contribution as well from monsoon rainfall. The mean winter temperature varies from -7.4 °C (Tmin) to 5.5 °C (Tmax); the total monsoon precipitation is very low (136 mm) as observed at Rakcham station. In socioeconomic terms, the Baspa Basin is important as many new mini- and micro-hydropower stations are being planned in this basin. For example, the glacier and snowmelt runoff plays an important role in the power generation of Karcham Wangtoo Hydro-electric plant, with a capacity of 300 MW (Gaddam et al., 2016). Therefore, it is important to have regular assessment of glaciers of the Baspa basin.
Satellite products are the major source of data for studying glacier changes all over the world (Barry, 2006; DeBeer and Sharp, 2007; Racoviteanu et al., 2008a; Paul and Andreassen, 2009; Shangguan et al., 2007). Multi-spectral and multi-temporal satellite data offer abundant potential to monitor the glaciers at regular intervals (Racoviteanu et al., 2009; Bolch et al., 2010; Paul et al., 2013; Mir et al., 2014a). Nonetheless, the conventional field surveys which are highly advisable for the assessment of glaciers are often laborious in high mountainous regions such as the Himalaya (Mir et al., 2014a, 2014b; Pratap et al., 2015a). In this study, the satellite data of the Landsat series, that is, MSS (60 m), TM/ETM+ (30 m) from the National Aeronautics and Space Administration (NASA) distributed through the University of Maryland's Global Land Cover Facility (GLCF), and the United States Geological Survey's (USGS's) Earth Resources Observation and Science (EROS) Center ( http://earthexplorer.usgs.gov/) have been used to investigate the glacier changes during 1976–1992, 1992–2000, 2000–2006, 2006–2011, and 1976–2011 periods. In addition, the topographic maps of the SOI (1:50,000 scale) and images of IRS-LISS-III (23.5 m) from National Remote Sensing Centre of Indian Space Research Organization have also been utilized as a source of additional information mainly for accurate mapping of glacier outlines. Because of a large uncertainty in their reliability (Bhambri et al., 2011; Mir and Majeed, 2016), the old SOI maps were excluded from the change detection studies and therefore, the MSS (1976) scene was used as the base image for glacier change studies. However, wherever required, the SOI maps (1962) were used to validate and support the mapping of the glaciers from the MSS scene. During the 1990s, two scenes of Landsat TM (i.e, 1992, 1994) were available in which the TM scene (1992) was used to extend the time span of observation and because of its little snow cover. But, wherever necessary, the TM (1994) scene was used to assist the mapping of glaciers from the 1992 scene. The ETM+ (2006) scenes are generally affected by scan line errors because of Scan Line Corrector (SLC) failure. However, our study area (Baspa basin) was fortunately free of scan line errors. Further, the available LISS III scene (2007) was used to check and maintain the accuracy of glacier extents delineated from the ETM+ scene through visual inspections. Using the Landsat TM (30 m) of 2011, a current and updated glacier inventory has been generated for the basin. The ASTER DEM (30 m spatial resolution) from the USGS was used for semiautomatic delineation of drainage basins and extraction of topographic parameters (such as elevation, aspect, slope, etc.). In addition, the glacier outlines available for Baspa basin at RGI v 5.0 generated for South Asia West ( http://www.glims.org/RGI/rgi50_dl.html) were used for comparative analysis. The details of data used are given in Table 1.
Data sources available for the study area.
The long-term meteorological data were used to investigate the changes in climatic conditions in the basin. The climatic data includes temperature (Tmin and Tmax) and precipitation procured from Bhakra Beas Management Board (BBMB). The precipitation data of snowfall as well as rainfall are available separately. In the Baspa basin, the snowfall is measured by melting the daily snow accumulations in the standard rain gauges and is represented as snow water equivalent (SWE) in millimeters. The snowfall was available from three stations (Sangla, ∼2742 m a.s.l.; Rakcham, ∼3154 m a.s.l.; and Chitkul, ∼3444 m a.s.l.) whereas the temperature and rainfall data were available from Rakcham station only. These data are available from 1985 to 2008 from Rakcham station, whereas from Sangla and Chitkul stations, the data are available from 1976 to 2008 (Fig. 1).
As mentioned above, this study is carried out using Landsat data series. But, additional information was obtained from the SOI topographic maps and IRS images to map the glaciers accurately. Using the GIS platform, the SOI maps were scanned, georeferenced, and primarily utilized for the registration of the 2011 Landsat TM image based on about >100 ground-control points evenly distributed within the whole basin. Subsequently, all other datasets were coregistered with the 2011 Landsat TM imagery. The root mean square error (RMSE) for MSS scene and other TM scenes was found to be 35 m (0.58) and 15 m (0.5), whereas the ETM+ scene registered an error of 20 m (0.68). The LISS III images were orthorectified based on the ASTER DEM before coregistration to the 2011 TM image and showed a registration error of 20 m (0.67). Furthermore, the images were geometrically rectified to the same projection system of WGS 84 UTM Zone 43.
From the Landsat MSS scene, the glaciers were delineated based on the false color composite (432), whereas, using the ETM+ and TM data images, the glaciers were first automatically delimited using band ratios of near infrared and shortwave infrared channels (i.e., TM4/TM5) of the images (e.g., Paul and Andreassen, 2009; Paul et al., 2015). In this way, a clean glacier-ice mask was generated in a binary image using threshold values of 2.5 (2000) and 2.2 (2006) for ETM+ images. For TM images, threshold values of 2.4 (1992) and 2.3 (2011) were used. These binary images were then converted into vector format. However, several misclassified glacier areas such as water bodies, shadows, and isolated rocks were present on the vector data, which were then easily eliminated by post-processing (i.e., manual exclusion/inclusion of wrongly classified/missing areas) to discern the glaciers prior to data retrieval. Nevertheless, ASTER GDEM has also been used for correction of glacier outlines through better visualization of the glacier extents and mapping of ice divides. Wherever possible, visual inspection of the glaciers, ice divides, and associated glacial features on Google Earth (≤⃒5 m) was used as an additional information source for mapping, and any mismatch found was eliminated manually (Paul et al., 2013). From SOI maps the glaciers were delineated by onscreen digitization.
Eventually, after deriving the information on glacier parameters, such as area and length, an overlay analysis of the vector data of glaciers for each study year was carried out by using the ArcGIS 10 spatial module. The minimum size of glaciers for inclusion into the inventories (e.g., 2011 inventory) was set at 0.01 km2. However, the comparison of glacier areas was restricted to those larger than 0.10 km2 in order to reduce the errors that might have been introduced because of the seasonal snow cover and the expected larger uncertainty of delineation of small glaciers, especially from the coarse-resolution Landsat MSS scene. A more meaningful approach to comparing glacier's dynamic states is to report area change rates in relative terms (i.e., percentage change in area per unit time). In this study, the percentage change in area per unit time (e.g., % change/year) was calculated by dividing the percentage change in area of a glacier by the number of years between two observations (Zemp et al., 2014).
Moreover, using the ASTER GDEM, the information on the terminus elevation, average/mid-elevation, median elevation, slope, and aspect were determined using Geographic Information System (GIS) tools as suggested by Paul et al. (2009). The mean slope data of each glacier were estimated as an average value of all pixels falling in their respective classes intersected by the glacier outlines. The mean aspect was estimated by dividing the mean sine by the mean cosine of aspect across all cells (Paul et al., 2009; Guo et al., 2015). The average mid-altitude was derived as the arithmetic mean of all glacier pixel altitudes, whereas the median elevation was extracted as the 50th percentile of the cumulative number distribution of pixel elevations as suggested by Guo et al. (2015). The glaciers were also classified into various altitude zones on the basis of the ASTER DEM. For this purpose, the glacier outlines of different periods under investigation were overlaid upon the DEM and subsequently clipped using ArcGIS software. The glacier area for each altitude zone was calculated on the basis of number of pixels falling in each zone.
The glacier length was manually derived along the central flow line from the snout to the highest point in the accumulation zone along the maximum length of the glaciers (Paul et al., 2009). Moreover, a band of stripes with a distance of 50 m between each stripe was drawn parallel to the main central flow line of the glaciers. The average length was calculated from the intersection of the stripes with the glacier outlines in terms of its retreat for comparison with average length change (Bhambri et al., 2012; Mir and Majeed, 2016).
The terminus/snout of glaciers was determined through visual interpretation technique, based on shape, size, pattern, tone, texture, and association (Kulkarni and Buch, 1991). The identification of the glacier snout was especially difficult in cases of heavily debris-covered glaciers. Therefore, to correctly map the snout location, several signs of movement (identified based on overlays of multi-temporal images), breaks in surface slope, spectral color differences, texture of glacier surface and periglacial area, emerging meltwater streams at the end of the terminus, braided streams, ice-wall shadow, presence of moraine-dammed supraglacial ponds/lakes, and proglacial morphological features were used as additional aids and diagnostic tools (e.g. Bhambri et al., 2013; Mir et al., 2014a). The debris-covered glaciers were delineated manually based on the above discussed indicators. The manual delineation of debris-covered glaciers is usually considered more accurate than automated methods (Raup et al., 2007), which is a major source of error in glacier mapping (Bhambri et al., 2011; Bolch et al., 2008)
Uncertainty in the Study
The quantification of error is crucial to ascertain the accuracy and significance of results. In general, the glacier outlines derived from satellite datasets at various spatial and temporal resolutions are subject to various degrees of uncertainty (Racoviteanu et al., 2009; Paul et al., 2013). Similarly, the positional accuracy of the glacier snouts is influenced by the sensor resolution (Williams et al., 1997) and coregistration error (Hall et al., 2003). However, in this study, the utilized Landsat (e.g., TM/ETM+) data sources were having similar spatial resolutions, except the Landsat MSS image. Therefore, errors of coregistration and glacier boundary delineations were considered to be the main factors that might have resulted in different levels of accuracy. The accuracy of the glacier outlines may have been further enhanced by the debris cover as it strongly influences the glacier outlines (Paul et al., 2013). The most accurate way to assess glacier outlines would be to use high-resolution images (Paul et al., 2013), but such data were not available for our study region. Still, wherever possible (e.g., toward western parts of the basin), the higher-resolution images from Google EarthTM were inspected in conjunction with the Landsat images to maintain accuracy.Bhambri et al., 2011; Shahgedanova et al., 2014). For instance, in this study basin, debris-covered glaciers cover a larger area than debris-free glaciers. Therefore, it is assumed that the uncertainty due to debris cover will be captured by taking an additional error of 7.5 m (Bhambri et al., 2011). Overall; the terminus accuracy was estimated to be ±64.92 m for TM images, ±69.92 m for the ETM+ image, and ±127.35 m for the MSS image. The terminus error of ±127.35 m for the MSS image is almost similar to a terminus uncertainty of ±134.13 m for the MSS image reported by Shukla and Qadir (2016).
The measurement of uncertainty of area extent (Uarea) for each glacier was determined by the buffer method (Granshaw and Fountain, 2006). It is a widely used and recommended approach to estimate the uncertainty using a buffer of up to one pixel size to the glacier margin (Granshaw and Fountain, 2006; Bolch et al., 2010; Paul et al., 2013). In this study, an area of buffer of each debris-covered and clean/debris-free glacier having a width equal to the digitizing error (RMSR) was created, and the uncertainty was calculated as an average ratio of the original glacier areas to the areas with a buffer increment. The overall uncertainty for each dataset was determined from the uncertainty of debris-covered and debris-free glaciers estimated separately. For debris-covered glaciers, a buffer distance of 15 m for TM/ETM+ images and 30 m for MSS image was used, whereas, for debris-free glaciers, buffer distances of 10 m and 20 m were used. Overall, an average uncertainty of ±0.020 km2 (2.0%) for TM (2011, 2000) and ETM+ (2006), ±0.022 km2 (2.2%) for TM (1992), and ±0.041 km2 (4.1%) for MSS (1976) images was determined. The uncertainty (4.1%) of the MSS image is similar to an uncertainty of 4.7% reported by Shukla and Qadir (2016). Furthermore, the uncertainty was also evaluated by comparing the glacier outlines derived from the Landsat ETM+ (30 m) of 2006 and those derived from an IRS image (23.5 m) of the same year for a sample of six selected large glaciers as suggested by Paul et al. (2013). The resulting uncertainty was found to be ±2.6%. The overall uncertainty values estimated are well within the previously reported acceptable ranges (Paul et al., 2002; Racoviteanu et al., 2008a; Bolch et al., 2010; Bhambri et al., 2011; Mir et al., 2014a).
To estimate uncertainty of glacier area change for all the time intervals of this study, the average uncertainty of area extent of all the glaciers for various data scenes was used. The area change uncertainty is estimated according to the law of error propagation using a given formula as
The trend estimation of the climatic variables such as temperature and precipitation (i.e., snowfall and rainfall) was carried out using the Mann—Kendall test (Mann, 1945; Kendall, 1975), Sen's slope estimator, and simple linear regression analysis. For this purpose, the daily climatic data were classified into annual, winter, pre-monsoon; monsoon, and post-monsoon seasons as per the prevailing climatic conditions of the region (Mir et al., 2015a, 2015c). The months of November, December, January, and February (NDJF) were classified as winter season. Similarly, the months of March, April and May (MAM) stand for pre-monsoon season, July and August (JA) stand for monsoon season, and September and October (SO) stand for post-monsoon season. Further, the months from October to March (ONDJFM) were defined as accumulation season followed by ablation season from April to September (AMJJAS). This particular division of the months is used because the nourishment of the glaciers takes place mainly during winter season from October to February, which in turn is because of the mid-latitude westerlies that have a dominant role in this part of the Himalaya (Benn and Owen, 1998). A significant accumulation also takes place during the month of March, but the ablation starts from April onward up to September and reaches to its peak during the month of August. The Zs (Z statistics) test statistic was used as a measure of significance of the trend. Positive values of standard normal variable Zs indicated an increasing trend, whereas the negative Zs valves showed decreasing trends. In this analysis, the trend significance was tested at 95% with Zs = 1.96 confidence level. The magnitude of slope (Qi), that is, change per unit time, was estimated using a simple non-parametric procedure developed by Sen (1968). Moreover, for the data without any trend, the probability values should be close to 0.5. The probability of the trend is represented by the p-value.
The inventory generated for the year 2011 in Baspa basin consists of 109 glaciers with an area of 187.0 ± 3.7 km2 (Fig. 1). The glacier cover is estimated to be 17% of entire basin area. The mean glacier size is about 1.7 km2 with the largest glacier covering an area of 32.3 ± 0.6 km2 (Baspa Bamak). Moreover, the glacier's sizes were divided into five categories as <0.5, 0.5–1, 1–5, 5–9, and >9 km2. The glacier size class >9 km2 has only 3 glaciers with an area of 69.8 ± 1.4 km2 in the basin constituting 37.3 ± 0.7% of the glacier area. Maximum number of glaciers (55) are in the size class of <0.5 km2 covering 10.7 ± 0.2 km2 (5.7 ± 0.1%) of area followed by 22 glaciers in the size class of 0.5–1 km2 covering an area of 13.9 ± 0.3 km2 (7.5 ± 0.1%). There are only 29 glaciers in the size classes from >1 km2 to <9 km2 that contribute 49.5 ± 1.0% to the total glacierized area of the basin.
The glacier termini are located around an altitude varying from 4244 to 5762 m a.s.l., whereas the mean/mid-altitude of the glaciers ranged from 4536 to 5842 m a.s.l. with an average mid-altitude of 5135 m a.s.l. The mid-altitude of 7 glaciers is greater than 5600 m a.s.l., and 1 (one) glacier reaches lower than 4500 m a.s.l. About 36 glaciers have debris-covered ablation zones where 20.5 ± 1.5% of their area is covered by debris. Moreover, the median-altitude of the glaciers, which is widely used for the estimation of long-term equilibrium line altitude (ELA) based on topographic data (Braithwaite and Raper, 2009), is slightly higher than the mid-altitude, with ranges from 4482 to 5872 m a.s.l. with an average of 5197 m a.s.l. Based on median-altitude of the glaciers, the distribution of the glacier area of total glaciers, debris-covered glaciers, and debris-free/clean glaciers against elevation (hypsometry) are shown in Figure 2. Most of the glaciers (78) have an area of 166.8 ± 3.3 km2, which is 89.2 ± 1.8% of the total glacierized area distributed between 4800 and 5500 m a.s.l. altitude. The clean glaciers cover a small part (29.5 ± 0.6 km2) of this range, whereas a major part (i.e., 137.7 ± 2.8 km2 which is 73.7 ± 1.5% of total glacierized area) is covered by debris-covered glaciers. Overall, the clean glaciers are located around higher elevation (above 5400 m a.s.l.) in the basin.
Most of the glaciers (74) covering about 75.0 ± 2.0% of the area have a northerly orientation varying from northwest to northeast as shown in Figure 3, part a. About 27 glaciers have a northwest (337.5°-22.5°) mean orientation followed by 25 glaciers with a northeast (22.5°–67.5°) and 22 glaciers with a north (337.5°-22.5°) mean orientation. The southeast (112.5°–157.5°) and west (247.5°–292.5°) orientations have 5 glaciers each. Most glacier areas, that is, 47.8 ± 1.0 km2 (25.6 ± 0.5%), are distributed on a northwest slope, whereas a small area of 14.4 ± 0.3 km2 (7.7 ± 0.2%) is located on the southeast slope. About 27.5 ± 0.6 km2 (14.6 ± 0.3%) of the glacier area is distributed on southern slope (157.5°–201.5°). Similarly, the mean slope of glaciers was found to be 18°, with a range from 5° to 35° (Fig. 3, part b). There are about 7 glaciers with a mean slope of 5° covering 47.4 ± 0.9% (88.2 ± 1.8 km2) of glacier area and 2 glaciers with a mean slope of 35°covering a small area of 0.6 ± 0.01% (1.1 ± 0.02 km2). Overall, 75.6 ± 1.5% (140.8 ± 2.8 km2) of the glacier area is covered by about 30 glaciers with a mean slope of <18°.
The glaciers of the basin were also classified into the different types as given in Table 2. The larger and compound valley glaciers (e.g., G-49, G-80) tend to extend down to lower altitudes, whereas the smaller glaciers (e.g., cirques) have high-altitude termini. The altitude also varies with glaciers size, that is, larger glaciers have a larger altitude range, and vice versa (Fig. 3, part c). The glacier length ranged from 0.30 ± 0.01 km (G-1) to 17.8 ± 0.8 km (G-49). Glaciers having a length of 2 km are prevalent in the basin (Fig. 3, part d).
Changes in Glacial Parameters
The glaciers have periodically lost area during the different time periods from 1976 to 1992, 1992 to 2000, 2000 to 2006, 2006 to 2011, and 1976 to 2011 (Table 3; Fig. 4). The glacier retreat map is shown in Figure 5. There were 103 glaciers in 1976 that covered an area of 227.4 ± 9.4 km2 (mean glacier area: 2.2 ± 0.1 km2) above a set size of 0.1 km2 selected for the study of glacier evolution in the basin. The area of the largest glacier (Baspa Bamak, G-49) was 35 ± 1.4 km2, whereas the smallest glacier size was 0.12 ± 0.005 km2 (G-53). Based on the 2011 inventory, the area of the observed 97 glacier units was 186.2 ± 3.7 km2 (mean glacier area: 1.9 ± 0.04 km2), in which the maximum size of the largest glacier was 32.3 ± 0.6 km2 (G-49) and the size of the smallest glacier was 0.10 ± 0.002 km2 (G-53). The glacier area has shrunk from 227.4 ± 9.4 km2 (1976) to 207.8 ± 4.6 km2 in 1992, followed by 198.1 ± 4.0 in 2000, 192.1 ± 3.8 km2 in 2006, and 186.2 ± 3.7 km2 by 2011. Overall, a total area loss of 41.2 ± 10.5 km2 (18.1 ± 4.1%) at a rate of decrease of 1.18 ± 0.3 km2 a-1 (0.51 ± 0.01%) was observed over a period of 35 years from 1976 to 2011.
Classification and other attributes of the glaciers based on Landsat-TM image (2011) and ASTER GDEM in the Baspa basin.
Total glacier areas, areas of debris covered, debris-free/clean glaciers, and corresponding changes during different time intervals from 1976 to 2011 in Baspa basin. The total number of glaciers and area for the studied years are based on a size of 0.10 km2 considered for the glacier evolution study in the basin.
During different time intervals, the percentage change in area shrinkage indicated a declining pattern from 1976 to 2011. For instance, a loss of 8.5 ± 4.1% was observed from 1976 to 1992 followed by 4.7 ± 2.2% from 1992 to 2000, 3.0 ± 2.0% from 2000 to 2006, and 3.1 ± 1.9% from 2006 to 2011. The rate of glacier loss almost indicated a similar pattern with minor variations. For instance, during 2000–2006, the rate of glacier change has been small, that is, 1.0 ± 0.9 km2 a-1 (0.50 ± 0.17%) than in the preceding time intervals of 1992–2000 and 1976–1992. During 1992–2000, a large glacier loss at a rate of 1.21 ± 0.76 km2a-1 (0.58 ± 0.13%) was found. The rate of glacier recession of 1.18 ± 1.0 km2 a-1 (0.61 ± 0.21%) was also higher during 2006–2011. Moreover, to understand the glacier changes for past three decades, that is, 1980–1990, 1991–2000, and 2001–2010, the rate of loss derived from available data sets of different dates was normalized and translated into decadal (10 years) changes. Depending upon the available data sets during a period of 10 years, the derived rate of loss for any date was supposed to be constant up to the next date and simply multiplied by 10 (i.e., number of years in a decade). For instance, a loss of 12.2 km2 during the 1981–1990 decade was derived by multiplying the rate of loss of 1.22 km2 a-1 obtained from MSS (1976) by 10. Overall, the analysis revealed a similar loss of 12.2 ± 6.5 km2 of glacier ice during last two decades, that is, from 1981 to 1990 and 1991 to 2000. However, the glacier ice loss has declined to 10.7 ± 9.4 km2 during a recent decade, from 2001 to 2010, which may possibly be because of the lower recession rate of 1 ± 0.91 km2 a-1 during 2000–2006 as discussed above. Among the glaciers, the loss ranged from 5.9 ± 0.25% (G-37) to 76.6 ± 0.20% (G-102) during 1976–2011. Overall, a mean decrease in glacial area at the rate of 0.51 ± 0.01% a-1 from 1976 to 2011 has been observed, which indicated a high and rapid glacial retreat in the basin.
Moreover, the area loss in the size classes of <0.5 km2, 0.5–1 km2, 1–5 km2, 5–9 km2, and >9 km2 corresponds to 53.0 ± 0.4%, 29.2 ± 0.3%, 20.8 ± 1.3%, 9.4 ± 0.8%, and 8.1 ± 1.5%, respectively. On an average, the small glaciers lost more of their surface area than large glaciers. Similar observations have been reported earlier in the other parts of Himalaya (Bolch et al., 2012; Mir et al., 2014b).
Out of 103 glaciers identified during 1976, about 10 glaciers melted out and disappeared during the intervening period until 2011. Out of the 10 disappeared glaciers, 5 glaciers (0.16 ± 0.003–0.8 ± 0.02 km2) vanished between the years from 1976 to 1992, whereas 3 glaciers (0.17 ± 0.003–0.65 ± 0.01 km2) vanished from 1992 to 2000, and 2 (0.12 ± 0.002–0.22 ± 0.004 km2) from 2000 to 2006. Simultaneously, four glaciers (i.e., G-22, G-25, G-95, and G-103) have fragmented into two individual glaciers during past four decades from 1976 to 2011. Overall, the number of glaciers decreased from 103 to 97 from 1976 to 2011 (Table 3).
An analysis of length changes in the present basin revealed a continued recession (Fig. 4), in common with most of the world's glaciers (Paul et al., 2004; Bhambri et al., 2011; Mir et al., 2014a; Mir and Majeed, 2016). Among the total glaciers mapped in this study, 33 compound valley type glaciers with well-defined, long, and clearly observable tongues, covering 85% areal extent were selected for the length change studies. Length of the selected glaciers varied from 1.35 ± 0.06 km (G-61) to 17.8 ± 0.8 km (G-41) with a mean length of 4.4 ± 0.01 km. The length change during the assessment period varied from 0.87 ± 0.06 km at the rate of 17.2 ± 1 m a-1 (G-34) to 0.60 ± 0.04 km at the rate of 24.8 ± 0.2 m a-1 (G-49). Overall, from 1976 to 2011, the mean length/snout retreat was found to be 615 ± 20 m at a rate of 19 ± 0.5 m a-1. The results are similar to previous reports, which suggest that glacier snouts are retreating between 5 and 20 m a-1 (e.g., Kulkarni and Bahuguna, 2002; Kulkarni et al., 2007; Kumar et al., 2007; Wagnon et al., 2007; Dobhal et al., 2008). However, the recession of the glaciers is found to be heterogeneous (Fig. 6, part a). This heterogeneous pattern of retreat probably may be because of complex topographic configuration, greater thickness, and the effects of thick supra-glacial debris on large glaciers. Moreover, the length change depends on various factors such as bed slope, length of the glacier, response time, and other factors. It is commonly seen that the flat-bedded glaciers show large area change as a result of vertical thinning with marginal frontal recession. Conversely, steep slope at the glacier front will lead to excessive frontal recession. It is also noticed that the greater width of the frontal part of the glacier terminus yields a higher areal loss with lower terminus recession (Chand and Sharma, 2015). In this basin, no apparent correlation exists between the glaciers length and length change (Fig.6, part a). In addition, a weak correlation (R2=0.24) was found between the annual retreat rate and glaciers size. Further, the rate of glacier area change also differs from that rate of length change (Fig. 6, part b).
Several studies have reported that glacial debris cover has increased over time, indicating an increased debris production (Bolch et al., 2008; Bhambri et al., 2011; Kamp et al., 2011) and glacier thinning. The debris-covered glaciers show lower recession rate as compared to clean, ice-covered glaciers in the Himalaya (Iwata et al., 1980; Bolch et al., 2008; Bhambri et al., 2011; Kamp et al., 2011). In this study, out of the total 97 glaciers, 36 glaciers were identified to have heavily debris-covered ablation zones. Most of the debris-covered glacier parts are situated in the southern section of the basin (Fig. 1). The total area of debris-covered glaciers for the studied years is given in Table 3. The debris-covered zones have a total area of 36.5 ± 0.7 km2, which accounts for 16.1 ± 0.4% of the total glacierized area during 1976, whereas, during 1992 and 2000, the debris cover has an area of 45.6 ± 0.9 km2 (21.9 ± 0.5%) and 50.5 ± 1.0 km2 (25.5 ± 0.5%), respectively. Similarly, during 2006 and 2011, the debris cover area has further increased to 53.7 ± 1.1 km2 (28.0 ± 0.6%) and 60 ± 1.2 km2 (32.2 ± 0.6%), respectively. The debris-cover area in the catchment has increased by 23.5 ± 1.4 km2 (16.3 ± 3.8%) at a rate of 0.67 ± 0.04 km2 from 1976 to 2011. The 16.3 ± 3.8% increase in debris cover in this basin is similar to 17.8 ± 3.1% in the Saraswati/Alaknanda basin and slightly higher than 11.8 ± 3.0% in the upper Bhagirathi basin reported from 1968 to 2006 in Garhwal western Himalaya (Bhambri et al., 2009). In this basin, the higher number of completely debris-free/clean glaciers as found was spread on a smaller area than the lower number of large debris-covered glaciers (Table 3).
Results of trend analysis of meteorological parameters based on the Mann-Kendall and Sen's slope estimator. The possibility of significant trend at a 95% confidence level is represented by Mann-Kendall (Zs) statistics. The slope of the trend is given by Sen's slope (Qi) and probability of the trend is represented by the p-value. In this table the values in bold represent the significant rising or decreasing trends (“—”sign reflects the decreasing trends).
Trend Analysis of Climatic Variables
Temperature (Tmax and Tmin)
The estimation of trends indicated that the mean annual temperature (Tmax and Tmin) for Rakcham station has increased from 1985 to 2008 with the value of Zs being more than 1.96, reflecting that the probability of positive trends is very high (Table 4). The trends for temperature (Tmax and Tmin) were identified as statistically significant at the 95% confidence level. The annual mean Tmax showed an increasing trend, at a rate of 0.071 °C a-1. The Tmax has also increased with significant rise for pre-monsoon period. However, during winter and monsoon seasons, the Tmax has declined insignificantly. The mean annual Tmin clearly showed an increasing trend at a rate of 0.076 °C a-1 (Table 4). The Tmin has increased significantly for all periods with an insignificant rise during the post-monsoon season. For instance, during the past two decades, the Tmin has increased highly at rate of 0.14 °C a-1 for the winter season followed by a rising rate of 0.071 °C a-1 during pre-monsoon season and 0.022 °C a-1 during the monsoon season. During pre-monsoon season, the Tmax has also increased significantly at a rate of 0.11 °C a-1. The mean value of annual Tmax and Tmin from June to August throughout the past two decades also showed an increasing trend, with a slope of 0.003 °C a-1 and 0.013 °C a-1, respectively. The most significant evidence for regional and global climate change is the increase in minimum temperature (Tmin). Previous observations over a large portion of the Earth's land area also suggested that the minimum temperatures have increased at a quicker rate, about three times as much as the corresponding maximum one (Karl et al., 1991). Furthermore, an analysis of temperature data for the accumulation (Oct–Mar) and ablation (Apr–Sep) seasons indicated a rise in Tmin at a rate of 0.013 °C a-1 and Tmax at a rate of 0.003 °C a-1 during the accumulation season, whereas during the ablation period, a rise in Tmin at a rate of 0.02 °C a-1 and Tmax at a rate of 0.009 °C a-1 was observed. The linear trends in temperature for different time periods are shown in Figure 7, parts a and b.
The changes in the precipitation amount and type may directly influence the process of ablation of the glacier. For instance, the percolating rainfall can affect the rain heat flux and thus modify the ice-surface melt rate (Röhl, 2008). The increase in the rainfall might also enhance the glacier melting through the process of advection of the ambient heat and heat from the warm debris to the ice as well as reduced snow accumulation. Warming can also lead to rainfall occurrence in the higher altitudes of the glacier during the monsoon period where snowfall used to occur previously and accelerate the glacier melt. In this study, the rainfall data are available from only one weather station (i.e., Rakcham station) located at 3154 m a.s.l. in the valley and the rainfall at this station does not always mean a rainfall on the glacier surface. The trend analyses of rainfall at this station from 1985 to 2008 revealed an increasing trend during the annual, pre-monsoon, monsoon, and post-monsoon seasons. The trends were significant for the annual and monsoon seasons, which indicated a rise of 0.018 mm a-1 and 0.05 mm a-1, respectively. Consequently, the number of rainfall days has also increased for the studied period. However, during accumulation season, the Zs values of increasing rainfall were found to be less than 1.96, and hence, demonstrated the presence of statistically insignificant trends at 95% confidence level. The linear trends in rainfall for different time periods are shown in Figure 7, part c, and the significance of these trends and Sen's slope are given in Table 4.
The changes in snowfall have direct influence on the process of accumulation of the glacier. For instance, a decline in the snowfall may expose the glacier ice early and reduce the net accumulation by the year end. In the Baspa basin, the snowfall data is available from three meteorological stations (i.e., Sangla, Rakcham, and Chitkul). Overall, the trend analysis of the snowfall data indicated a decreasing pattern during the studied periods. For instance, the snowfall from 1976 to 2008 at Sangla station (2742 m a.s.l.) showed a reduction during all the studied seasons along with a decline in number of snowfall days. The decline was significant for annual and winter seasons, during which the snowfall has deceased at a rate of -6.28 mm a-1 and -0.11 mm a-1, respectively. Moreover, the value of Zs was found significantly higher than a standard value of 1.96, meaning the presence of a statistically significant negative trend at 95% confidence level around lower altitudes. Similarly, the snowfall from 1985 to 2008 at Rakcham station (3154 m a.s.l.) showed a reduction during annual, pre-monsoon, accumulation, and ablation seasons. During accumulation season, the decline was statistically significant at 95% confidence level and has decreased at a rate of -0.03 mm a-1. Contrary to the decline around low-altitude stations, the snowfall at Chitkul station, which is located at a relatively higher altitude (3444 m a.s.l.), showed an insignificant rise for annual, winter, and pre-monsoon seasons except for a decreasing trend during the post-monsoon season. At this station, the number of snowfall days has also increased during the period from 1976 to 2008. Overall, this nature indicated that the temperature in this basin still remains below freezing during the cold winter months, especially at high-altitude areas despite the increasing pattern around low-altitude areas. The linear trends in snowfall for different time seasons are shown in Figure 7, parts d, e, and f, and the significance of these trends with Sen's slope is given in Table 4.
Comparison with Other Studies
The present study generated five glacier inventories based on various temporal satellite datasets of Landsat series. For 1976, the inventory consists of 115 glaciers followed by 113 glaciers for 1992, 111 glaciers for 2000, and 109 glaciers for 2006. Last in the series, the 2011 glacier inventory also had 109 glaciers with a total area of 187.0 ± 3.7 km2 for the entire Baspa basin based on the Landsat-TM image. Previously, the glacier inventories for the basin have been generated by the GSI (Raina and Srivastava, 2008) and ESA (Paul and Andreassen, 2009). The ESA has generated the glacier outlines of this basin under the GlobGlacier project, which were later utilized in the compilation of RGI v 5.0 (Pfeffer et al., 2014). But, it is important to note that the main purpose of the GlobGlacier project was to derive the glacier extents for the whole northwestern Himalayan region rather than the Baspa river basin.
While comparing the present inventories with earlier inventories, a noticeable variation in the delineation of glaciers and reporting of the total number and area of glaciers have been observed. For the present Baspa basin, the GSI glacier inventory (2008) has reported about 89 glaciers in it. During an over-analysis of these glaciers upon the glaciers of SOI maps (1962), it was clearly confirmed that these glaciers have been derived from SOI maps for this basin. The SOI maps of 1962 are close to the Landsat MSS scene of 1976, and hence the GSI glacier inventory was compared with the glaciers derived from the MSS scene in this study. During analysis, a significant variation in the mapping and number of glaciers (22% or 26 glaciers) has been observed. About 19 glaciers as reported by the GSI inventory did not exist on the Landsat MSS and other data images. Additionally, 9 glaciers reported as individual glacier units by the GSI inventory were found disjointed/separated into more than 9 glaciers. Importantly, there are about 17 glaciers in the basin that are not reported by the GSI glacier inventory. For example, there are about 10 glaciers in western parts of the basin that are completely missed by the GSI inventory. There is also a significant variation between the glacier areas of the GSI inventory (∼235.6 km2) and the MSS-derived glacier area (227.4 km2) of the present study. The difference between the GSI glacier inventory as compared to current glacier inventories is probably because of the scale limitations of topographic maps of 1:50,000 scale (Raina and Srivastava, 2008). Moreover, the aerial photographs/satellite imageries of the months of November and January frequently used for the generation of SOI maps might have resulted in the inaccurate glacier identification because of seasonal snow and debris cover (Bhambri and Bolch, 2009). Furthermore, the smallest size fixed for identifying glaciers is found to be a major limiting factor for the exclusion of many glaciers by the GSI inventory. For instance, the smallest size of the glaciers included into the GSI inventory is 0.09 km2, which is larger than the size (0.01 km2) accounted for by the present inventory. Figure 8, part a, illustrates the variation between our analysis and the GSI inventory data (Raina and Srivastava, 2008).
Similarly, in comparison to the RGI v 5.0 databases derived from Landsat ETM+ 2001, there is a significant difference in the number of glaciers in the basin. The RGI v 5.0 inventory shows 140 glaciers as compared to the 109 glaciers (2011) of the present inventory. Further, a comparison of glacier areas indicated that the glacier area in the Baspa basin has been underestimated by RGI v 5.0 (∼159.7 km2) in comparison with the present study (187 km2). The difference of 27.3 km2 constitutes 14.6% of the total glacier area in the basin in 2011. This difference can be attributed to inadvertent mapping of mainly large debris-covered glaciers as fragmented, which arises out of automated mapping. There are about 18 glaciers shown as fragmented into more than 3 glaciers by the RGI v 5.0 databases. Figure 8, part b, illustrates the discrepancy between our analysis and RGI v 5.0 data. At a few places, the seasonal snow/ice patches and/or snow/ice aprons are also included in the GlobGlacier/RGI v 5.0 databases. Overall, the variation in the total number and area of glaciers in this basin can be attributed to the temporal differences in terms of acquired images and mapping period and differences in classification/definitions of glacier areas/boundary. Moreover, the small glacier size included in the RGI v 5.0 database was found to be 0.05 km2 as compared to the 0.01 km2 set for the present inventory, and the maximum area of the largest glacier (G-49) was observed to be 30.3 km2 compared to an area of 32.3 km2 in this study.
In this study, the glacier area loss of 41.2 ± 10.5 km2 (18.1 ± 4.1%) at a rate of 1.18 ± 0.3 km2 a-1 (0.51 ± 0.01%), and fragmentation as well as the disappearance of the glaciers was similar in comparison to other studies carried out previously in the Indian Himalaya. For example, based on 19 glaciers in the same Baspa basin, a loss of 19% glacier area has been reported from 1962 to 2001 (Kulkarni and Alex, 2003). Based on a study of 466 glaciers in the western Himalaya and archival records from 1962 to 2001, Kulkarni et al. (2005) has reported a loss of 21% of the glacierized area. Similarly, a number of tributary glaciers have detached from the trunk glaciers resulting in an increase in their number with reduced extent in the western Himalaya (Kulkarni et al., 2007; Chand and Sharma, 2015). Considering a smaller/insignificant loss from 1962 to 1976, the present study suggests a glacier area loss of 12.8 ± 4.4% (29.02 ± 10.2 km2) from 1976 to 2000 in comparison to a higher loss of 19% between the years 1962 and 2001 as reported by Kulkarni and Alex (2003). More importantly, the present glacier area loss of 18.1 ± 4.1% from 1976 to 2011 is lower than a loss of 24% (41.1 ± 8.4 km2) reported for 19 glaciers from 1962 to 2014 recently in the basin (Gaddam et al., 2016). The noticeably higher loss (e.g., 24%) is because of the use of SOI topographic maps as base maps for change detection studies in the above-discussed works. There are many studies which have reported problems of glacial extent demarcation in the SOI topographic maps (e.g., Bhambri and Bolch, 2009; Chand and Sharma, 2015). In this study, also, we have observed that besides mapping error of glacier boundaries, the SOI maps have ignored many glaciers in the basin as discussed above.
Impact of Glacier Size, Topography, and Debris on Glacier Changes
The glacier size, topographic features (e.g., elevation, slope, and aspect) and debris cover vary from glacier to glacier and therefore may affect their long-term response to climate as well. However, to understand the glacier area changes within different elevation, slope, and aspect ranges in the Baspa basin, the glacier's areas were estimated by counting the number of glacier DEM cells in different altitude, aspect, and slope ranges (Paul and Andreassen, 2009; Guo et al., 2015).
Present study shows that, on an average, the small glaciers (<0.50 km2) lost more of their surface area (53.5 ± 0.4%) than large glaciers, that is, a loss of 8.1 ± 1.5% (Fig. 9, part a) and indicated a significant control by glacier size on glacier loss. The correlation between glacier area and absolute area change for all glaciers was highly significant, that is, R2 = 0.74 (Fig. 10, part a). There are a number of studies suggesting an inverse relationship between the glacier size and rates of glacier area loss (Mark and Seltzer, 2005; Chueca et al., 2007; Chand and Sharma, 2015). Statistical analysis between size and percentage area loss of individual glaciers also shows a lower percentage of area loss for bigger glaciers in general, with wider distribution among smaller glaciers (Fig. 10, part b), which is attributed to the larger number of small glaciers in the basin as well as diverse aspect and altitude ranges of these glaciers.
In addition to size, the glacier area losses are also dependent on the elevation/altitude of the glaciers (Fig. 9, part b). A decreasing trend of glacier loss with rise in altitude was observed and thereby indicated that the lower-altitude glaciers are prone to disappear or lose more area. For instance, the glacier loss varied from 93.4 ± 4.2% at the 4100–4300 m a.s.l. zone to 2.6 ± 1.8% at the 6300–6500 m a.s.l. zone from 1976 to 2011. Similarly, in another following low-altitude zone of 4300–4500 m a.s.l., the glacier area loss of 45.6 ± 4.3% was observed. Up to the altitude zone of 5100–5300 m a.s.l., the glacier area loss has been very large and significant (e.g., 16.3 ± 2.3 km2, that is, 16 ± 4.5% up to this zone). The higher loss of 93.4 ± 4.2% of glacier area below an altitude of 4300 m a.s.l. was attributed to the complete disappearance of many small glaciers around this altitude zone during different time periods from 1976 to 2011 (Fig. 5). As observed, about 10 glaciers have completely vanished in the basin during the past three decades (Table 3).
Similarly, the exposure of the glaciers can influence the observed shrinkage trends. In this study, the orientation distribution showed that the location of glaciers is also dependent on local topographical constraints (Andreassen et al., 2008). The effect of exposure was quantified by grouping the values of extent into aspects having exposures ranging from north through east, south, and west. In this basin, the glaciers have, in general, northerly preferred aspect varying between northwest and northeast. However, the shrinkage of glaciers was higher in south to eastern orientations than in other aspects (Fig. 9, part c). For instance, the northwest and southwest aspects of the glaciers lost an area of 11.0 ± 4.5% (5.9 ± 2.4 km2) and 11.1 ± 4.4% (2.1 ± 0.84 km2), respectively, which is almost two times lower than a decrease of 29.6 ± 4.4% (11.6 ± 1.7 km2) of area in southern aspects. In addition, the southeast and east glacier aspects also revealed relatively higher glacier shrinkages of 27.7 ± 4.4% and 23.7 ± 4.4%, respectively. Despite the large number, the northern orientations of the glaciers also revealed smaller area loss of 11.2 ± 4.3% (3.0 ± 1.2 km2) in the basin. Previously, it has been reported that the south-facing slopes are exposed for longer time to solar radiation and therefore receive higher radiation resulting in higher melting of ice than the northern slopes (Ahmad and Rais, 1999, Bhambri et al., 2011; Mir et al., 2014b; Mir and Majeed, 2016). Less radiation makes north-facing glaciers react more slowly to climate change. The higher shrinkage of south-facing glaciers was also attributed to their smaller size (Mir et al., 2014b).
The glacier's surface slopes play a crucial role in determining the glacier recession, that is, the steeper the glacier, the larger the glacier area losses. In this basin, a direct relationship between the glacier area changes and slope was observed and thereby indicated that steep glaciers are more prone to losing more area. The glaciers were classified into five slope classes as shown in Figure 9, part d. As observed, the gentle glacier slope class of 0°–15° revealed a lower area loss of 9.0 ± 4.5% (8.7 ± 4.3 km2) than a steep slope class of 60°–80°, which revealed a higher area loss of 59.9 ± 4.2% (2.2 ± 0.15 km2). Overall, an increasing trend of glacier loss with rise in slope was observed. In the basin, small size glaciers generally exhibit steeper slopes because of low altitudinal ranges than do large glaciers.
The thick debris cover is considered one of the main factors in controlling glacier dynamics. It is suggested that the presence of supra-glacial debris cover, with thicknesses exceeding a few centimeters, lowers the rate of melting (Østrem, 1975; Pratap et al., 2015b), as these glaciers are less sensitive to climatic change (Kulkarni and Bahuguna, 2002; Pratap et al., 2015b), which leads to less glacier loss (Bolch et al., 2008; Nainwal et al., 2008b; Racoviteanu et al., 2008a, 2008b; Bhambri et al., 2011; Scherler et al., 2011). Moreover, the thick debris cover lowers the steady-state Area Accumulation Ratio (AAR), and thus a larger ablation area is required to balance accumulation, which further influences the glacier hypsometry (Benn and Evans, 2010; Chand and Sharma, 2015). In this study, based on the studies of 36 debris-covered glaciers, the effect of debris cover on glacier shrinkage was assessed by classifying glaciers into debris-covered and debris-free/clean categories for each year, respectively. During the study period from 1976 to 2011, clean glaciers lost an area of 26.6 ± 4.4% (0.59 ± 0.10 km2 a-1), whereas debris-covered glaciers lost only 13.8 ± 4.5% (0.58 ± 0.19 km2 a-1) (Fig. 9, part e). The large-sized glaciers with steep accumulation areas may have resulted in the production of heavy debris cover on their glacier ablation zones, which thereby protected these glaciers from greater melting. In addition, an analysis of the slope of the surrounding hills was observed to be >30°, which might be responsible for the increased hill-slope erosion rates and greater flux of rocky debris to the glacier surfaces. In this study basin, an analysis of debris-covered glaciers indicated that the debris has increased progressively during different time periods. For instance, the debris has increased by 9.1 ± 1.8 km2(5.9 ± 4.1%) from 1976 to 1992, by 4.9 ± 1.4 km2(3.6 ± 2.5%) from 1992 to 2000, by 3.3 ± 1.4 km2(2.6 ± 1.8%) from 2000 to 2006, and by 6.2 ± 1.6 km2(4.2 ± 3.0%) from 2006 to 2011. Overall, an increase of 16.3 ± 3.8% has been found between the years 1976 and 2011.
Impact of Climatic Factors on Glacier Changes
There are a number of studies reporting a significant effect of temperature and precipitation changes on glacier dynamics (Yao et al., 2012; Wiltshire, 2014; Mir et al., 2014a, 2014b; Yang et al., 2014). Overall, significant variations in temperature have a vital impact on glacier loss as well as the other controlling factors of glacier behavior (e.g., precipitation type). The precipitation type is a crucial variable because solid precipitation is contributing to accumulation, and liquid precipitation is contributing to runoff and additionally enhancing ablation via heat transfer and other processes. In this study, we analyzed the changes in air temperature and different precipitation types and their impacts on glacier changes. As already discussed, the meteorological data for the study area were available from three stations (i.e., Sangla, ∼2742 m a.s.l.; Rakcham, ∼3154 m a.s.l.; and Chitkul, ∼3444 m a.s.l.). However, it is essential to note that rainfall at meteorological stations at lower elevations in valleys than glacier elevations did not always mean simultaneous rainfall on the glacier surface. But, as per the temperature lapse rate, if it is snowfall at the base station it will necessarily be snowfall at the higher altitudes including glacier surface, whereas the rainfall at valley bottom may convert into snowfall in the higher reaches. Therefore, the variations in the solid precipitation type, that is snowfall, may be considered significant for the glacier health in this basin.
In this study, the trend analysis indicated an increasing pattern in temperature, especially a significant rise in Tmin. During the past two decades from 1985 to 2008, the mean annual Tmin and Tmax has increased at a rate 0.076 °C a-1 and 0.0071 °C a-1, respectively. Contrary to other seasons, the winter Tmin rose at a higher rate of 0.14 °C a-1, thereby suggesting the occurrence of warmer winters. But, in comparison to Tmin, the winter Tmax showed an insignificant decline. However, during the pre-monsoon season, the Tmax observed a significant increase of 0.11 °C a-1. The increase in Tmin is also significant for the studied accumulation and ablation seasons as well (Fig. 11, parts a and b). Thus, the warming of Tmin is considered to be one of the dominant factors responsible for glacier changes in this basin. The continuous increase in temperature (Tmin) for the winter months is more significant in terms of glacier fluctuation as it may be one of the causative factors for simultaneous ablation during the winter season (Kulkarni et al., 2011). In general, the rise in temperature leads to a rise in energy available for ice and snowmelt, followed by decline in snow accumulation and lower albedo of the glacier surface (Fujita and Ageta, 2000; Wang et al., 2014), all of which results in more ice melts and thus larger negative mass balance on the glaciers. The larger negative mass balance leads to a reduction in glacier area and retreat of glacier termini (Yao et al., 2004).
The observed warming is likely to lead to a decline in proportion of solid precipitation and an increase in the liquid precipitation particularly around low altitude areas even during the winter season when the glaciers receive most of their nourishment. Overall, the changes in the precipitation may cause shifting in snowline toward higher areas and that may have a significant effect on glacier loss. The changes in the form of precipitation may further enhance and reduce the accumulation and increase ablation. In this basin, snowfall has declined for all the studied seasons at Rakcham and Sangla stations, located around low altitudes, but not for Chitkul station, which is located at a comparatively higher altitude. The decrease in snow precipitation is more significant and pronounced in winter months and could be related to the significantly increasing Tmin during winter season as observed in the basin. The winter snowfall on glaciers contributes not only to glacial accumulation but also influences glacial ablation, as prolonged snow cover shields glaciers from radiation influx in the subsequent summer (Thayyen et al., 2007). Contrary to this, rainfall has increased during all the studied periods except the accumulation period, thereby indicating the conversion of form of precipitation from solid into liquid because of the warming pattern (Fujita, 2008; Shekhar et al., 2010; Mir et al., 2015a, 2015c), especially around low-altitude zones. The number of snowfall days has also decreased and is associated with an increase in the number of rainfall days. The larger increase in air temperature at lower elevations results in a large upward shift of the 0 °C isotherm as well, thereby causing more precipitation to fall as rain than as snow in the glacier elevations. The reduction in snowfall along with a rise in rainfall because of rising temperatures results in reduced accumulation and accelerated ablation along with an upward shift in the 0 °C isotherm line, particularly during summer season. The low snow deposition may also reduce the shielding effect of the high snow albedo during summers, leading to an increase in radiation absorption and consequently a higher melting. Consequently, the low-altitude glaciers face higher rates of retreat (Racoviteanu et al., 2008b; Bolch et al., 2012).
In the present basin, from 1976 to 2011, the reduction in glacier area has been mostly caused by the disappearance of low-altitude glaciers (∼10 glaciers disappeared from 1976 to 2011) and higher rates of loss and recession of small (<0.5 km2), low-altitude, south-facing, and debris-free glaciers in association with a significant warming temperature trend, especially Tmin, and therefore reflects their marked sensitivity to changes in climatic factors. The small glaciers are generally more sensitive to climatic changes because of their shorter response time to climate change (Bahr et al., 1998; Ye et al., 2001) and are likely to face more melting than larger ones (Armstrong et al., 2009). Additionally, the formation and progressive expansion of a Proglacial Lake at the Baspa Bamak glacier terminus (i.e., G-49) (Fig. 12, parts a–d) and presence of supra-glacial lakes on the glaciers such as G-66, G-12,G-80, G-49, and others indicated a clear response of glaciers to climate change in the basin. In general, the presence of lakes on the glacier's surface may further accelerate the glacier melting by modifying the stress regime of the glacier ice in contact and hence enhance the glacier loss (Komori, 2008; Röhl, 2008; Sakai et al., 2009; Gardelle et al., 2011).
Changes occurring to the Himalayan glaciers have wider implications for the region. Baspa river is a major tributary of River Satluj with a number of hydroelectric power projects downstream and sustains the agricultural activities of the fertile lands of Punjab and Haryana. In this study, remote sensing and GIS techniques have been used to generate five glacier inventories to assess the evolution of the glaciers and their controlling factors in the Baspa basin.
Using the satellite data from the Landsat series, the glacial parameters have indicated continuous and accelerated changes from 1976 to 2011. The glacier area has reduced by 41.2 ± 10.5 km2(18.1 ± 4.1%) at a rate of 1.18 ± 0.3 km2a-1 from 1976 (227.4 km2) to 2011 (186.2 km2). Similarly, the glacier retreat studied among 33 selected glaciers varied from 3.3 ± 0.03%, that is, 0.87 ± 0.06 km at the rate of 17.2 ± 1 m a-1 to 30 ± 6.6%, that is, 0.60 ± 0.04 km at the rate of 24.8 ± 0.2 m a-1. Furthermore, an analysis of 36 debris-covered glaciers indicated that the debris cover has increased by 23.5 ± 1.4 km2(16.3 ± 3.8%) from 1976 to 2011.
The present investigations revealed that the glacier evolution is controlled by various factors such as glacier size, topography (altitude, aspect, and slope), and debris cover. Factors assisting higher glacier loss in the basin are found to be (1) debris-free surface, (2) small size, (3) low altitude,(4) south to east orientation, and (5) steep surface slope. In addition to the local factors, glacier evolution is significantly controlled by climatic factors, particularly air temperature and precipitation. For instance, the statistical trend estimation suggested a significant rise in temperature (mainly Tmin) at 95% confidence level during the 1985 to 2008 period. Based on Sen's slope (Qi), the mean Tmin and Tmax have increased at a rate of 0.076 °C a-1 and 0.0071 °C a-1. During the accumulation period, Tmin and Tmax have increased at a rate of 0.013 °C a-1 and 0.003 °C a-1, whereas during the ablation period, Tmin and Tmax have increased at a rate of 0.021 °C a-1 and 0.009 °C a-1, respectively.
Tandem with the temperature increase, the snowfall has shown a reducing trend. However, rainfall has shown an increasing trend at the valley bottom, which may be due to the conversion of solid precipitation to liquid precipitation in response to the warming of the region. It is suggested that increasing temperatures and the changes in the form of precipitation (from snow to rain) may be the main factors responsible for the glacier recession in the basin. Moreover, the higher glacier loss (53.0 ± 0.4%) of small glaciers (<0.5 km2) highlighted the marked sensitivity of small glaciers to climatic changes. However, depending on glacier size, topography, and characteristics of clean/debris-covered glaciers, the response of the glaciers in the basin is found to be heterogeneous. Thus, the warming pattern and the corresponding glacier shrinkage as observed in the basin are reducing the overall water storage, both in the form of glaciers as well as seasonal snow in the basin.
First, thanks are due to BBMB for providing hydrometeorological data, NRSC ISRO for the LISS III data set, and NASA for making the Landsat datasets freely available under the umbrella USGS web server. The presented work is a part of Ph.D. thesis of R. A. Mir. The authors are grateful for detailed reviews by anonymous reviewers, which greatly helped to improve the content and structure of the manuscript.