Mapping debris-covered glaciers and identifying factors affecting the accuracy

Mapping debris-covered glaciers and identifying factors affecting the accuracy

Cold Regions Science and Technology 106–107 (2014) 161–174 Contents lists available at ScienceDirect Cold Regions Science and Technology journal hom...

4MB Sizes 0 Downloads 45 Views

Cold Regions Science and Technology 106–107 (2014) 161–174

Contents lists available at ScienceDirect

Cold Regions Science and Technology journal homepage: www.elsevier.com/locate/coldregions

Mapping debris-covered glaciers and identifying factors affecting the accuracy Anshuman Bhardwaj a,b,⁎, Pawan Kumar Joshi a, Snehmani b, Mritunjay Kumar Singh b,c, Lydia Sam d, R.D. Gupta c a

Department of Natural Resources, TERI University, 10 Institutional Area, Vasant Kunj, New Delhi 110 070, India Snow and Avalanche Study Establishment-Research and Development Center (SASE-RDC), Him Parisar, Plot No. 1, Sector 37A, Chandigarh 160036 (UT), India Motilal Nehru National Institute of Technology (MNNIT), Allahabad, 211004 UP, India d Defence Research and Development Organisation, New Delhi, India b c

a r t i c l e

i n f o

Article history: Received 16 July 2013 Accepted 19 July 2014 Available online 27 July 2014 Keywords: Supraglacial debris Morphometry Remote sensing Hamtah Glacier Patsio Glacier

a b s t r a c t Supraglacial debris significantly hampers the mapping of glaciers using remote sensing data. A semi-automated approach for the mapping of debris-covered glacier was applied, which combined the inputs from thermal and optical remote sensing data and the Digital Elevation Model (DEM) derived morphometric parameters. A thermal mask that delineates the supraglacial debris extent was generated by the thresholding of surface temperature layer obtained from Landsat TM/ETM+ thermal band satellite data. The extent of clean glacier ice was identified by band ratioing and thresholding of TM/ETM+ 4 and TM/ETM+ 5 bands. Morphometric parameters like slope, plan curvature and profile curvature were rearranged in similar surface groups using the technique of cluster analysis. All these masks were vectorized and final classification maps were generated using geographic information system (GIS) overlay operations. The areal extent of semi-automated outlines of Hamtah and Patsio Glaciers derived from cluster analysis varied from manually derived outline using pan-sharpened Landsat ETM+ September 2000 image by −1.3% and −1.6%, respectively. Year 2011 classification map for Patsio Glacier was compared with the field observations and a high correlation and overall accuracy (~91%) were observed. The same classification methodology was adopted for images of years 2000 and 1989 for Patsio Glacier to observe the effects of varying snow cover patterns on adopted methodology. Also the methodology was adopted and verified for Hamtah Glacier, with different geometry and terrain conditions as compared to Patsio Glacier. Although the spatial resolution limitation of ASTER GDEM and Landsat TM/ETM+ thermal band limits the automated mapping of small debris-covered glaciers, the outcomes are still favorable enough to apply such methodologies for mapping different types of debris-covered glaciers in the future. © 2014 Elsevier B.V. All rights reserved.

1. Introduction The Himalayan valley glaciers are often covered with the thick layer of debris toward the ablation zone and the snout (Shroder et al., 2000). The debris cover is considered as a part of an efficient sediment transport system in cold and high mountains (Kirkbride, 1995). This transport can be supraglacial, englacial or subglacial (Small, 1987). The debris externally originates from rockfalls, debris flows, rock avalanches and snow/ice avalanches (Snehmani et al., 2013) from adjacent slopes on the glacier ice (Fort, 2000; Hewitt, 2009; Shroder et al., 2000). The areal and volumetric assessments of supraglacial debris are important for studies related to glacier dynamics and mass and energy balance.

⁎ Corresponding author at: Snow and Avalanche Study Establishment-Research and Development Center (SASE-RDC), Him Parisar, Plot No. 1, Sector 37A, Chandigarh 160036 (UT), India. Tel.: +91 8146388210. E-mail addresses: [email protected] (A. Bhardwaj), [email protected] (M.K. Singh).

http://dx.doi.org/10.1016/j.coldregions.2014.07.006 0165-232X/© 2014 Elsevier B.V. All rights reserved.

Many studies have shown that thick debris cover on glaciers can greatly reduce their melting rate, thus causing their delayed response to climate change (Benn and Evans, 1998; Mattson, 2000; Pelto, 2000). Debris cover is also believed to promote formation of supraglacial lakes which can over a period lead to catastrophic floods (Benn and Warren, 2000; Komori, 2008; Reynolds, 2000). Because of the high altitudinal location of Himalayan glaciers and associated inclement weather conditions in most of the months of a year, field based monitoring of glaciers is risky, time-consuming and very costly. A coupled approach involving both, field based studies and remote sensing is therefore important to better study the debris-covered Himalayan glaciers (Vincent et al., 2013). Mapping of debris-covered glaciers based on remote-sensing techniques is very challenging as the use of visible (VIS) and near-infrared (NIR) bands is insufficient to detect the supraglacial debris cover (Bhambri et al., 2011a; Frey et al., 2012). But the recent progress in semi-automated techniques of monitoring large glacier coverage is promising (Bhambri and Bolch, 2009; Bhambri et al., 2011b; Bolch et al., 2010; Karimi et al., 2012).

162

A. Bhardwaj et al. / Cold Regions Science and Technology 106–107 (2014) 161–174

Fig. 1. Hypsometric distribution: a) Hamtah Glacier, b) Patsio Glacier.

