Simulation of glacial lake outburst flood hazard in Hunza valley of upper Indus Basin

The UIB (Upper Indus Basin) is prone to GLOFs (Glacial Lake Outburst Floods). Physical monitoring of such a large area on a regular basis is a challenging task, especially when the temporal and spatial extent of the hazard is highly variable. The purpose of this study was to map the potentially dangerous glacial lakes and simulate the associated hazard in the downward settlements using HEC-RAS in the GIS environment using Landsat 7 remote sensing data. The study was conducted in Hunza Valley of UIB, where there are several human settlements that are endangered due to the GLOF hazard. Sudden breaches in the unstable moraine dams adjoining receding glaciers may occur because of the rapid and huge accumulation of turbulent water in the glacial lakes. The ASTER GDEM (Digital Elevation Model) is utilized to detect flow accumulation of glacial hazard involving slope, elevation, and orientation of the mountain glaciers. The study results revealed that settlements of Hunza Valley are threatened by the GLOFs hazard. Keeping in view the seasonal growth of the potentially dangerous glacial lakes of Hunza Valley, a low discharge of 3500m 3 /s from a potentially dangerous glacial lake can affect 40%, whereas a moderate discharge of 5000m 3 /s can affect 60%, and a high discharge of 7000m 3 /s can affect 80% of the Shimshal village habitat. The results of the study can provide a platform for the establishment of an early warning and monitoring system to minimize the impact of future GLOFs. Accurate and comprehensive knowledge of potentially dangerous GLOFs is of utmost importance for risk management. A digital repository of GLOFs can enhance the ability to inform policymakers on the vulnerability, risk mitigation, and action/adaptation measures.


