The cryosphere is one of the main factorsthat affect changesin globalsea level.Alpine glaciers are the most sensitive to climate change and are used as indicator for climate change (Zhou et al., 2018). Glacier retreat is accelerating with global warming as the amounts of CO2 and other greenhouse gases increase. To obtain information regarding the distribution and changing trends of glaciers, glacier inventory programs have been launched by the National Snow and Ice Data Center (NSIDC), theWorld Data Center in Cambridge (WDC-C) for Glaciology, the National Oceanic and Atmospheric Administration (NOAA), the National Aeronautics and Space Administration (NASA) and other agencies (Li et al., 1999). There are over 36,000 glaciers in the western China, only a few of them are under long-term observations (Shi and Liu, 2000). Glacier No. 1 in Tianshan, the headwater of Urumqi River, is the only one with nearly 30 years of observation data (Liu et al., 1999). Chinese Glacier Inventory (CGI) was built over 20 years and completed in 1999. This inventory mainly reflects the status of glaciers during the first aerial mapping period (late 1950s to late 1970s) in China and is a comprehensive collection of aerial photographs, topographic maps and field survey data forindividual glaciers and forChina’s basic glacier data census (Liu et al., 2000). However, the CGI does not reflect present status of glacier changes.
With the development ofsatellite imaging technology, remote sensing images can be used to identify glaciers for real-time, dynamic and large-scale monitoring. Remote sensing technology can be used to monitor alpine glaciers that have no field survey data. Thus, it can be used to properly evaluate and predict glacier changes and to analyze the impacts of the observed glacier changes on water resource and river runoff in arid areas in western China. In addition, it provides basic information on the use of water resources and ecological environmental protection.Traditional glacier identification methods include the use of threshold statistics, the ratio method, the Normalized Difference Snow Index (NDSI)method, and so on (Yan andWang, 2013). Shangguan et al. (2004) used the NDSI to extract glacier information of headwater region of the Yulongkashi River for examining the glacial changes. However, those methods are unable to distinguish debris-covered glacier, which has an area of 26,000 km2 at a global scale (Scherler et al., 2018). The presence of debris can not only reduce glacier ablation, but also strongly altersthe surface boundary condition and thus heat exchanges with the atmosphere (Collier et al., 2015)
Generally, the terminus of a glacier in the western China is covered by a layer of debris (Shi and Liu, 2000). The spectral characteristics of glaciers depend on their material, glacial moraine and snow and ice compositions. The proportions of glacial moraines, especially surface debris, and snow depth on glaciers affect their spectral features. Methods using optical remote sensing (e.g., NDSI method) can only identify relatively clean glacier but not debris-covered glacier (Shukla et al., 2010). The spectral information for this debris is similar to that of the bare rocks around the glacier.Thus, these debris-covered glaciers are difficult to distinguish, and this problemaffects glacier boundary extraction. However, the surface temperature of debris islower than that of the surrounding objects due to the underlying glacier (Bhambri et al., 2011; Karimi et al., 2012). Therefore, the surface temperature can be used to identify debris-covered glacier (Paul et al., 2004; Lambrecht et al., 2011). Remote sensing images with improved spatial resolution, and surface temperature extracted from thermal infrared image in the support of other auxiliary data such asthe digital elevation model (DEM) can be used to more accurately extract debriscovered glacier boundaries (Bhambri et al., 2011). Brenning et al. (2012) used day and night ASTER images to measure the near-surface temperature to extract the boundaries of the debris-covered glaciers, and the reuslt showed that the surface temperature of the Chilean Andes debris-covered glacier was significantly lower than that of the surrounding rocks and vegetation. This finding indicated that the method could effectively identify debris-covered glaciers when combined with optical and thermal infrared images. Shukla et al.(2010) used IRS-P6-AWiFS andTERRAASTER optical and thermal infrared images to map the boundaries of the Samudra Tapub glacier, and the result indicated that this method was suitable and very accurate for the identification of debris-covered glaciers.
Therefore, the free optical and thermal infrared images of Landsat are used to extract the boundaries of debris-covered glaciers in the Yigong Zangbo basin, and our study aims to improve the accuracy of glacier identification and to analyze dynamic glacier changes between 1990 and 2015. The impacts of the glacier change on moraine lakes are examined, and one of impact factors of glacier change, air temperature, is explored as well.
2. Study area
The study area is a sub-watershed of the Yigong Zangbo basin (94°32′ - 95°12′E, 30°15′ - 30°38′N) in the eastern region of the Nyainqentanglha mountains and at the eastern end of the Himalayas. This area is under the administration of Bomi County, Tibet Autonomous Region (Fig. 1). The area is 1,662 km2 and the elevation of the terrain is high in the north and low in the south with an altitude of 2,200 to 6,500 m. The terrain is characterized by high mountains and deep valleys and is dotted with numerous glaciers and snow-capped mountains. Yigong Zangbo River is a tributary of the Brahmaputra River.There isthe largest group of glaciersin theTibetan Plateau. Qiaqing glacier, one of China’sthree major glaciers and China’slargest temperate-glacier, is located in the Yigong Town of Bomi County. Southwest monsoon from the Indian Ocean brings plenty of precipitation, consequently, temperate-glaciers form in this area (Shi and Liu, 2000). The moisture-laden warm air current from the Indian Ocean enters the Yigong Zangbo valley along the Brahmaputra River and results in an apparent altitudinal climatic zonation from the bottom of the valley to the ridge of the watershed. The climate is subtropical under 2700 m, warm and semi-humid highland at 2700-4200 m, and cold and wet temperate highland above 4200 m. The records of Bomi station, the closest meteorological station to the study area, show that annual mean air temperature in the Yigong Zangbo basin is 8.5°C, with a mean air temperature of 17.2°Cin July and -1.7°Cin January, and annual mean precipitation is 958 mm.
Fig. 1. The location of Yigong Zangbo basin, overlying SRTM DEM
As a sub-watershed, the study area accounts for 12.3% of theYigong Zangbo basin (Fig. 1).According to records from the Global Land Ice Measurements from Space (GLIMS) (https://www.glims.org/), there are 207 glaciers in the study area, with most glaciers situated at an altitude of more than 3500 m. On the one hand, seasonal changes of temperature and precipitation in this basin are significantso that glacier melt is very evident and the equilibriumline islow (Yao et al., 2010). On the other hand, glaciersin the Yigong Zangbo basin are affected by warm and humid airflow from the Indian Ocean, these glaciers flow rapidly and have active glacial processes, and are extremely sensitive to global warming (Ke et al., 2013). Thus, it is necessary to examine glacier change in the Yigong Zangbo basin.
3. Data and methods
Generally, the most significant snowmelt occurs in summer.Therefore, it isthe least interference for glacier identification and remote sensing images selected from this period can improve the accuracy of glacier boundary extraction (Xu et al., 2013). However, due to cloudy weather and coarse temporal resolution, it is difficult to find suitable Landsat images in the same season, especially in summer. Here, Themetic Mapper (TM) images of Landsat5 (Path 135-Row 39) in 1990 and 2000, as well as Operational Land Imager (OLI) and Thermal Infrared Sensor (TIRS) of Landsat8 (Path 135-Row 39) in 2015 (Table 1), downloaded from the United States Geological Survey (USGS) (https://earthexplorer.usgs.gov/), are used for glacier identification and change analysisin the study area.TM images in 1990 and 2000 are in spring, but May is last month of spring, and glacier condition is close to that in summer, therefore, in our opinion, the images can be used for glacier change comparison with that of Landsat 8 in July 2015. Certainly, we have no other choice but to use them because of cloudy weather.
Table 1. Information for Landsat images used in this study
The Shuttle Radar Topography Mission (SRTM) has a resolution of 90 m, and this resolution can be improved to 30 m by smoothing with 11 * 11 Neighborhood Statistics and resampling.The processed DEM data are used to extract the boundaries of the study area, and to analyze changes in the glacier terminus elevation. The temperature data of Bomi meteorologicalstation (1961-2015), obtained from the National Meteorological Information Center of China, are used to analyze the impacts oftemperature changes on glaciers. The CGI and GLIMS database are used to obtain basic information of glaciers, for example, international uniformed name and code of glaciers, and also used as references to validate and correct the classification results.
2) Data preprocessing
All bands of TM and OLI/TIRS are calibrated, and atmospheric correction is performed with Fast Line-ofsight Atmospheric Analysis of Spectral Hypercubes (FLAASH), and topographic correction isimplemented with the improved C modification method (Huang et al., 2005). The brightness temperature obtained from the thermal bands is resampled to a resolution of 30 m in ordertomatch the optical band resolution.In addition, the brightnesstemperature data are normalized.All data are converted to the same projection (WGS 1984, UTM 47N), and cropped with vector boundary of the study area.
For debris-covered glaciers, the pre-processed images are classified using the spectral reflectance of optical bands and the brightnesstemperature ofthermal infrared bands(Fig. 2). Because optical characteristics and the thermal infrared radiation values ofilluminated areas differ from those of the shaded areas, the images are divided into illuminated and shaded areas according to D values calculated with formula (1) (Shukla et al., 2010).
\(D= \varphi -A\) (1)
Fig. 2. Glacier identification method based on optical and thermal infrared bands of Landsat TM and OLI/TIRS.
where φ is the sun azimuth at the time of image acquisition and A is the aspect angle derived from the DEM. When D is in the range of [-90°, 90º], covering area of the image is an illuminated area at that time, otherwise in a shaded area.The classified typesinclude clean glacier, debris-covered glacier, moraine lake and non-glacier area. The training samples are selected on the basis of the ground object characteristics of the optical band and the temperature characteristics of the thermal band. For the illuminated areas of the images, training samplesfor the optical band and thermal band are selected for each class,respectively. Forthe shaded areas, training samples are selected with the same method as the illuminated areas. Table 2 shows us the detail information of training samples for each class. The illuminated and shaded areas are classified with the maximum likelihood classification (MLC) method, respectively. In addition, the DEM data are used to facilitate classification and interpretation. Finally, the areas classified asilluminated orshaded are merged to complete the classification.
Table 2. The number of classification samples in the illuminated and shaded areas of the images
Due to distinguishing snow and glacier difficultly in the thermal infrared and optical bands if there is snow in the study area in May and July, a semi-automatic post correction is performed on the classification results by visual interpretation and referring to the CGI and GLIMS database (https://www.glims.org/). SRTM DEM was used to extract mountain ridges, then the boundary of each glacier was quantitatively defined. The imagesin 1990, 2000 and 2015 are classified using this method to obtain glacier boundaries and to distinguish clean glacier from debris-covered glacier, and our final result is close to the two glacier inventories. Next, statistical analyses of glacier area, moraine lake area, typical glacier length and glacier terminus elevation are carried out.
4) Accuracy assessment
In order to evaluate the accuracy of the classified results derived from remote sensing images covering the study area, we selected 20 same glaciers (about 10% of total glaciers) from the images in 1990, 2000 and 2015. For considering all kinds of conditions for glacier presence as much as possible and their representativeness of samples, these selected glaciers have different size, shape, average elevation, aspect, etc., moreover, some of them are debris-covered glaciers. Then, glacier experts with field experiences and remote sensing knowledge identified these glaciers on the images by visual interpretation, and this interpreted result can be regarded as the ‘true values’ for accuracy assessment of the classified results. The overall accuracies ofthe classified resultsin 1990, 2000 and 2015 were found to be 92.4%, 89.6% and 94.1%, respectively. On the one hand, snow is a possible impact factor for errors. It is difficult to distinguish snow and glacier because of similar spectral features. Particularly, images in 1990 and 2000 are in May, and there is possible snow cover at this time, their accuracies are therefore lower than that of images in July 2015. On the other hand, debris is also a possible impact factor for errors. If debris on a glacier is too thick, thermal remote sensing can not accurately identify the underlying glacier because of little or no surface temperature difference between debris and other surrounding objects (Brenning et al., 2012). Consequently, field validation is needed to reduce the error and improve the identification of debris-covered glacier.
The final classified results are shown in Fig. 3. Glacier areasin 1990, 2000 and 2015 were 928.76 km2 , 918.46 km2 and 901.51 km2 (Table 3),respectively.The total glacier area decreased by 27.25 km2 over 25 years, and this indicated that glacier retreat is very evident in theYigong Zangbo basin. However, the debris-covered glacier area slightly increased, and their areas in 1990, 2000 and 2015 were 63.39 km2 , 66.24 km2 and 71.16 km2, respectively. The debris-covered glaciers are mostly located at the glaciertongues and sides(Fig. 3), and debrisisresulted from freeze-thaw action of rocks under global warming, which can prevent glacier from melting to a certain extent (Vincent et al., 2016).
Fig. 3. Glacier and moraine lake distribution in the Yigong Zangbo basin in 1990, 2000 and 2015, 1-11 (figure a) denote index number of measurable glacier (Table 5), and A-I (figure c) denote index number of moraine lake (Table 6)
Table 3. Glacier areas in the Yigong Zangbo basin in 1990, 2000 and 2015 (unit: km2)
In an average condition, the rate of glacier retreat was 1.09 km2 yr-1 between 1990 and 2015 (Table 4). The glacier area decreased by 10.30 km2 with a retreat rate of 1.03 km2 yr-1 from 1990 to 2000, and it decreased by 16.95 km2 from 2000 to 2015 with a retreat rate of 1.13 km2 yr-1 . The retreat rate increased by 9.71% between 2000 and 2015 relative to that occurred between 1990 and 2000. The glacier shrinkage in the southeastern Tibetan Plateau (such as Yigong Zangbo basin) wasthe most pronounced in the entire Tibetan Plateau, and it also shows accelerating since 2000 (Azam et al., 2018; Yao et al., 2012).
Table 4. Area and rate of glacier retreat in the Yigong Zangbo basin from 1990 to 2015
Within 11 large glaciers examined (Table 5), the largest glacier 5O281B0729 (Gongpu) retreated from 186.53 km2 in 1990 to 182.03 km2 in 2000 and 181.17 km2 in 2015, its length shortened 665 m from 1990 to 2015, and its terminus elevation of ice tongue risen from 2777 m in 1990 to 2911 m in 2015. Meanwhile, the glacier 5O281B0768 (Nalong) showed the most significant area reduction of 3.37% from 1990 to 2000 and 2.56% from 2000 to 2015. All glaciers had the same trend, their areas reduced, resulting in shorter lengths and higher terminus elevations. This is consistent with the overall trends for glacier retreat in theTibetan Plateau (Brun et al., 2017). Liu et al.(2005) drew a similar conclusion and indicated that glaciersin the southeastern Tibetan Plateau retreated faster after 1980 due to global warming.
Table 5. Measurable glacier area, length and terminus elevation in the Yigong Zangbo basin from 1990 to 2015
The classified results ofimagesshowed that moraine lakes emerged at the glacier terminus and gradually increased in size (Fig. 3), and 9 moraine lakes exist in the study area (Table 6).The total area ofthese moraine lakes was 1.43 km2 in 1990 and increased by 38% (0.55 km2 ) in 2000, and it was 3.41 km2 in 2015, which was 2.38 times of its area in 1990. The area of the largest moraine lakeAincreased from 0.69 km2 in 1990 to 1.61 km2 in 2015. The area of moraine lake E in 2000 was 4.00 times of that in 1990, and expanded to 0.15 km2 in 2015, which is 15.00 times of that in 1990.
Table 6. Changes in area of moraine lakes in the Yigong Zangbo basin from1990 to 2015 (unit: km2)
All moraine lakes mainly expanded towards the glacier terminus, the glacier terminus retreated and its location was eventually occupied by the moraine lake (Fig. 3).The accelerated expansion ofthe moraine lake area after 2000 indirectly indicates the accelerated melting and retreat of the glaciers. Che et al. (2014) demonstrated that glacial lakes of the Pumqu River basin in the Tibetan Plateau expanded faster in the period of 2001-2013 than in the period of 1970s-2000. This is verified by our present study in the Yigong Zangbo basin. Due to global warming, the moraines formed depressions and blocked the rivers as the glaciers retreated, which resulted in a moraine lake at the glacier terminus and caused the area to expand gradually, and the increasing number and area ofmoraine lakes became more widespread in the Tibetan Plateau after 2000 (Veh et al., 2018). In addition, outburst floods of moraine lakes often occurred and threatened downstream residents(Harrison et al., 2018), therefore, monitoring moraine lake area expansion due to glacier retreat is important for preventing such disasters. Knowledge regarding the expansion of moraine lakes can provide a basis for making decisions to prevent moraine lake outburst and ensure the safety ofresidents living in the downstream.
The analysis of data from Bomi station indicates an overall increasing trend in air temperature of the Yigong Zangbo basin from 1961 to 2015 (Fig. 4), with a growth rate of approximately 0.325 °C/decade. The air temperature increased significantly faster between 1990 and 2015, with +0.603 °C/decade. The rapid air temperature rise in the recent decades enhanced freeze-thaw action, leading to more debris covering on glaciers and increase in area of debriscovered glaciers(Lambrecht et al., 2011). On the other hand, the rapid airtemperature rise could have been the main reason for the glacier retreat, which resulted in moraine lake expansion. The study by Shi and Liu (2000) also indicated that maritime glaciers had poor sensitivity to changes in precipitation, however, increasing temperaturessignificantly impacted changes in maritime glaciers and accelerated glacier melting.
Fig. 4. Change in annual mean air temperature (ºC) at Bomi meteorological station from 1961 to 2015.
The optical and thermal infrared band images of Landsat in 1990, 2000 and 2015 are used to identify and extract the glacier boundariesin theYigongZangbo basin. These measurements are based on the surface temperature differences among the clean glaciers, debris-covered glaciers and surrounding objects. Changes in the glacier area, length and terminus elevation are analyzed. In addition, changes in the moraine lake areas and the impacts of temperature rise on glaciers are further discussed.
Overall, the glaciers were receding and melting, with an area reduction of 27.25 km2 over 25 years. However, the debris-covered glacier area slightly increased, and their areas in 1990, 2000 and 2015 were 63.39 km2, 66.24 km2 and 71.16 km2 , respectively. In an average condition, the glaciersretreated at a rate of 1.03 km2 yr-1 between 1990 and 2000, at a rate of 1.13 km2 yr-1 between 2000 and 2015. The glaciers generally became shorter, and the terminus elevations of the glacier tongues risen. The total area of the moraine lakes was 1.43 km2 in 1990 and 1.98 km2 in 2000. The moraine lakes expanded drastically after 2000, reaching 3.41 km2 in 2015 (2.38 times ofits area in 1990).The annual mean air temperature risen significantly from 1961 to 2015, and the greater airtemperature increase accelerated the glacier retreat after 1990.
The data for glacier identification, especially for the parameters of the debris-covered glaciers on complex terrains, are obtained through the classification of optical and thermal bands of Landsat images. These data can be used as a supplement to the glacier inventory. Even in summer, the interference ofsnow in the high elevation areas could not be eliminated completely due to the similarspectral features ofsnow and glacier. In addition, thick debris is also an impact factor, and results in certain error occurrence in the glacier identification. Therefore, high-resolution SAR image should be used to overcome these drawbacks to obtain more accurate identification information of debris-covered glacierin the future study (Huang et al., 2018).
This research was supported by the National Nature Science Foundation of China (No. 41830105) and also funded by the International Scholar Exchange Fellowship(ISEF)programatKFAS(KoreanFoundation of Advanced Studies). The Landsat data are obtained from the United States Geological Survey (USGS) (https://earthexplorer.usgs.gov/).