For mapping clean glacier ice and snow, which is free of supraglacial debris cover, a number of methods are available which work well with multispectral satellite images. For example thresholding of the ratio images (Landsat Thematic Mapper (TM) band 4/TM band 5) is helpful to distinguish clean glacier ice (Sidjak and Wheate, 1999; Williams and Hall, 1998). But in case of the supraglacial debris-covered valley glaciers, the debris shows the same spectral characteristics as the surrounding terrain and is therefore very difficult to distinguish only by using optical images (Bishop et al., 2001; Raup et al., 2007). Manual heads-up digitization of debris-covered glacier is very time and effort taking for studying a large number of glaciers. Therefore in the past two decades, several semi-automated approaches for such classification have been proposed to overcome the challenges. Bishop et al. (1995, 1999) applied an ANN (Artificial Neural Network) based approach for approximating the supraglacial debris loads of Himalayan glaciers, using pre-defined glacier outlines. Paul et al. (2004) combined pixelbased multispectral classification with morphometric features, neighborhood relations (debris cover/clean ice) and change detection to attain a good accuracy. Taschner and Ranzi (2002) investigated the opportunities of Landsat TM and Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER) thermal band for debriscovered glacier monitoring. They concluded that slightly better results obtained using ASTER thermal band were because of its higher spatial resolution. Recently, Paul et al. (2013) have concluded that automated clean ice mapping is to be preferred over manual digitization with corrections for debris-covered and shadow areas. Rastner et al. (2013) compared pixel-based classifier with object-based classifier for glacier classification and reported a 12% higher accuracy in the case of objectbased classification of debris-covered glaciers. They attributed this better performance to post-processing possibilities like neighborhood analysis in the object-based classification. Several researchers (Bhambri et al., 2011a; Bolch et al., 2007; Shukla et al., 2010) have used morphometry based approach to map the debris-covered glacier with thermal inputs from ASTER data.

Racoviteanu et al. (2008) have reviewed the applicability of ASTER sensor for monitoring of Himalayan glaciers in detail. But ASTER data are not freely available for the Himalayas. There was a need to verify the approach on other Himalayan glaciers using freely available remote sensing data such as Landsat series. The knowledge gap still exists for developing a semi-automated classification methodology using Landsat bandwidths. The present study is an attempt to utilize freely available Landsat data for mapping two debris-covered glaciers in Western Himalayas: Patsio and Hamtah. Cloud-free Landsat images of distinct years and with varying seasonal snow cover were taken. Landsat images are apt for such studies and subsequent validation as they are freely available with considerably better resolution, availability of the Panchromatic band at 15 m spatial resolution (ETM +) and a thermal band (at 60 m spatial resolution for ETM+ and 120 m for TM). The conventional pixel value based classifiers of remote sensing images often do not perform appropriately in case of glacial terrain because of the heterogeneity and complexities involved. The objective of the study was to develop a semi-automated workflow using Landsat data for classifying different surface types (supraglacial debris, snow, ice) of a glacier. The existing methods in the literature for different kinds of data were refined accordingly and proper thresholds were found out. The outcome was properly validated using field sample points and the problem areas were discovered and discussed in detail. 2. Study area and data used The selection of the two glaciers in this study was not random and it was based on the fact that they represented different geometries, terrain conditions and extent of snow cover. Hamtah is a typical valley glacier with very simple geometry and side boundaries marked by prominent lateral moraines. Patsio, on the other hand, is a combined basin type glacier with more complex geometry and branched accumulation zones. The glaciers have different hypsometric distribution along with north-facing slopes (Fig. 1). The elevation layers for these glaciers

A. Bhardwaj et al. / Cold Regions Science and Technology 106–107 (2014) 161–174

163

Fig. 2. Study area: September, 2000 Landsat ETM+ FCC (False Color Composite) pan-sharpened image (4–3–2 combination).

were divided into 8 classes of equal range (100 m for Hamtah and 150 m for Patsio) for hypsometric analysis. The altitudinal setting of both the glaciers is widely non-overlapping. The last elevation band

of Hamtah constitutes the first elevation band of Patsio (Fig. 1). Almost all the elevation classes show significant differences in hypsometric profile. Thus, identifying these differences and selecting these two

Table 1 Details of satellite images used. S. no.

Data

Date of acquisition

Observed RMSE after co-registration (in pixel)

Hamtah Glacier 1. 2. 3.

Landsat TM (Path/Row: 147/38) Landsat ETM+ (Path/Row: 147/38) Landsat TM (Path/Row: 147/38)

09.10.1989 29.09.2000 03.10.2010

0.32 Master image 0.41

Patsio Glacier 1. 2. 3.

Landsat TM (Path/Row: 147/37) Landsat ETM+ (Path/Row: 147/37) Landsat TM (Path/Row: 147/37)

09.10.1989 13.09.2000 21.09.2011

0.45 Master image 0.42

164

A. Bhardwaj et al. / Cold Regions Science and Technology 106–107 (2014) 161–174

Fig. 3. Spatial distribution of GPS observations for Patsio Glacier (26 September 2011): Background image is Landsat TM (21 September 2011) FCC (4–3–2).

Fig. 4. Field photographs: a) Patsio Glacier, b) Stream coming out of snout, c) Glacier discharge meeting Bhaga River, d) Lateral moraine, e) Crevasse, f) Accumulation zone of the glacier, g) Ice mixed with debris near the snout (black circle represents mixed pixels in satellite image, red circle represents supraglacial debris), h) Supraglacial debris.

A. Bhardwaj et al. / Cold Regions Science and Technology 106–107 (2014) 161–174 Table 2 Details of Landsat bands used. S. no.

Band

Thematic Mapper (TM) 1. 4 2. 5 3. 6

Wavelength (in micrometers)

Spatial resolution (in meters)

0.76–0.90 1.55–1.75 10.40–12.50

30 30 120

Enhanced Thematic Mapper Plus (ETM+) 1. 4 0.77–0.90 2. 5 1.55–1.75 3. 6 10.40–12.50

30 30 60