Introduction
*The alpine glaciers of the Hindukush-Karakoram-Himalaya region are a renewable natural freshwater storehouse that benefits hundreds of millions of people downstream (Immerzeel et al., 2010). Owing to global warming acceleration, the glaciers of the mid-latitude region of Pakistan are retreating since the second half of the 20 th century (Das and Meher, 2019). On the retreating glacier terminus, this phenomenon has accounted for the accumulation of many disastrous glacial lakes (Shah and Kanth, 2012). The damming by unstable moraines has caused several glacial lakes. The disastrous GLOFs containing debris and a large quantity of turbulent water lead to the sudden breaches of these unstable moraines, which hold huge quantities of water (Bhambri et al., 2013). Glaciation and interglaciation are natural processes that have occurred several times during the last 10,000 years. A situation that provides a large space for retaining meltwater and leads to the formation of morainedammed lakes. Glacier-connected lakes have likely accelerated the glacial retreat via thermal energy transmission and contributed to over 15% of the area loss in their connected glaciers (Che et al., 2014). On the other hand, significant glacial retreats led to disconnections from their proglacial lakes, which appeared to stabilize the lakes in the Himalayas. Continuous expansions in the lakes connected with debris-covered glaciers, therefore, need additional attention due to their potential outbursts (Nie et al., 2018). Glacier retreat is an indication of glacial lake formation. The glacial hazard of GLOFs can cause loss of life, livestock, property, valuable forests, costly mountain infrastructures, farmlands, and pasture resources. Damages to settlements and farmland can take place at very great distances from the outburst source. For example, in Pakistan, damage occurred 1,300 Km from the outburst source (Che et al., 2014). Much of the damage created during GLOF events are associated with the large amounts of debris that accompany the floodwaters (Budhathoki et al., 2010). In the past 20 years, glaciers in the Himalayas have retreated and thinned rapidly as a response to regional climate warming, leading to the formation of new glacial lakes and the expansion of existing glacial lakes (Kaushik et al., 2019). These areas are located in the border belt and the Eurasian plates, where tectonic seismic activity is frequent and intense. Earthquakes have often compromised the stability of mountain slopes, glaciers, and moraine dams, resulting in an imbalance in the state of glacial lakes (Wang and Zhou, 2017). Another emerging hypothesis of more GLOF events is the change in the pattern of rainfall (Khan et al., 2019;Harrison et al., 2018). A historic glacial flood burst had a depth of around 30m at the junction of Shimshal and Hunza, (about 40Km from the assumed position of the lake) and destroyed the village of Passu near the Hunza river (Goudie et al., 1984). Glacier-fed lakes are dominant in both quantity and area and exhibit an overall faster expansion trend compared to the nonglacier-fed lakes in the Himalayas (Mir et al., 2018). Glaciological characteristics of the ablation zone of the glaciers of the Shigar River basin have changed in the recent past because of a host of factors (Mayer et al., 2006). The diversity of glacial material is the prime reason for the diverse behavior and peculiar dynamics of the glaciers. Heterogeneity in the Karakoram glacier surges is observed because of the peculiar dynamics of the glaciers (Quincey et al., 2015). The glacier surges are propagated once coupled with glacial lakes. Glacier changes in the Karakoram region are mapped temporally in order to observe the diverse behavior of the glaciers (Rankl et al., 2014). On debris-covered glaciers, glacial lake formation is observed at a faster pace in the alpine region of Pakistan (Hambrey et al., 2008;Raup et al., 2007). A conceptual analysis model of supra-glacial lake formation on debris-covered glaciers is based on GPR (Ground Penetrating Radar) (Mertes et al., 2017). The analysis has put forth the argument of increased melting observed in glaciers of northern Pakistan. The risk factor increases exponentially with the presence of supra-glacial lakes. The trend is observed through modeling and risk assessment of GLOFs (Lala, 2018). The frequency of glacier-dammed lakes and outburst floods in the Karambar Valley of Hindukush-Karakoram has increased the risk to infrastructure and living organisms in this region (Iturrizaga, 2005). The glacier surge is a seasonal phenomenon owing to the extreme flow velocities resulting in the formation of a dammed lake (Steiner et al., 2018). The balance in accumulation and ablation zones of a glacier is very vital for its stability. A hydrometeorological perspective on the anomaly of glacier dynamics has originated the argument of heavy accumulation zones, thus disturbing the mass balance (Bashir et al., 2017). The natural stability and behavior of the glacier are very much dependent on slope, elevation, aspect, and geomorphology of the vicinity. Glacier expansion is very much related to the elevation from the mean sea level in the Karakoram region (Hewitt, 1998). Himalayan glaciers are a focus of public and scientific debate. Prevailing uncertainties are of major concern because some projections of their future have serious implications for water resources. Most Himalayan glaciers are losing mass at rates similar to glaciers elsewhere, except for emerging indications of stability or mass gain in the Karakoram (Bolch et al., 2012). Rising global temperature is the major factor in the glacial lake formation, which is caused by the glacial retreat in mountainous regions. In the era from 1550 to 1850, the glaciers were quite significant in length in comparison with today. With the inception of global warming, moraines formation adjacent to glaciers blocks the glacial lakes (Bhutiyani, 1999). Since the Little Ice Age, it is said that the glaciers of the Himalaya have experienced a retreat of approximately one kilometer in length (Mool et al., 2001). A situation leads to the formation of moraine-dammed lakes with the provisioning of a large space for meltwater retention. The proximity analysis of the settlements with respect to glacial lakes is very vital with respect to geospatial analysis and modeling of GLOFs hazard in Hunza Valley of UIB. Nearly 35 devastating GLOF events have occurred during the last 200 years in Gilgit Baltistan (Din et al., 2014). The frequency and the intensity of GLOF events have risen over the past few years, according to available records. During the last three decades, the glacier cover has been decreased on an average of 10%, which has caused many GLOF events in Hunza Valley (Ali et al., 2019). Moreover, during the year (2008-2009), five GLOF events took place in Hunza Valley (Kreutzmann et al., 2011). A study of the GLOF events shows that such an event has been associated with weather conditions in terms of temperature increase, precipitation, and heatwaves.
Keeping in view the fact that watersheds of Pakistan are covered by major glaciers, which are quite susceptible to disastrous outbreak/flooding hazards, the objectives of the study were (i) to map potentially dangerous glacial lakes in Hunza Valley and (ii) to simulate the flood extents of the potentially dangerous glacial lakes using HEC-RAS model and do damage assessment to the downstream settlements.

