Datasets
Glacier Data
We analyzed the glacier changes within seven periods: 1990–1994, 1994–1999, 1999–2004, 2000–2004, and 2004–2009, 2009-2015 and 2015-2018. The main sources of data were Landsat Thematic mapper (TM), Enhanced thematic mapper (ETM), Enhanced thematic mapper Plus (ETM+), and Operational land imager (OLI), which are available from the U.S. Geological Survey (http://earthexplorer.usgs.gov/). Alaska Satellite Facility (ASF) Digital elevation model (DEM) having a spatial resolution of 12.5m is downloaded from https://earthdata.nasa.gov/eosdis/daacs/asf. All remote sensing tiles of ablation season (July to August) having minimum cloud cover (<10%) were acquired from 1990-2019 (Table 1). All remote sensing images were preprocessed (image enhancement, co-registration and image stacking, image mosaicking) and were projected into World Geodetic System 1984 (WGS84), Universal Transverse Mercator (UTM) zone 43 projection, using Environment for Visualizing Images (ENVI), Erdas Imagine and ArcMap.
Table 1 list of Landsat images used in the study
Sensor
|
Path/Row
|
Date acquired
|
Cloud cover
|
RMSE
|
Landsat 8 OLI
|
149/35
|
4/8/2018
|
2.57%
|
3.9m
|
150/35
|
14/08/2018
|
2.57%
|
Landsat 7 ETM+
|
149/35
|
3/7/2015
|
2.00%
|
4.1m
|
150/35
|
23/07/2014
|
2.00%
|
149/35
|
26/07/2009
|
5.00%
|
150/35
|
25/07/2009
|
5.00%
|
150/35
|
28/06/2005
|
1.00%
|
149/35
|
5/8/2004
|
1.00%
|
Landsat 5 TM
|
149/35
|
16/08/1999
|
8.00%
|
4.5m
|
150/35
|
20/08/1998
|
10.00%
|
150/35
|
8/7/1994
|
4.00%
|
149/35
|
17/07/1994
|
5.00%
|
149/35
|
7/8/1990
|
1.00%
|
TM
|
149/35
|
26/07/1989
|
Less than 10%
|
6.3m
|
Note: Root mean square error (RMSE) is obtained during the image registration process
Meteorological data
In the upper Indus basin, around 20 meteorological stations are operated by different data collection organizations such as the Pakistan meteorological department (PMD). As PMD operated their stations from 1960s, after that in 2005 EvK2-CNR operated three stations in high altitude. The third automatic meteorological stations known as Data Collection Platforms (DCPs) are operated by Water and Power Development Authority (WAPDA) these stations are providing Data since 1995. There are around three high altitudes meteorological stations in Hunza river basin with one valley-based station namely Karimabad installed Pakistan met department in 2006. These stations provide daily climate data (temperature maximum and minimum, precipitation, solar radiation etc.). For this study, we have collected daily data of Minimum temperature, maximum temperature and precipitation from three stations (Khunjerab, Ziarat and Naltar) which WAPDA operates between 1995 and 2018 (Table 2).
Table 2. Overview of meteorological data of Hunza river basin (HRB) from 1995-2018
Station
|
Longitude
oE
|
Latitude
oN
|
Elevation
m.a.s.l
|
Tmax
oC
|
Tmin
oC
|
Ts
oC
|
Tw
oC
|
P
mm/yr
|
Khunjerab
|
74.397
|
36.838
|
4730
|
-0.88
|
-10.97
|
6.43
|
-13.44
|
202.21
|
Ziarat
|
74.459
|
36.817
|
3669
|
6.56
|
-2.65
|
9.46
|
-4.06
|
228.59
|
Naltar
|
74.177
|
38.157
|
2810
|
10.01
|
1.403
|
13.78
|
-1.01
|
648.89
|
Note: Tmax-mean maximum temperature, Tmin-mean minimum temperature, Ts-mean summer temperature, Tw-mean winter temperature and P- total annual precipitation
Method
Variations in glacier coverage
The automated mapping of glaciers from satellite/remote sensing imagery is hindered by debris cover of the glacier (Racoviteanu et al., 2008). Though automated methods are intended to these debris-covered glaciers, the results are less reliable, requiring post manual corrections (Paul et al., 2004 and Bhambri et al., 2011). Therefore, manual digitization was performed to ensure the interpretation of glacier boundaries as accurate as possible and control the error of <2% (Bolch et al., 2008). We performed visual and manual interpretation of remote sensing data for the identification and delineation of glacier masses. For the delineation of glacier outlines, we used the Normalized Difference Snow Index (NDSI).
NDSI measures the relative magnitude of the reflectance difference between visible (green) and shortwave infrared (SWIR). Firstly, we created our glacier inventory then we compared this with Randolph glacier inventory (RGI) version 6.0. After extraction of glaciers, we calculated the slope and Aspect of each glacier using ASF DEM. Finally, area and length of each glacier are calculated using geometry toolbox of ArcMap. All of the work is conducted in a system research institute (ESRI) ArcGIS 10.8 software. Additionally, it is for preprocessing remote sensing data like co-registration, image enhancement, and image mosaicking ENVI software is used. We have also used Google Earth is also used for identification of delineation of glacier boundaries.
Aspect and slope of glaciers
Using the approaches used in Manley (2008) and Paul et al. (2009), each glacier aspect is calculated. Using statistical tool of ArcMap tool, the mean gradient is estimated for each glacier, and then used dominant aspect is used to name an aspect of each glacier. Firstly, for slope calculation, using median filtering removed potential outlier altitudes then extracted ASF DEM for each glacier, and the elevation range is calculated based on their outlines.
Climate Trend
Daily climate data taken from three stations of HRB (Khunjerab, Ziarat and Nalter) were processed. For this study 3 climate indices (Tmax, Tmin, precipitation) were selected to explore changes in the climate of HRB from 1995-2018 using the 'Rclimdex' package (Zhang and Yang 2004). The RClimDex program uses linear regression for trend calculations (Powell and Keim, 2015). Using the package as mentioned above of RStudio mean monthly and mean annual Temperature and Precipitation is calculated.
2.3 Errors and Accuracy Assessment
The error in glacier delineation comes mainly from the operator's experience who identifies and demarcates the glacier boundaries Xiang et al., (2014). For the current study, the uncertainty is estimated based RS uncertainty approach proposed by Ye et al. (2006) and Li et al. (2015).
the area uncertainty is calculated byes
Where, UA is Area uncertainty of glacier area, and Ust is linear uncertainty which is given in eq.3. Area uncertainty ranges between 0.004 and 0.01 km2.
While the linear uncertainty is described as
Where Ul is linear uncertainty of glacier terminus, is the original pixel resolution of the satellite image, and is a co-registration error. Linear uncertainty was 6.3m for TM (1990) is, 4.5m for ETM (1994), 4.1m for ETM+ and 3.9m for OLI.