glaciers were crucial steps in order to assess the effectiveness of the methodology and to prove that it is independent of hypsometric distribution of the glaciers. Hamtah is a north-facing glacier in the Chandra River basin of Lahaul and Spiti district of Himachal Pradesh, India. Hamtah is a simple basin glacier. According to WGMS (1989), a simple basin glacier is usually small sized with moderate thickness and has a single well defined accumulation zone. Detailed decade long glaciological studies by GSI (Geological Survey of India), Lucknow (India) since the year 2000, have been undertaken to know the various parameters of glaciers and the variety of climatic and hydrological inputs (GSI, 2007). This glacier is approximately 6 km long with an average observed retreat of ~ 16 m/year during the period 2000–2007 (GSI, 2007). This glacier falls in SoI (Survey of India) topographical sheets numbered 52H/7 and 52H/8. Patsio Glacier is located in the middle Himalayas (Sharma and Ganju, 2000) of Himachal Pradesh, India (Negi et al., 2013). Total length of this glacier is nearly 3.2 km as per measurement based on SoI topographical sheets. The glacier terminates approximately at an altitude of approx. 4609 m asl (Fig. 1). A high rate of retreat (22 m/year) and increasing Equilibrium Line Altitude (ELA) have been observed for this glacier during 1971–2011 (Negi et al., 2013). As observed in field visits, this glacier has mostly steep slope with a prominent bulge in the middle flow line. The discharge from this glacier contributes to the

165

Bhaga River. Fig. 1 shows the altitudinal ranges for both of these glaciers. Fig. 2 shows the geographical location and extent of the glaciers. Freely available satellite and elevation data were used for this study (source: earthexplorer.usgs.gov). The list of images with dates of acquisition is given in Table 1. The images were of mid-September to early October months of different years with varying seasonal snow cover. The acquisition of cloud and seasonal snow cover free images is a big issue in carrying out such study. Therefore, the images were chosen based on the presence of least cloud and snow cover over the glaciers for a particular year. A sufficient temporal gap of almost 10 years was kept in between the selected images to observe if the methodology was able to detect varying supraglacial debris cover. The presence of varying seasonal snow cover was similarly helpful to observe its effects on efficiency of the classification methodology. The Digital Elevation Model (DEM) used was ASTER GDEM version 2 product. Well distributed 23 point observations, made for Patsio Glacier on 26 September 2011 using handheld GPS, were used for accuracy analysis of the semiautomated classification image of 21 September 2011. Snow and Avalanche Study Establishment (SASE) has got a manned observatory near the snout of the glacier. Meteorological records from this observatory showed no evidence of snowfall or major deviation in temperature conditions between 21 and 26 September 2011, thus suggesting unchanged field conditions on the glacier. The field observation points are shown in Fig. 3 while Fig. 4 shows some field photographs of Patsio Glacier. A limitation of such handheld GPS device is its accuracy as it is not differentially accurate and thus not very admissible for precise point location based studies. But in context of the present study, for monitoring surrounding land cover and subsequent validation of classification maps, these devices could be considered appropriate sources of ground control points (GCPs). Such GPS is reported to provide vertical accuracies in order of ± 15 m in mountainous terrain (Racoviteanu et al., 2007) with horizontal accuracies up to ± 3.9 m (Ackerman et al., 2001). Here we were more concerned about horizontal accuracies. Description of the Landsat bands used is given in Table 2. The bands used in the study were TM/ETM + bands 4 and 5 and thermal band (band 6).

Fig. 5. Methodology adopted for semi-automated classification.

166

A. Bhardwaj et al. / Cold Regions Science and Technology 106–107 (2014) 161–174

Fig. 6. Morphometric layers for a) Hamtah Glacier, b) Patsio Glacier.

Fig. 7. Results of cluster analysis and vectorized boundaries. a) Hamtah Glacier, b) Patsio Glacier.

A. Bhardwaj et al. / Cold Regions Science and Technology 106–107 (2014) 161–174

3. Methodology Concurrently, three workflows were applied which were further integrated using GIS overlay operations. A thermal mask for the supraglacial debris extent was generated by thresholding of surface temperature images obtained from Landsat TM/ETM + thermal band. The ratioing and thresholding of TM/ETM+ 4 and TM/ETM+ 5 bands provided extent of clean glacier ice. Morphometric parameters like slope, plan curvature and profile curvature were rearranged in similar surface groups using cluster analysis. Overlay analysis of generated extents for the various glacier surface types (supraglacial debris, snow, ice) gave us final classification maps. The manual boundary delineation and classification of panchromatic bands of ETM+ images of the year 2000 of both the glaciers were done for validation of the classification approach. The flowchart of the methodology followed is shown in Fig. 5. The methodology can be broadly classified into the following steps: 3.1. Georectification ETM + images of the year 2000 were taken as master images because of the fact that ETM + multispectral bands could be pan-

167

sharpened to 15 m spatial resolution (Zhang, 2002) which further increased the ease and accuracy of the georectification with respect to GCPs as well as with other used images and of manual glacier boundary delineation. ETM + images of the year 2000 were rectified (RMSE in pixels = 0.43) with respect to 15 GCPs selected out of 23 GPS points for Patsio Glacier for subsequent accuracy assessment and error matrix generation. These 15 GCPs were particularly selected because they were clearly visible in both images of the year 2000, irrespective of the presence of seasonal snow. The year 2000 ETM+ images were taken as master images for further image-to-image co-registration. All the satellite images of both the glaciers were co-registered with respect to master images with considerably low RMSE (Root Mean Square Error) (Table 1). ASTER GDEM V2 was also co-registered with the master images before performing morphometric analysis (RMSE in pixels = 0.43). Landsat ETM+ images were co-registered with respect to the ASTER GDEM using the projective transformation algorithm (Bhambri et al., 2011b). To assess positional accuracy, 15 common location points such as river junctions were first used to rectify ETM+ images. Since ASTER GDEM already fits over the Landsat scenes available on USGS website, the changed transformation coefficients were applied to ASTER GDEM in order to directly co-register it with the newly rectified ETM+ scenes.