Study area-Hunza valley of UIB
The Hunza River sub-basin drains the Karakoram Mountains, which are comprising of a large glaciated area in the north (Ashraf et al., 2012). The Karakoram highway linking Pakistan to China passes across this Valley. Part of the road runs along the Hunza River and ends near Khunjerab Pass (Geerken and Bräker, 2017). The tributaries joining the Hunza River are Chabursan, Khunjerab, Ghujerab, and Shunsha River. The sub-basin of the Hunza River comprises potentially dangerous Khurdopin glacial lakes on the course of the Shimshal River (Fig. 1). Karimabad, the capital of the Hunza Valley, is stretched over miles and miles of terraced fields and fruit orchards. The Valley encompasses Rakaposhi, Ultar, and Balimo peaks. Gulmit is shining white and deeply crevassed-just as you would expect a glacier to look. Above this glacier to the left is the jagged line of the Passu and Batura peaks, seven of which are over 7,500m. Passu is the setting-off point for climbing expeditions up the Batura, Passu, Kurk, and Lupgar groups of peaks, and for trekking trips up the Shimshal Valley and Batura Glacier (Singh, 2015;Kreutzmann, 2018).

Geomorphology
The high mountains of Pakistan comprise the western end of 2,400km long Himalayan range and some parts in the Hindukush and Karakoram ranges. Northern areas spread over 72,496km 2 with a midst towering snow-clad peaks having heights varying from nearly 1,000 to over 8,000 meters above sea level. Of the 14 over 8,000m peaks on the earth, 4 occupy an amphitheater at the head of Baltoro glacier in the Karakoram Range. These are K-2 (Mount Godwin Austen), which is 8,611m, and is the world's second-highest peak, Gasherbrum-I (8,068m), Broad Peak (8,047m), and Gasherbrum II (8,035m). In addition to these, there are 68 peaks over 7,000m and hundreds, which are over 6,000m high. Generally, because of their rugged topography and the rigors of the climate, the northern highlands and the Himalayas to the east have been formidable barriers to movement into Pakistan throughout history (Isserman et al., 2010).

Climate
Pakistan is basically a dry country of the warm temperate zone. The climate of the area is transitional between that of central Asia and the monsoonal region of South Asia, which varies considerably with latitude, altitude, aspect, and local relief. There is not only high spatial variability, but temporal variability is quite high as well. Except for a small strip of sub-tropical terrain in Punjab and the wet zone on the southern slopes of the Himalayan and Karakoram mountain ranges, most of the country is arid or semi-arid steppe land. The snowmelt run-off constitutes a substantial part of the water resources of the rivers of Pakistan (Singh et al., 2011). The Indus River, primarily supplied by glaciers in its upper reaches, and subject to the least seasonal variation, still has a maximum flow more than fifty times its minimum. Alpine glaciers contribute 50% of the Indus water flow. The Indus River is about 2,800km long, and 62% of its catchment lies in Pakistan (Singh et al., 2011). The swelling of the Indus and its tributaries is subjected to a volumetric decrease of glaciers and, if coupled with heavy monsoonal rains, can cause floods during summer.

Dataset
Landsat ETM+ Images of Hunza Valley, within the substantial time span from May to September, have been acquired from the USGS (United States Geological Survey) using the Earth Explorer interface (earthexplorer.usgs.gov), Digital Elevation Model of Hunza Valley is used for obtaining the elevation, aspect, and slope of the glaciers hosting the glacial lakes. ASTER interpolated data at 15m is used for this purpose, Geomorphologic data of Hunza Valley of Pakistan acquired from the Geological Survey of Pakistan.

Methodology
The Study encompasses the acquisition of Satellite Images, the performance of geospatial analysis, and identification of GLOFs to assess the glacial hazard-prone areas. By utilizing height information obtained through DEMs, orientation, and slope maps are formed (Nabi et al., 2018). The Landsat images of different time spans were downloaded and studied in detail for quality input. Capturing Digital Data of Glacial Lakes from Imagery Landsat images were used for the identification of glacial lakes by applying the Normalized Difference Water Index (NDWI), taking advantage of the low water reflectance in the NIR band.
Thereafter, mapping of Glacial Lakes of Hunza Valley of UIB in contact with glaciers and upstream of settlements was carried out. In this connection, a direct hydrological connection and lake dam type was determined. The lakes' volume was calculated based on the surface area. Finally, the simulation and modeling of potentially dangerous glacial lakes in HEC-RAS were conducted. The HEC-RAS model (Fig.  2) contains several river analysis components for steady flow water surface profile computations and one and two-dimensional unsteady flow simulation, including velocity and water surface depth analysis. The release of Version 5.0 introduced twodimensional modeling of flow as well as sediment transfer modeling capabilities. The program was developed by the US Department of Defense, Army Corps of Engineers in order to manage the rivers, harbors, and other public works under their jurisdiction; it has found wide acceptance by many others since its public release in 1995 (Osti and Egashira, 2009).