Fig. 8. Masks for supraglacial debris cover superimposed on Landsat FCC (4–3–2) images: Red rectangles show the pixels which were misclassified because of seasonal snow and were corrected manually in final classification maps (Fig. 11).

168

A. Bhardwaj et al. / Cold Regions Science and Technology 106–107 (2014) 161–174

3.2. Geomorphometry and cluster analysis for delineating glacier boundary This step was important for extracting the boundaries for both the glaciers. Using ASTER DEM, geomorphometric layers of plan curvature, profile curvature and slopes were generated. For the generation of slope and curvature, maximum triangle slope method (Tarboton, 1997) and degree polynomial approach (Zevenbergen and Thorne, 1987) were applied, respectively. The formed morphometric layers for both the glaciers are shown in Fig. 6. Cluster analysis was used in the study as a tool of classifying a raster image by clustering classes of different morphometric layers, representing similar surface characteristics together, in one layer. A combination of two different clustering processes, Iterative Minimum Distance (Forgy, 1965) and Hill-Climbing (Rubin, 1967), was applied simultaneously to all the three morphometric parameters (slope, profile and plan curvatures) and 3 cluster classes were given as output. The generated glacier boundaries were easily distinguishable, and free of adjoining misclassified grids (Fig. 7). It was vectorized to generate boundary polygon (Fig. 7). The areal extent of semi-automated outlines of Hamtah and Patsio Glaciers derived from cluster analysis was also compared with the manually derived outline on pan-sharpened Landsat ETM + September 2000 image. 3.3. Thermal mask for supraglacial debris Thermal band (band 6) of TM/ETM+ proved to be useful for mapping supraglacial debris cover. Taschner and Ranzi (2002) had performed a comparative study on Landsat and ASTER thermal bands for

mapping supraglacial debris. Lower spatial and spectral resolution of Landsat TM in the thermal region limits its capability, but still in the present case, the thresholding applied proved to be good enough to give the final classification of the glacier images with good accuracy. For both ETM+ images, the thresholds came out to be between 3.6 °C to 6.1 °C and for the rest of the TM images, it came out to be between 3.9 and 5.6 °C. This variation may be attributed to the varying dates of image acquisition. For all the TM images, the dates were of late September (21 September 2011) and October. But one important observation was that in any case, all the supraglacial pixels got covered between the threshold range of 3.6 and 6.1 °C. The thermal masks generated are shown in Fig. 8. Few of the misclassified polygons derived from TM thermal bands were manually corrected in the final step of overlay analysis.

3.4. Mask for clean glacier ice/snow Previous studies have shown that thresholded ratio images (Landsat Thematic Mapper (TM) band 4/TM band 5) are helpful to distinguish clean glacier ice (Sidjak and Wheate, 1999; Williams and Hall, 1998). The clean glacier ice/snow mask was generated in the binary image format based on TM/ETM+ 4 to 5 ratio. A threshold value of clean glacier ice ≥5 was taken for generating the binary mask. It was further converted into a vector polygon map. This step was repeated for both the glaciers for all the three years. The generated vector masks are shown in Fig. 9. These masks were mainly representative of the accumulation zones of the glaciers which were free of any impurity or debris cover.

Fig. 9. Masks for clean glacier ice/snow superimposed on Landsat FCC (4–3–2) images: Green rectangles show the pixels which were misclassified because of shadow of glacier head wall and were corrected manually in final classification maps (Fig. 11).

A. Bhardwaj et al. / Cold Regions Science and Technology 106–107 (2014) 161–174

169

Fig. 10. Land cover classes for glaciers in different years.

3.5. GIS overlay operations This is very crucial step in removing the misclassifications of the above mentioned steps. Actually the masks shown in Figs. 8 and 9 are the refined outcome of overlay operations. First of all, separate intersections of glacier boundary vector obtained from cluster analysis were performed with supraglacial debris masks and clean glacier ice/snow masks. The outcomes gave the supraglacial debris cover and clean glacier ice/snow lying within glacier boundaries only. Now union operation was applied to glacier boundary and to both outcomes of the intersection operation, which provided the final classification maps. After union, very few polygons of conflict (results of misclassifications) were observed. They were manually merged into appropriate land cover classes.

snowfall just before the image acquisition, which covered the glacier valley in the lower reaches, thus almost doubling the area under clean snow category and reducing the debris-covered part. The mixed pixels of snow and debris also increased because of this snowfall. In case of Patsio Glacier, the area under clean ice/snow category remained almost

Table 3 Observed areal extents of different land cover classes in different years: decreased supraglacial debris in case of Hamtah in the year 2010 is because of seasonal snowfall just before the day of image acquisition. It is shifting the areal extent of debris cover equally to other two classes, i.e., mixed pixels and clean snow. The double fold increase in the areal extent of supraglacial debris in the year 2000 in case of Patsio Glacier is because of reported rockfall and landslide activities in lower reaches (just above the snout) of the glacier. This debris covered a part of area coming under mixed pixels (Fig. 4g). Hamtah Glacier Land cover

1989 (area in km2)

2000 (area in km2)

2010 (area in km2)

Clean ice/snow Supraglacial debris Ice/snow mixed with debris

0.75 2.13 1.08

0.90 2.14 0.92

1.84 0.71 1.41

Land cover

1989 (area in km2)

2000 (area in km2)

2011 (area in km2)

Clean ice/snow Supraglacial debris Ice/snow mixed with debris

1.80 0.44 1.40

1.72 0.84 1.08

1.91 0.48 1.25

4. Results and validation The final classification maps are shown in Fig. 10. Areal extents of different land cover classes are given in Table 3. Year 2011 classification map for Patsio Glacier showed a higher overall accuracy of about 91.3% when compared to field observations. The same classification methodology was adopted for years 2000 and 1989 for Patsio Glacier to observe the effectiveness of the methodology in case of varying snow cover. In case of Hamtah Glacier, the area of supraglacial debris cover remains unchanged from 1989 (2.13 km2) to 2000 (2.14 km2) but in 2010 it dips down to 0.71 km2. This might be primarily because of fresh

Patsio Glacier

170

A. Bhardwaj et al. / Cold Regions Science and Technology 106–107 (2014) 161–174

Fig. 11. The observation points falling in different land cover classification categories: Black circles show the misclassified points. The blue rectangles highlight the misclassified regions between supraglacial debris cover and mixed pixels.

unchanged. But supraglacial debris cover almost doubled up in 2000 as compared to other two years (1989 and 2011). This was because of the fact that the image of 2000 was acquired in mid-September (13 September) when this glacier was still to receive fresh snowfall of the coming winter season. Most of the glacier surface was snow free and exposed. Therefore the extent of supraglacial debris cover increased and it took out a large portion from mixed pixels, converting them to the supraglacial debris cover category. Accuracy assessment was done at two different levels of data processing. Initially it was performed to check the result of morphometry and cluster analysis. The extracted glacier boundaries, when compared to manually digitized boundaries using pan-sharpened September 2000 Landsat ETM + images and ASTER GDEM, gave very less differences in areal extents. The measured areas of Hamtah and Patsio Glaciers, derived from cluster analysis, varied from manually derived outlines by only −1.3% and −1.6% respectively. The automatic extraction of glacier boundaries was most crucial step as further overlay operations required it as input. Afterward, as the next part of accuracy assessment, the accuracy of classification maps was assessed using available ground truth points (Fig. 3). Field observations for Patsio

Glacier were used for accuracy analysis of the semi-automated classification map of 21 September 2011 (Fig. 10). The observation points falling in different glacier surface types are shown in Fig. 11 to better represent their distribution on the glacier. The misclassified points are shown within the black circle in Fig. 11. The inaccuracies in the classification techniques usually are prominent in the areas where the land cover changes. Therefore, most of the field sample points were selected at the transition zones of the land cover classes to observe the accuracy of classification technique. Maximum efforts were involved to make sure that no biases were present while collecting field sample points and they were randomly distributed. An error matrix was generated to get an overall accuracy for the classification output (Table 4). It was found to be ~91%. In addition, producer's accuracy and user's accuracy were also estimated. The lowest user's accuracy (80%) was observed for mixed pixels because the misclassified pixels of debris cover went into this class. Subsequently, the lowest producer's accuracy of 75% was observed for supraglacial debris with 100% accuracy for the other two classes. All the GPS points were found to be accurately classified in case of classes representing clean glacier ice and mixed pixels (debris mixed with ice/snow). As expected, two misclassified GCPs out of total

Table 4 Error matrix. Reference data

Map data

Land cover

Clean ice/snow

Supraglacial debris

Ice/snow mixed with debris

Total

User's accuracy

Clean ice/snow Supraglacial debris Ice/snow mixed with debris Total Producer's accuracy Overall accuracy

7 – – 7 100% 91.3%

– 6 2 8 75%

– – 8 8 100%

7 6 10 23

100% 100% 80%

A. Bhardwaj et al. / Cold Regions Science and Technology 106–107 (2014) 161–174

eight were observed in the case of supraglacial debris. This can be attributed to 120 m resolution of the Landsat TM thermal band. 5. Discussion 5.1. Error estimation An error matrix is very effective to represent accuracy as the accuracies of each category are simply described along with both the errors of inclusion (commission errors) and errors of exclusion (omission errors) in the classification (Congalton, 1991). Other normally used ways of accuracy assessment are normalization of error matrix (Feinberg, 1970) and use of Khat (Kappa) statistics (Cohen, 1960). In the present case, there were only 3 classes out of which two were giving 100% producer's accuracy. Also, we were keen to identify the problem areas in the applied methodology rather than going into the robustness of accuracy assessment. This purpose was solved by a simple error matrix which clearly denoted that classification of supraglacial debris cover was compromised up to an extent.

171

Gangotri Glacier. They initially performed the cluster analysis for layers of plan and profile curvatures into 10 clusters and the outcome they further combined with slope layer through another cluster analysis with 10 clusters as output. They manually selected 3 out of these 10 categories and reclassified to get a glacier boundary. But the result also covered part of the lateral moraine near the accumulation zone of another glacier which they removed using information extracted from the ASTER thermal band. In the present case, following the same methodology (Bhambri et al., 2011a) gave a similar kind of output with some misclassification which did not give actual glacier boundary. But refining this result using Landsat thermal band did not work out well because of poorer spectral as well as spatial resolution (in case of TM images of 1989, 2010 and 2011). Therefore a different approach was taken here for cluster analysis. A combination of two different clustering processes, Iterative Minimum Distance (Forgy, 1965) and Hill-Climbing (Rubin, 1967), was applied simultaneously to all the three morphometric parameters (slope, profile and plan curvatures) and 3 cluster classes were given as output instead of 10. It excluded the lateral moraine in the lower reaches of the glacier (shown in Fig. 4d) and the results are shown in Fig. 12.