Results and discussion
The potentially dangerous glacial lakes which are concentrated at the headwaters of the UIB can affect settlements, infrastructure, and agricultural fields situated in the downstream River Valleys (Stäubli et al., 2018). The ability of decision-makers on the adaption of risk mitigation measures and reduction in vulnerability will be enhanced with a detailed digital data repository of glacial lakes and GLOFs occurrences. This forms the basis for global warming studies and future climate change research in Pakistan, as the irrigation network is primarily dependent on summer season snowmelt (Mukhopadhyay and Khan, 2015).

GLOFs hazard assessment
By using Landsat-7 images, the study of glaciers and glacial lakes is carried out coupled with field investigations of potentially dangerous GLOFs. Using remote sensing satellite images, the monitoring of the glaciers as per created inventories, and the impact assessment of the GLOFs extent is done precisely. The accuracy is achieved with the remote sensing data and techniques for the evaluation of geophysical conditions of the terrain with the help of satellite images. The ability and precision of the analysis performed is increased with the multistage approach of field investigation coupled with remote sensing dataset. The study involving glaciers and GLOFs becomes reliable once visual image interpretation techniques are integrated with GIS analysis.
For this research, the identification of glaciers and glacial lakes has been made by utilizing Landsat-7 ETM+ images. The bandwidth of TM and ETM+ is slightly different, ranging from the blue to the farinfrared wavelength. For feature identification, Landsat-7 band combinations and indices are utilized. The glacial lakes can be easily identified in the band combination of RGB (Red-Green-Blue) (Pan-7-6b) due to their better contrast with the surrounding features (Gilany and Iqbal, 2020). In this FCC, the fresh snow and ice of the glaciers appear in light to dark red color. In the image of the winter season, the glacial lakes with a smooth texture and varying gray tone due to their semifrozen ice surface are easily identified (Figs. 3a-3c). The identification/monitoring of glaciers and glacial lakes is made with the integration of remote sensing technique and GIS (Geographic Information System) data analysis. It played a major role in decision making and application of rules of land cover types and features discrimination in GIS analytical techniques, which enabled better presentation and perspective views. DEM is utilized to create slope and aspect data sets of the study area. Even though the glacial lake surfaces are flat and covered by snow, the glaciers ice and snow ice create slope angles (Gilany and Iqbal, 2016). The probability of a glacial lake outburst is a function of the basic susceptibility of the dam to fail and the potential for external trigger processes (Huggel et al., 2002). Antecedent, decision rules of integrated GIS analysis are applied, that if the surface texture is smooth and the slope is not pronounced, then such areas are recognized as glacial lakes.

GLOFs analysis of Hunza valley of UIB
There are 110 glacial lakes covering an area of about 3.22Km 2, out of which 47 are major glacial lakes in Hunza Valley of UIB. The high relief and unstable deposits along the Valley sides have made the slopes prone to mass movements. The Hunza Valley provides an ideal and easily accessible location for the study of ice-dammed and mass movement-dammed lakes. The largest-sized lake in this category has a 0.38Km 2 area and a length of 5Km. It is oriented towards the North West and is associated with the Khurdopin glacier. This potentially dangerous glacial lake is hazardous to the settlement of Shimshal Valley (Fig. 4a).

HEC-RAS model simulation of glacial lake outburst flood hazard risk in Hunza valley of UIB
The HEC-RAS model simulation is utilized for the identification of GLOF hazard risk extents to Shimshal village in Hunza Valley. Shimshal village is located in Gojal Tehsil of Hunza District, in the Gilgit-Baltistan region of Pakistan. It lies at an altitude of 3,100m amsl and is the highest settlement in Hunza Valley. It is a border village that connects the Gilgit-Baltistan with China. The total area of Shimshal is approx 3,800km 2, and there are around 2000 inhabitants with a total of 250 houses. The input parameters for the HEC-RAS model are listed in Table 1.