5.2. Misclassification of supraglacial debris and effect of seasonal snow 5.4. Unavailability of temporal elevation data for this part of Himalayas The problem areas were the pixels in the transition zone of supraglacial debris cover and mixed pixels consisting of ice/seasonal snow along with debris (Fig. 4g). The blue rectangles highlight the misclassified regions between supraglacial debris cover and mixed pixels. Coarse spatial resolution of the single thermal band of Landsat TM can be a possible explanation for this. But this much of misclassification has been reported by earlier studies also even if they used ASTER 90 m resolution thermal band (Bhambri et al., 2011a). A significant limitation of thermal remote sensing data for detecting debris-covered parts of the glacier is that the cooling effect of the underlying ice is supposed to be detectable in case of thin debris or in the presence of enough exposed ice. In case of very thick debris, the temperature of the debris is more affected by the lithology. In the present case, both the glaciers were quite smaller in size with moderate adjoining slopes in ablation zones as compared to larger glaciers like Gangotri with steeper adjoining slopes. This is the reason that the two glaciers studied have thinner debris cover. Several studies (Anders et al., 2006) have shown that effect of precipitation (especially monsoonal) and hence weathering and thickening of debris cover is much more prominent in the glaciers of Garhwal Himalayas. Since the methodology has worked well for even huge glaciers like Gangotri with thicker debris cover (Bhambri et al., 2011a,b), it was adopted with some modifications for the concerned glaciers and outcome was further assessed by field observations. Table 3 shows observed areal extents of different land cover classes in different years. Decreased supraglacial debris in case of Hamtah in the year 2010 can be attributed to seasonal snowfall just before the day of image acquisition, thus shifting the areal extent of debris cover equally to the other two classes, i.e., mixed pixels and clean snow. Another interesting observation, as can be seen in Fig. 10 in year 2000 classification map of Patsio Glacier, is the presence of increased supraglacial debris cover near the snout. This double fold increase in the areal extent of supraglacial debris was because of reported rockfall and landslide activities in the lower reaches (just above the snout) of the glacier. This debris covered a part of the area coming under mixed pixels (Fig. 4g). Over the period of 10 years, this extra debris got cleared because of glacier movement and the surface became as it was prior to the year 2000 (Fig. 10, 1989 classification map of Patsio).

Another point of discussion about the study was that temporal variations in land cover were observed, assuming that glacier geometry and boundary conditions remained same for all the 3 years. This was because of unavailability of freely available time series elevation data for this part of the Himalayas. Therefore, this study should be taken as a methodology development for glacier classification and glacier boundary extraction. Updated elevation data are expected to provide similar accuracy in boundary extraction and further classification as the slope is the key morphometric parameter in such kind of studies (Bhambri et al., 2011a). 5.5. Effect of shadow It was observed that shadow areas affected the supraglacial debris cover mapping as well as pure ice/snow pixel mapping. Especially in case of the Hamtah Glacier (Fig. 13) the effect was very clearly seen

5.3. Semi-automated glacier boundary delineation The significance of geomorphometric parameters in the delineation of glacial boundary is because of the fact that they are best representative of sudden changes in terrain conditions. Bhambri et al. (2011a) applied iterative minimum distance statistical technique (Forgy, 1965) for

Fig. 12. The variation in glacier boundaries near the snout, derived for Patsio Glacier, using iterative minimum distance statistical technique (yellow) (Forgy, 1965) and a combination of Iterative Minimum Distance (Forgy, 1965) and Hill-Climbing (Rubin, 1967) (red): 1st technique did not give actual glacier boundary as it included the lateral moraine in lower reaches within glacier boundary. The combined approach removed this issue and very closely matched with manually derived boundary as well as field observations.

172

A. Bhardwaj et al. / Cold Regions Science and Technology 106–107 (2014) 161–174

Fig. 13. The effect of mountain shadow on the classification of clean ice can be clearly seen here in case of accumulation zone of Hamtah Glacier: The blue boundary encompasses all the pure pixels of ice/snow whereas the portion in red rectangle shows the left out area because of partial shadow of glacier headwall.

on a large area in the accumulation zone, which was under shadow in every scene of Landsat or other widely used sensors (Fig. 2). It was easier to identify or remove because there was no possibility of debris cover or mixed pixels in the accumulation zone. Such instances were also reported by Bhambri et al. (2011a) in case of the Gangotri Glacier. These errors can be manually identified and removed during the final overlay step. 5.6. Developing thermal mask Bhambri et al. (2011a) applied thresholding directly on DN (Digital Number) image of thermal band. They made thermal mask by taking threshold values between 92 and 115. But if we are willing to observe temporal changes, the thresholding based on DN cannot be generalized for images of all the years. Therefore, we first converted the DN into LST (Land Surface Temperature) values and then we applied further

thresholding based on visual interpretation and literature. Taschner and Ranzi (2002) modeled surface temperature using the energy balance model for bare rock and supraglacial debris, located at the same elevation (2000 m above sea level). The range of surface temperature for debris cover, given by them in September at around 10.45 am (time of pass for Landsat TM and ETM+ over present study area) varied from 4 to 6 °C. In the present study also, similar ranges of surface temperature for masking out debris cover were obtained while thresholding. Few of the misclassified polygons derived from the thresholding of TM thermal bands were observed. The coarse resolution of the TM thermal band was the reason behind the existence of such polygons. 6. Conclusions This study was an effort to establish the usefulness of free of cost Landsat series data for the purpose of a GIS based semi-automated

A. Bhardwaj et al. / Cold Regions Science and Technology 106–107 (2014) 161–174

glacier classification. The outcome was properly evaluated with respect to field sample points. 1. The variables such as seasonal snow cover, changing light conditions and satellite bandwidth variations, which could cause the errors in the classification were also taken into account by incorporating temporal remote sensing images. Since the sun elevation angle is always changing with respect to time and a particular location on the Earth, changing sun elevation leads to changing light conditions and thus shadows. 2. In order to observe the area–elevation control over the applied thresholds and the classification outcomes, two glaciers with different hypsometric profiles were selected for the study. The present study is different from earlier such studies in its approach of comparing the same methodology for two different glaciers quite far apart from each other and in different elevation zones with different geomorphology. A hypsometric comparison was thought out to be a proper way of prominently looking into these differences while mapping these glaciers. These differences were not found to have much controlling effects because same thresholds and methods properly mapped both the glaciers. 3. The spatial resolution limitation of ASTER DEM and low spectral (only 1 band) and spatial resolution of Landsat TM/ETM+ thermal bands were main bottlenecks in this study especially because the glaciers were having comparatively smaller areal extents. Still the outcomes are favorable enough to apply such methodologies for mapping debris-covered glaciers in the future. 4. Analysis of surface temperature profiles of supraglacial and periglacial debris/rocks, given by Taschner and Ranzi (2002) clearly revealed the considerable temperature differences existing between them, which was the basis for attempting such methodology using Landsat data. 5. The methods of morphometry, cluster analysis and overlay analysis applied in this study were different from earlier studies which were done on different glaciers with different lithology, climatic setting, large size and different kind of remote sensing data in terms of spectral and spatial resolutions. For instance, Bhambri et al. (2011a) used ASTER data for Gangotri Glacier which is many times bigger than Patsio or Hamtah. The approach they adopted for morphometric and cluster analyses to get good result for Gangotri, did not entirely apply here (as mentioned in Sections 3.2 and 5.2). 6. The replication of the suggested methodology for different types of glaciers involves a deep analysis of the data available and the adjustments in the thresholds, based on the problem areas discussed above. Keeping in mind the recent Landsat-8 data with availability of two thermal bands (for debris cover mapping), freely available data were used in the entire study to establish a specific methodology as proper thresholding and method selection are crucial parts of such studies. References Ackerman, T.,Erickson, T.,Williams, M.W., 2001. Combining GIS and GPS to improve our understanding of the spatial distribution of snow water equivalence (SWE). Proceedings of the 2001 ESRI User Conference, 10 July 2001, San Diego, CA. Anders, A.M., Roe, G.H., Hallet, B., Montgomery, D.R., Finnegan, N.J., Putkonen, J., 2006. Spatial patterns of precipitation and topography in the Himalaya. Spec. Pap. Geol. Soc. Am. 398, 39–53. Benn, D.I., Evans, D.J.A., 1998. Glaciers and Glaciation. Arnold, New York. Benn, D.I., Warren, C.R., 2000. Rapid growth of a supraglacial lake, Ngozumpa Glacier, Khumbu Himal, Nepal. Debris-Covered Glaciers 2, 177–185. Bhambri, R.,Bolch, T., 2009. Glacier mapping: a review with special reference to the Indian Himalayas. Prog. Phys. Geogr. 33, 672–704. Bhambri, R., Bolch, T., Chaujar, R.K., 2011a. Mapping of debris-covered glaciers in the Garhwal Himalayas using ASTER DEMs and thermal data. Int. J. Remote Sens. 32 (23), 8095–8119. Bhambri, R., Bolch, T., Chaujar, R.K., Kulshreshtha, S.C., 2011b. Glacier changes in the Garhwal Himalayas, India 1968–2006 based on remote sensing. J. Glaciol. 57 (203), 543–556.

173