Khurdopin glacial lake inlet and Shimshal village domain area
First, the domain area surrounding the flow accumulation of stream flowing out of the potentially dangerous glacial lake is drawn. The habitat of Shimshal village is included in the domain area to ascertain the damage extends to the settlements. The two-dimensional domain area is assigned the pixel value of 15x15m. The inlet to the domain area is drawn at the outflow of potentially dangerous Khurdopin glacial lake (Fig. 4b).
The peak seasonal discharge from the potentially dangerous lake of Khurdopin glacier is calculated using the empirical formula of peak discharge (Costa, 1988).
where, Qmax=Peak Discharge (m 3 /s); Vo=Volume (m 3 ). Basing on the parameters (Table 2) of the potentially dangerous Khurdopin glacial lake (Fig.  4c), the peak seasonal discharge flow is calculated as 3500m 3 /s, which can generate low flood damage extent to the settlements of Shimshal village ( Table  2). The peak simulated scenario-1 discharge flow is calculated as 5000m 3 /s, which can generate moderate flood damage extent. The peak simulated scenario-2 discharge flow is calculated as 7000m 3 /s, which can generate high flood damage extent to Shimshal village (Table 2).

HEC-RAS simulated hydrograph max depth, max velocity, and max water surface elevation at x-sec of Shimshal village
Basing on a peak seasonal discharge of 3500m 3 /s and a data input interval of 10mins, the twodimensional hydrograph profiles are generated from the HEC-RAS model. The max depth of flow hydrograph at Shimshal village x-sec is calculated as 40m, which has generated a low flood to the settlements (Fig. 4d). The max velocity of flow is calculated as 7m/s (Fig. 4e). The reference surface elevation is 3070m amsl, and the max water surface elevation at Shimshal village x-sec is 3110m amsl (Fig. 4f).

HEC-RAS output parameters at Shimshal village x-sec (low, moderate, and high discharge)
Basing on the simulated 2D hydrograph profiles generated from the HEC-RAS model, the output parameter values of potentially dangerous Khurdopin glacial lake at Shimshal village x-sec generating low, moderate, and high discharge are as shown in Table 3.

Conclusion
The GLOFs play a vital role in sedimentation and erosion in UIB of Pakistan. Their significance cannot be denied, especially, which lies in very exceptional risk to the infrastructure and human installations. The historically recorded floods gain height well beyond peak discharge estimated values for the seasonal precipitation. The erosion capacity and competence are immensely enhanced by the active dynamic character of the GLOFs. In the context of erosion in these Valleys and the sedimentation of the reservoirs in the downstream area, the vital importance is of the GLOFs happening. The passage of this dam burst involving turbulent floods has contributed to huge numbers of landslides that have occurred in Valley sides and on the terraces of Hunza Valley. Keeping in view the seasonal growth of potentially dangerous Khurdopin glacial lake of Hunza Valley and simulated scenarios, the low discharge of 3500m 3 /s from the glacial lake can affect 40% of the Shimshal village habitat, the moderate discharge of 5000m 3 /s from the glacial lake can affect 60% of the Shimshal village habitat, and the high discharge of 7000m 3 /s from the glacial lake can affect 80% of the Shimshal village habitat. The Hunza Valley of UIB is prone to glacial lake outburst floods hazard based on the proximity of glacial lake with respect to infrastructure, geomorphology of underneath surface, geo-cover of the vicinity, crevasses, ice melt, and anthropogenic activities. Therefore, continuous monitoring through physical gauge stations and satellite images is very vital for the streams nearing settlements of Hunza Valley. Knowing the extent of damages beforehand can help in mitigating the impact of GLOF surges. Antecedent, the most vital mitigation step to reduce flood risk, is to gradually reduce the volume of the glacial lake to decrease the dynamic peak surge of glacial lakes containing a huge volume of water. In order to protect the infrastructure in downstream areas against the destructive/dynamic forces of surging GLOFs, pre-disaster mitigation measures must be taken. An early warning and monitoring system should be placed in advance in order to safeguard against such catastrophic events. While choosing the appropriate method or starting any mitigation measure, precise evaluation involving detailed analysis studies of lakes, mother glaciers, surrounding conditions, and damming materials are the foremost requirements. The measures adopted must be such that those must not increase the risk of a GLOF event during or after the placement of mitigation measures. At different stages of the mitigation process, i.e., during or after, the onsite monitoring gadgets at the mother glaciers, the lake, the dam, and the surroundings are very vital. cooperation and help in this research work. We acknowledge PMD (Pakistan) for the provision of data support for this paper.

Conflict of interest
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.