Bishop, M.P., Shroder Jr., J.F., Ward, L.J., 1995. SPOT multispectral analysis for producing supraglacial debris-load estimates for Batura Glacier, Pakistan. Geocarto Int. 10, 81–90. Bishop, M.P., Shroder Jr., J.F., Hickman, B.L., 1999. SPOT panchromatic imagery and neural networks for information extraction in a complex mountain environment. Geocarto Int. 14, 19–28. Bishop, M.P., Bonk, R., Kamp, U., Shroder Jr., J.F., 2001. Terrain analysis and data modeling for alpine glacier mapping. Polar Geogr. 25, 182–201. Bolch, T., Buchroithner, M.F., Kunert, A., Kamp, U., 2007. Automated delineation of debriscovered glaciers based on ASTER data. Geoinformation in Europe (Proc. of 27th EARSel Symposium, 04–07 June 2007), Bozen, Italy, pp. 403–410. Bolch, T., Menounos, B., Wheate, R., 2010. Landsat-based glacier inventory of western Canada, 1985–2005. Remote Sens. Environ. 114, 127–137. Cohen, J., 1960. A coefficient of agreement for nominal scales. Educ. Psychol. Meas. 20, 37–46. Congalton, R.G., 1991. A review of assessing the accuracy of classifications of remotely sensed data. Remote Sens. Environ. 37 (1), 35–46. Feinberg, S.E., 1970. An iterative procedure for estimation in contingency tables. Ann. Math. Stat. 41, 907–917. Forgy, E., 1965. Cluster analysis of multivariate data: efficiency vs. interpretability of classifications. Biometrics 21, 768. Fort, M., 2000. Glaciers and mass wasting processes: their influence on the shaping of Kali Gandaki Valley, Nepal. Quat. Int. 65, 101–119. Frey, H., Paul, F., Strozzi, T., 2012. Compilation of a glacier inventory for the western Himalayas from satellite data: methods, challenges, and results. Remote Sens. Environ. 124, 832–843. GSI, 2007. Detailed glaciological studies on Hamtah Glacier, Lahaul and Spiti district. In: Himachal Pradesh: Rec. Geol. Surv. India v. 139, pt. 8, 136–139. Hewitt, K., 2009. Rock avalanches that travel onto glaciers and related developments, Karakoram Himalaya, Inner Asia. Geomorphology 103, 66–79. Karimi, N., Farokhnia, A., Karimi, L., Eftekhari, M., Ghalkhani, H., 2012. Combining optical and thermal remote sensing data for mapping debris-covered glaciers (Alamkouh Glaciers, Iran). Cold Reg. Sci. Technol. 71, 73–83. Kirkbride, M.P., 1995. Processes of transportation. In: Menzies, J. (Ed.), Modern Glacial Environments: Processes, Dynamics and Sediments. Glacial Environments, vol. 1. Butterworth-Heinemann, Oxford, pp. 261–292. Komori, J., 2008. Recent expansions of glacial lakes in the Bhutan Himalayas. Quat. Int. 184, 177–186. Mattson, L.E., 2000. The influence of a debris cover on the mid-summer discharge of Dom Glacier, Canadian Rocky Mountains. IAHS Publ. 264, 25–33. Negi, H.S.,Saravana, G.,Rout, R.,Snehmani, 2013. Monitoring of great Himalayan glaciers in Patsio region, India using remote sensing and climatic observations. Curr. Sci. 105 (10), 1383–1392. Paul, F.,Huggel, C.,Kääb, A., 2004. Combining satellite multispectral image data and a digital elevation model for mapping debris-covered glaciers. Remote Sens. Environ. 89 (4), 510–518. Paul, F.,Barrand, N.,Baumann, S.,Berthier, E.,Bolch, T.,Casey, K.,Frey, H., et al., 2013. On the accuracy of glacier outlines derived from remote sensing data. Ann. Glaciol. 54 (63), 171–182. Pelto, M.S., 2000. Mass balance of adjacent debris-covered and clean glacier ice in the North Cascades, Washington. IAHS Publ. 264, 35–42. Racoviteanu, A.E., Manley, W.F., Arnaud, Y.,Williams, M., 2007. Evaluating digital elevation models for glaciologic applications: an example from Nevado Coropuna, Peruvian Andes. Global Planet. Chang. 59, 110–125. Racoviteanu, A.E.,Williams, M.W.,Barry, R.G., 2008. Optical remote sensing of glacier characteristics: a review with focus on the Himalaya. Sensors 8 (5), 3355–3383. Rastner, P., Bolch, T., Notarnicola, C., Paul, F., 2013. A comparison of pixel- and objectbased glacier classification with optical satellite images. IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens. http://dx.doi.org/10.1109/JSTARS.2013.2274668. Raup, B.H., Kääb, A., Kargel, J.S., Bishop, M.P., Hamilton, G., Lee, E., Paul, F., Rau, F., Soltesz, D., Khalsa, S.J.S., Beedle, M., Helm, C., 2007. Remote sensing and GIS technology in the Global Land Ice Measurements from Space (GLIMS) project. Comput. Geosci. 33, 104–125. Reynolds, J.M., 2000. On the formation of supraglacial lakes on debris-covered glaciers. AHS Publ. 264, 153–161. Rubin, J., 1967. Optimal classification into groups: an approach for solving the taxonomy problem. J. Theor. Biol. 15 (1), 103–144. Sharma, S.S., Ganju, A., 2000. Complexities of avalanche forecasting in Western Himalaya — an overview. Cold Reg. Sci. Technol. 31, 95–102. Shroder Jr., J.F., Bishop, M.P., Sloan, V., Copland, L., 2000. Debris-covered glaciers and rock glaciers in the Nanga Parbat Himalayas, Pakistan. Geogr. Ann. 17–31. Shukla, A., Arora, M.K., Gupta, R.P., 2010. Synergistic approach for mapping debriscovered glaciers using optical and thermal remote sensing data with inputs from geomorphometric parameters. Remote Sens. Environ. 114 (7), 1378–1387. Sidjak, R.W., Wheate, R.D., 1999. Glacier mapping of the Illecillewaet Icefield, British Columbia, Canada, using Landsat TM and digital elevation data. Int. J. Remote Sens. 20 (2), 273–284. Small, R.J., 1987. Englacial and supraglacial sediment transport and deposition. In: Gurnell, A.M., Clark, M.J. (Eds.), Glaciofluvial Sediment Transfer: An Alpine Perspective. John Wiley & Sons, Chichester, pp. 111–145. Snehmani, Bhardwaj, A., Pandit, A., Ganju, A., 2013. Demarcation of potential avalanche sites using remote sensing and ground observations: a case study of Gangotri Glacier. Geocarto Int. http://dx.doi.org/10.1080/10106049.2013.807304. Tarboton, D.G., 1997. A new method for the determination of flow directions and upslope areas in grid digital elevation models. Water Resour. Res. 33 (2), 309–319.

174

A. Bhardwaj et al. / Cold Regions Science and Technology 106–107 (2014) 161–174

Taschner, S., Ranzi, R., 2002. Comparing the opportunities of Landsat-TM and Aster data for monitoring a debris-covered glacier in the Italian Alps within the GLIMS project. Geoscience and Remote Sensing Symposium, 2, pp. 1044–1046. Vincent, C., Ramanathan, Al, Wagnon, P., Dobhal, D.P., Linda, A., Berthier, E., Sharma, P., Arnaud, Y., Azam, M.F., Jose, P.G., Gardelle, J., 2013. Balanced conditions or slight mass gain of glaciers in the Lahaul and Spiti region (northern India, Himalaya) during the nineties preceded recent mass loss. Cryosphere 7 (2), 569–582. Williams Jr., R.S., Hall, D.K., 1998. Use of remote-sensing techniques. In: Haeberli, W., Hoelzle, M., Suter, S. (Eds.), Into the Second Century of Worldwide Glacier Monitoring: Prospects and Strategies. UNESCO Publishing, Paris, pp. 97–111.

World Glacier Monitoring Service (WGMS), 1989. WGMS (1989): World Glacier Inventory — Status 1988. In: Haeberli, W., Bösch, H., Scherler, K., Østrem, G., Wallén, C.C. (Eds.), IAHS (ICSI)/UNEP/UNESCO, World Glacier Monitoring Service, Zurich, Switzerland (58 pp.). Zevenbergen, L.W., Thorne, C.R., 1987. Quantitative analysis of land surface topography. Earth Surf. Process. Landf. 12 (1), 47–56. Zhang, Y., 2002. Problems in the fusion of commercial high-resolution satellite as well as Landsat 7 images and initial solutions. ISPRS, Vol. 34, Part 4, GeoSpatial Theory, Processing and Applications, Ottawa, Canada, vol. 34.