Debris-covered glacier systems and associated glacial lake outburst flood hazards: challenges and prospects

Glaciers respond sensitively to climate variability and change, with associated impacts on meltwater production, sea-level rise and geomorphological hazards. There is a strong societal interest in understanding the current response of all types of glacier systems to climate change and how they will continue to evolve in the context of the whole glacierized landscape. In particular, understanding the current and future behaviour of debris-covered glaciers is a ‘hot topic’ in glaciological research because of concerns for water resources and glacier-related hazards. The state of these glaciers is closely related to various hazardous geomorphological processes which are relatively poorly understood. Understanding the implications of debris-covered glacier evolution requires a systems approach. This includes the interplay of various factors such as local geomorphology, ice ablation patterns, debris characteristics and glacier lake growth and development. Such a broader, contextualized understanding is prerequisite to identifying and monitoring the geohazards and hydrologic implications associated with changes in the debris-covered glacier system under future climate scenarios. This paper presents a comprehensive review of current knowledge of the debris-covered glacier landsystem. Specifically, we review state-of-the-art field-based and the remote sensing-based methods for monitoring debris-covered glacier characteristics and lakes and their evolution under future climate change. We advocate a holistic process-based framework for assessing hazards associated with moraine-dammed glacio-terminal lakes that are a projected end-member state for many debris-covered glaciers under a warming climate.

Glaciers respond sensitively to climate change and variability, with associated impacts on meltwater production (Kaser et al. 2010; Huss 2011; Immerzeel et al. 2012), sea-level rise (Berthier et al. 2010; Leclercq et al. 2011; Marzeion et al. 2012) and geomorphological hazards (Kääb et al. 2005; Benn et al. 2012; Harrison et al. 2018). Glacier behaviour has potential knock-on effects for valley-scale sediment fluxes, surface energy balance, water storage and geomorphological hazards. Therefore, there is a keen societal interest in understanding how the different types of glacier systems are currently responding to climate change and how they will evolve in the context of the whole landscape. In particular, the role of debris-covered glacier landsystems and associated lakes in water supply and related hazards is less well understood. Although they account for only c. 4 to 7% of the global glacierized area (Scherler et al. 2018; Herreid and Pellicciotti 2020), debris-covered glaciers are a prominent feature of high-relief orogenic belts where high denudation rates supply abundant rock debris to the glacier surface, often producing debris-covered glacier tongues up to tens of kilometres long such as Baltoro Glacier (62 km) in the Karakoram (Mihalcea et al. 2006) or Ngozumpa Glacier (c. 25 km) in the Nepal Himalaya (Casey and Kääb 2012) (Fig. 1).They are an especially well-developed feature in the Hindu Kush Himalayan region, where c. 13% of the total glacierized area is debris-covered, ranging from 9% in the Karakoram to 15% in the eastern Nepalese and Bhutanese Himalaya (Kääb et al. 2012). They are also found in the Tien Shan (Hagg et al. 2008), Caucasus (Stokes et al. 2007), Alaska (Berthier et al. 2010), New Zealand (Anderson and Mackintosh 2012), parts of the Andes (Racoviteanu et al. 2008) and the Alps (Deline 2005), Greenland, and the Dry Valleys of Antarctica. Supraglacial debris varies in thickness from several centimetres up to 2 m or more (Benn and Evans 1998; Anderson and Anderson 2018).

Satellite-derived inventories show that glacier-wide ice mass loss from debris-covered glacier tongues over recent decades is substantial and increasing (Bolch et al. 2011; Kääb et al. 2012; Brun et al. 2019; Maurer et al. 2019; King et al. 2020b). As mountain glaciers continue to diminish in the coming decades, an increasing proportion of the remaining ice is expected to become debris-covered (Herreid and Pellicciotti 2020). This makes it critical to understand how debris cover impacts glacier meltwater production in order to make projections of regional water resources and global sea-level rise. Furthermore, mass loss from debris-covered glaciers in particular is closely associated with the formation of ice-contact and moraine-dammed lakes (Reynolds 2000; Benn et al. 2012; Sakai 2012; King et al. 2019). Their impact of lake evolution on local hazard potential in the context of future climate projections is still unclear (Harrison et al. 2018).

Over the last decade, debris-covered glaciers have become a ‘hot topic’ in glaciological research following concerns about the fate of glaciers, particularly in High Mountain Asia (Cogley et al. 2010; Bolch et al. 2012). During this time, several satellite remote sensing studies showed that thickly debris-covered glaciers display high rates of surface lowering, comparable to those of clean ice glaciers (Kääb et al. 2012; Nuimura et al. 2012; Gardelle et al. 2013; Pellicciotti et al. 2015; Brun et al. 2019), even though a thick debris mantle has been conclusively shown to locally reduce ablation compared to that of clean ice (Østrem 1959; Kayastha et al. 2000; Nicholson and Benn 2006). This mass loss has been attributed to modified ice dynamics (Vincent et al. 2016; Brun et al. 2018; Anderson et al. 2021; Rounce et al. 2021) and to localized ice ablation rates related to ice cliffs and ponds (Sakai et al. 2000b; Miles et al. 2018b; Buri et al. 2021). The complex surface topography of debris-covered tongues exhibits exposed ice cliffs (Steiner et al. 2015; Buri and Pellicciotti 2018), surface ponds of various sizes (Sakai and Fujita 2010; Watson et al. 2016; Miles et al. 2018b), debris cones/hummocks (Moore 2018; Bartlett et al. 2021), medial moraines (Anderson 2000), lateral and terminal moraines (Hewitt and Shroder 1993; Owen et al. 2003; Benn et al. 2004), supraglacial streams (Fyffe et al. 2019a; Miles et al. 2020), surface depressions (Mertes et al. 2016; Benn et al. 2017; Miles et al. 2017a), relict englacial conduits (Gulley and Benn 2007; Gulley et al. 2009b), base-level lakes (terminal, proglacial or supraglacial, and proto-lakes) (Benn et al. 2012) and supraglacial vegetation (Fickert et al. 2007; Tampucci et al. 2016; Anderson et al. 2020) (Fig. 2). Certain supraglacial features act as ‘hot spots’ for ice melt, particularly ice cliffs (Sakai et al. 2002; Han et al. 2010; Steiner et al. 2015; Buri et al. 2016, 2021) and supraglacial ponds (Sakai et al. 2000b; Miles et al. 2016, 2018b; Salerno et al. 2017). In spite of the prevalence of these features, the insulating effect of a thick debris cover is still predominant on many glaciers (Vincent et al. 2016; Brun et al. 2018; Anderson et al. 2021), and is evidenced at a mountain-range scale by debris-covered glaciers having lower terminus elevations than clean ice glaciers.

Local, regional and global patterns of glacier thinning and mass loss are coupled with an increase in debris cover extent as the upper limit of the debris cover migrates upglacier with the equilibrium line of the glacier (Deline 2005; Anderson and Mackintosh 2012; Herreid and Pellicciotti 2020; Xie et al. 2020a). Debris thickness increases due to cumulative exposure of englacial debris as glaciers thin due to surface ice ablation. While progress has been made in understanding local glacier-surface dynamics related to these supraglacial features, the extent to which the evolution of the debris-covered glacier surface influences overall glacier behaviour remains uncertain.

Insights into the surface characteristics of debris-covered glaciers and the enhanced local ablation rates gave rise to another concern relevant to both the scientific community and local communities, namely the accelerated growth of supraglacial and proglacial lakes associated with glacier thinning and recession documented around the world (Paul et al. 2007; Komori 2008; Gardelle et al. 2011; Thompson et al. 2012; Wang et al. 2014; Nie et al. 2017; Shukla et al. 2018; Chen et al. 2020; Shugar et al. 2020). Proglacial lakes exhibit an effect on glacier ice dynamics through enhanced ablation at the glacier margins via mechanical and thermal stresses; they modify meltwater routing and sediment fluxes through sedimentation (Carrivick and Tweed 2013). While proglacial lakes and their hazard potential have been addressed in several studies (Reynolds 1999; Watanabe et al. 2009; Aggarwal et al. 2017; Haritashya et al. 2018; Wilson et al. 2019), the link between supraglacial lakes and glacial hazards is less well studied.

There has been an increased interest in understanding the conditions for the formation and evolution of supraglacial ponds on debris-covered glaciers (Sakai and Fujita 2010; Sakai 2012). Some are highly dynamic, quickly evolving and growing, while others are persistent; some are short-lived (Miles et al. 2018a), while others may coalesce through time to create a larger moraine-dammed lake at the terminus of the glacier, adjoining a calving glacier front (Benn et al. 2012). The formation of moraine-dammed lakes is associated with glacier retreat and downwasting patterns (Korup and Tweed 2007; Carrivick and Tweed 2013) and is favoured by low glacier slope and velocities (Quincey et al. 2007) and/or changes in supraglacial debris flux (Benn et al. 2012). Water stored behind a weak moraine has the potential to breach the moraine dam, resulting in glacier lake outburst floods (GLOFs). These involve the rapid release of large volumes of water and sediments, with disastrous consequences for communities downstream. The Dig Tsho flood event in the Khumbu Himalaya in 1985 was one of the largest such events recorded (Vuichard and Zimmermann 1987; Mool et al. 2002).

The evolution of moraine-dammed lakes associated with debris-covered glaciers has been addressed in few studies (e.g. Bolch et al. 2008b; Thompson et al. 2012, 2017; Harrison et al. 2018). Assessing the hazard potential of these lakes presents significant challenges both in the field due the highly dynamic environment and in remote sensing due to limited multi-temporal, high-resolution data needed to estimate glacier surface evolution in some areas (Wang et al. 2020). Another important challenge is posed by the lack of systematic approaches for classifying and ranking proglacial and supraglacial glacier lakes in terms of their hazard potential. This results in significant gaps in glaciological and geomorphological characteristics of debris-covered glaciers and associated lakes, hindering the understanding of the future evolution of these glaciers and its implications for glacier hazards and water resources. A deeper and broader, contextualized understanding is prerequisite to identifying and monitoring the geohazards and evaluating the hydrologic implications associated with changes in the debris-covered glacier system under future climate scenarios.

Addressing this research gap requires a more complete estimate of glacier responses to climate change and their impacts, notably a better understanding of the debris-covered glacier landsystem, its components and the interplay of various factors such as local geomorphology, ice ablation patterns, debris characteristics, glacier lake growth and development. In order to achieve this, a systems approach is needed. Traditionally, debris-covered glacier processes are addressed by a single discipline such as geomorphology, glaciology or hydrology; very few studies have adopted the required holistic, systemic approach. Although glaciology is an increasingly interdisciplinary field, most scientists are driven to specialize and only study one or few aspects of the glacier system (e.g. melt processes, hydrology or ice dynamics). Benn et al. (2012) is one example of a holistic study that links climate, mass balance, ice dynamics, topographic evolution and hydrology in the Everest region of Nepal and explores how observed behaviour and hazard potential emerge from interactions between these process domains. Lake hazard assessments are often conducted from a remote sensing or geophysical perspective (Pant and Reynolds 2000; Rana et al. 2000; Richardson and Reynolds 2000; Reynolds 2006; Hambrey et al. 2008; Thompson et al. 2012, 2017). A systems approach combines various perspectives to provide a more complete picture of the debris-cover–lake system, i.e. the interplay between the glacier surface topography, the lake dynamics and the ablation patterns related to supraglacial debris cover and its characteristics. Furthermore, this type of holistic perspective allows a better understanding of the past, present and future states of debris-covered glaciers, as well as their position, role and consequences within landscapes. This is key for various applications, but in this paper we focus on its value for estimating the hazards related to rapidly evolving moraine-dammed glacier lakes and their impacts on populations. We identify misunderstandings related to these concepts and failures in the way they have been communicated, and suggest ways to bridge these gaps to develop an understanding of resilience to climate change.

In this section we consider the question ‘what is a debris-covered glacier?’ from a landsystem perspective. In doing so, our aim is to widen the focus to understand how the glacier and landscape interact. In light of this, we must consider a debris-covered glacier as a system with a particular assemblage of features, associated with the higher rate of debris loading and exhumation than a typical valley glacier system. For example, while a simple and commonly used definition for identifying a debris-covered glacier is one where debris material covers the full width of part of its ablation zone (Kirkbride 2011), here we define a debris-covered glacier as one whose surface mass balance is sufficiently affected by supraglacial debris as to alter the glacier geometry, ice dynamics, surface features and hydrology (or some subset of the above) compared to that of a clean-ice glacier. The relative debris-richness of the glacier system controls its propensity to become debris-covered (Kirkbride 1989). Placing debris-covered glaciers on a continuum of mixed ice/debris landforms in this way ties the degree of debris cover to the relative abundance of snow v. debris supply. From this viewpoint, negative glacier mass balance conditions in the presence of abundant debris are expected to lead to development of a debris cover.

Rock debris is supplied to the glacier surface by gravitational mass movements from the surrounding terrain (e.g. Nakawo et al. 1986; Hambrey et al. 2008; Nagai et al. 2013), which generally occur as isolated events in space and time, and are poorly represented by the scant available measurements (Deline 2009; Hewitt 2009; Reznichenko et al. 2011) or by estimates derived from long-term headwall retreat rates (Heimsath and McGlynn 2008; Seong et al. 2009). Debris deposits onto the glacier surface within the accumulation zone are buried and transported englacially until they are exhumed to the surface by ice melt within the ablation zone (Kirkbride and Deline 2013). Deposits directly onto the ablation zone remain at the surface. Theoretical considerations and modelling studies (e.g. Rowan et al. 2015; Anderson and Anderson 2016; Wirbel et al. 2018; Scherler and Egholm 2020) highlight that the specific location of debris inputs strongly influences the spatial pattern of supraglacial debris. A constant rate over a longer time period cannot produce a localized, high debris concentration within the glacier, but will lead to an extended zone of lower concentration, which will produce distinctly different surface debris patterns compared to an initially localized zone of high concentration (Wirbel et al. 2018). The emergence of debris in the ablation zone is governed by the debris input location, englacial transport and melt-out rates. Thus, to accurately compute the point of emergence and thickness of debris at this melt-out location, the englacial transport pathways and deformation must be captured in some way (Kirkbride and Deline 2013; Wirbel et al. 2018). Once on the glacier surface, debris is transported by advection with underlying ice flow, where gradients in the ice surface velocity and resulting zones of compressive or extensional ice flow will thicken or thin the debris cover layer. In addition, as soon as debris emerges on the glacier surface, it is affected by other processes such as gravitational reworking (Anderson 2000; Kirkbride and Deline 2013; Moore 2018; Fyffe et al. 2020). Irregular englacial debris concentration and subsequent surface reworking causes local debris thickness variability (Moore 2018; Nicholson et al. 2018; Westoby et al. 2020), leading to strong small-scale variability in ablation rates and the formation of pronounced surface relief and features (Benn et al. 2001; Mertes et al. 2016).

The supraglacial debris mantle has a profound influence on the underlying ice ablation rate. Early studies (Østrem 1959) showed that glacier surfaces with patchy or very thin cover of supraglacial debris experience accelerated ablation, whereas a continuous debris cover of thickness greater than a few centimetres inhibits the underlying ice ablation. The specific effect of debris on ice melt is influenced by its individual characteristics such as debris thickness, debris material, porosity, grain size, moisture and liquid water content and the prevailing meteorological conditions. However, field studies (e.g. Fujii 1977; Mattson et al. 1993), laboratory experiments (Reznichenko et al. 2010) and modelling studies (e.g. Nakawo and Young 1982; Nicholson and Benn 2006; Reid and Brock 2010) demonstrate that debris thickness is the primary determinant of how sub-debris ice ablation rates differ to clean-ice melt rates, with the properties of the debris layer playing secondary roles (Reznichenko et al. 2010; Nicholson and Benn 2012; Collier et al. 2014). As surface debris is continuously conveyed downglacier with ice flow, debris cover thickness increases towards the glacier terminus (Rowan et al. 2015; Anderson and Anderson 2018). This profoundly alters the glacier-scale ablation regime, in principle causing an inversion of the ablation gradient toward the terminus such that maximum ablation occurs some distance upglacier instead of at the terminus as is the case for clean-ice glaciers (Benn and Lehmkuhl 2000; Bisset et al. 2020). This in turn has consequences for ice dynamics as the ablation gradient influences the development of the glacier surface longitudinal profile and thereby the driving stresses through the ablation zone.

Debris can be removed from the system by marginal deposition or by surface meltwater. Some debris-covered glaciers form large, impounding latero-terminal moraine complexes (e.g. Benn and Owen 2002; Hambrey et al. 2008) while other debris-covered glacier termini end in outwash plains without substantial terminal moraines (e.g. Mayer et al. 2006). Large terminal moraines affect the englacial water table, and increase the potential for water to be stored behind this impounding moraine (e.g. Benn et al. 2012). In the absence of impounding latero-terminal moraines, terminal lakes can only form by external geomorphological processes, for example where water courses are impounded by advances of neighbouring glaciers, or slope failures damming the valley downstream (Rashid et al. 2020). Large latero-terminal moraines also inhibit the evacuation of debris from the glacier surface by gravitational processes, and exert a physical constraint on the glacier terminus position and upstream ice dynamics. The sedimentological, geomorphological and dynamic context of debris-covered glaciers has been discussed by Hambrey et al. (2008). They presented a conceptual model for the eastern Himalaya applicable to other glaciers to explain the development of large, lateral-terminal moraine complexes and associated moraine dams. The presence or absence of confining moraine dams may play a decisive role in determining the end-member of glacier development under declining ice content. If present, they facilitate the formation of terminal lakes, while their absence may allow the transition of debris-covered glaciers to rock glaciers or other ice-debris landforms (Whalley and Martin 1992; Jones et al. 2019).

Understanding the future states of the debris-covered glacier landsystem requires the knowledge of the location of such glaciers, as well as their current extent and state. This has been the subject of previous mapping efforts, including regional and global estimates of supraglacial debris cover (e.g. Scherler et al. 2018; Herreid and Pellicciotti 2020). However, a complete understanding of the system also requires information on the fundamental debris surface characteristics in order to understand the ice ablation processes, the velocity and dynamics, the evolution of ice cliffs and ponds and their importance to hydrology and hazard events. In this section we focus on state-of-the-art techniques for mapping and monitoring these characteristics, focusing in particular on supraglacial debris cover extent, debris thickness, physical properties and associated surface features. For each, we present both remote sensing and field methods, we outline current advances, remaining gaps and challenges, and offer recommendations to overcome these.

Mapping of debris-covered glaciers received considerable attention in the late 2000s and early 2010s, with important improvements in monitoring capacity as satellite imagery improved in both spatial and temporal resolution and coverage. The release of the Randolph Glacier Inventory (RGI; Pfeffer et al. 2014; RGI_Consortium 2017) and subsequently the global supraglacial debris cover datasets constructed on the basis of the RGI (Scherler et al. 2018; Herreid and Pellicciotti 2020) enabled a step-change in understanding the distribution of debris-covered glaciers at the large scale. However, while these global databases provide an initial and global perspective of glacier and supraglacial debris extent, they suffer from several limitations. The RGI's composition from distinct sources means that both datasets suffer from inconsistent methods between and often within regions, varying representative dates, user-subjective post-processing and manual delineation (e.g. Paul et al. 2013), geolocation or projection errors, and occasional inclusion of spatially descriptive but not explicit sources (e.g. World Glacier Inventory). These problems were partially mitigated in a manual revision limited to glaciers larger than 1 km² (Herreid and Pellicciotti 2020), but such an effort is laborious for repeat application at the global scale. Consequently, although debris-covered glaciers may be accurately represented in available databases for well-known sites that have been mapped carefully, their representation may be inconsistent and may occasionally commit and/or omit entire features for areas that are less well surveyed (Racoviteanu et al. 2021). Thus, while acknowledging that the mapping of supraglacial debris within the bounds of the RGI (e.g. Scherler et al. 2018; Herreid and Pellicciotti 2020) constitutes an important advance in high-level understanding of debris-covered glacier distribution globally, we consider that the current representation of debris-covered ice within the RGI is not sufficiently robust or consistent for understanding debris-covered glacier processes and change. Accurate, large-scale mapping of debris-covered glacier tongues at multi-temporal resolution remains a gap.

In optical remote sensing, identifying the glacier boundary is difficult due to the spectral similarity of supraglacial debris to the surrounding moraines (Racoviteanu and Williams 2012). Previous remote sensing studies have used a combination of terrain information, spectral information and terrain curvature (Bishop et al. 2001; Paul et al. 2004; Bolch et al. 2007; Shukla et al. 2010a; Kamp et al. 2011) to map debris-covered glaciers. Recent studies combined these criteria in machine-learning algorithms in order to automate the mapping process (Robson et al. 2015; Zhang et al. 2019; Xie et al. 2020b; Holobâcă et al. 2021). Opportunities for method development have greatly improved in the past decade with the increased availability of new operational, rapid-repeat, and public satellite imagery (e.g. Sentinel). In addition to optical remote sensing, several other proxies offer promise and overcome its limitations, among which we mention the following.

• oSurface motion derived from pairs of satellite optical or radar images helps identify active debris-covered areas (Gardner et al. 2018; Dehecq et al. 2019).

• o Sequential synthetic aperture radar (SAR) coherence images which indicate changes in surface backscatter between repeat observations; SAR helps to identify the active parts of the debris-covered glaciers due to their motion-related decorrelation compared to the highly coherent surrounding areas (Strozzi et al. 2010; Frey et al. 2012; Robson et al. 2015; Lippl et al. 2018; Holobâcă et al. 2021). While SAR coherence images are widely available and overcome the limitations posed by cloud cover in optical remote sensing, wide application of SAR techniques is hindered by complex processing.

• oSatellite thermal imaging helps distinguish the debris underlined by glacier ice and the surrounding non-ice moraines based on the brightness temperature difference (Taschner and Ranzi 2002; Shukla et al. 2010b; Bhambri et al. 2011; Racoviteanu and Williams 2012; Alifu et al. 2015).

• oDigital elevation model (DEM) differencing, derived from topographic maps for example (GLAMOS 2020; Linsbauer et al. 2020) serves to identify surface lowering. This is based on the concept that even where debris cover is thick, some heat is generally transmitted through the debris layer seasonally, leading to a small amount of ice melt and thus resulting in surface lowering.

• oSurface roughness and characteristics including pitted, hummocky topography with sharp breaks in slope, incised channels, etc. from a high-resolution DEM (King et al. 2020a) help to identify debris characteristics that may differ quantitatively from other land surfaces;

• oObject-oriented and machine-learning techniques (OBIA; Robson et al. 2015; Khan et al. 2020) based on shape is a complement to debris-cover mapping procedures.

The challenge in remote sensing mapping of debris cover is how to capitalize on new monitoring tools to develop targeted repeat monitoring in a systematic manner. Currently, most of the methods presented above still need manual post-processing of glacier shapes to ensure that glacier outlines conform to our understanding of ice flow and known landforms. These studies remain limited to specific regions; there remains a need for an automated tool, workflow or set of methods to delimit debris-covered ice independent of the bounds of given glacier polygons. In order to achieve this, the above proxies should be used in addition to optical remote sensing to construct a method that should be transferable to different sites, should rely on freely available, global coverage data sources, and should be validated against the best available field data.

Field-based delineation of debris cover extent is generally difficult. Many debris-covered glaciers occupy remote, rugged domains, with limited field access and the precise boundary of a debris-covered glacier is challenging to identify in the field. Validating the above remote sensing methods is particularly difficult as field methods for debris-covered glacier extent mapping are also not in an advanced state and detailed field studies are site-specific. Promising field-based methods based on ground temperatures, geophysics, drone-deployed optical and thermal imaging and time-lapse photography are time consuming, site-specific and difficult to extrapolate over larger areas. Simple recognition of debris-covered glacier surface features in the field (cliffs, ponds, etc.) can often be helpful to identify the presence of sub-debris ice, but it is often easier to (subjectively) identify debris-covered glacier extent from a planimetric perspective (e.g. high-resolution optical satellite data, ideally with multitemporal data) than in the field.

Surface debris thickness modulates the surface temperature and exerts a physical control over sub-debris melt rates (Østrem 1959; Nicholson and Benn 2006); it is probably the most crucial but also the most difficult debris cover property to quantify and monitor; field measurements of debris cover thickness are difficult to obtain in the field due to the rugged terrain. Therefore, satellite remote sensing approaches have been increasingly used in recent decades to overcome this challenge.

• oThermal imagery: a variety of approaches of varying complexity have used satellite thermal data to estimate debris thickness (Boxall et al. 2021). These range from simple band thresholding (Ranzi et al. 2004) to exponential curve fitting based on the empirical relationship between surface temperature and thickness (Juen et al. 2014; Kraaijenbrink et al. 2018) or energy-balance inversion, often requiring model spin-up (Mihalcea et al. 2008; Zhang et al. 2011; Foster et al. 2012; Rounce and McKinney 2014; Schauwecker et al. 2015; Rounce et al. 2021; Stewart et al. 2021). An intercomparison of these methods is needed and has been identified as a research target for the IACS Debris-covered Glaciers Working Group (https://cryosphericsciences.org/activities/wgdebris/).

• oElevation change/surface mass balance: since the thickness of debris moderates energy transfer to the ice, it also controls ice melt rates (Østrem 1959; Nicholson and Benn 2012). Surface mass balance data can thus be used to invert an energy mass balance model for debris thickness (Ragettli et al. 2015; Rounce et al. 2018), although this often requires careful consideration of ice dynamics to estimate surface mass balance from elevation change, and the long-duration melt modelling is computationally expensive (Rounce et al. 2021).

• oPolarimetric SAR: as certain wavelengths of radar can penetrate into the debris surface, the attenuation of radar signals is indicative of debris thickness. Debris cover thickness can be estimated based on inversion of the volume scattering power and other parameters after target decomposition. This method is in its infancy, but shows promise for an independent assessment of debris thickness (Huang et al. 2017).

Debris thickness can be inferred from proxy data in remote sensing studies. Regional and global applications of these methods (e.g. Rounce et al. 2021) represent a key advance towards including explicit inclusion of the effects of debris cover in global glacier simulations. However, as with debris-covered glacier extent, methods for determining debris thickness require validation, which is usually achieved through field measurements. We thus recommend the continued investigation of local-scale debris thickness, and the compilation of a database of available debris thickness measurements. This is another important aspect of the IACS Debris-covered Glacier Working Group, which has begun to assembly of a repository of debris-related measurements via the Zenodo data repository (https://zenodo.org/communities/iacswgondcgs). Each dataset is assigned a unique DOI to ensure that the responsible parties receive appropriate credit.

Field measurements of debris cover thickness are all spatially limited and labour-intensive to obtain. Manual excavation (e.g. Zhang et al. 2011) is practically limited to debris cover thickness < 0.8 m; high-frequency ground-penetrating radar (GPR) (e.g. McCarthy et al. 2017) is hindered by the difficulty of deploying GPR equipment at remote, high-elevation sites; estimates based on oblique terrestrial photography (e.g. Nicholson and Mertes 2017) require a number of crude geometric assumptions; estimates based on unmanned aerial vehicle (UAV) or terrestrial thermal imagery (Steiner and Pellicciotti 2016; Kraaijenbrink et al. 2018) have large uncertainties associated with image processing, due to local surface temperature variations as a result of shading or moisture and signal saturation at thickness >∼30 cm (Steiner et al. 2021); calculation from field measurements of surface lowering by any terrestrial data acquisition (e.g. structure from motion DEMs, TLS/LIDAR, terrestrial radar) are difficult to upscale to the glacier extent. In general, while these methods all have great potential to complement regional remote sensing studies, they currently have limited application due to logistical difficulties and limited spatial extent.

Remote sensing methods are key for estimating ice velocities, and can be applied to the surface of debris-covered glaciers. Several studies used remote-sensing surface velocities to show the deceleration of ice flow downglacier leading to stagnation at the snout (Quincey et al. 2007; Hambrey et al. 2008; Haritashya et al. 2015). Flow velocities can be derived by feature tracking using satellite imagery such as ASTER, Landsat series or Sentinel (e.g. Berthier et al. 2003, 2005; Kääb 2005; Scherler et al. 2008; Dehecq et al. 2015, 2019; Millan et al. 2019) and established image coregistration methods (Leprince et al. 2007). Methodological advances and data availability led to the globally comprehensive and temporally dense multi-sensor record of land ice velocity from the Inter-mission Time Series of Land Ice Velocity and Elevation project. However, the spatial resolution of these data (120–240 m) remains an issue, as the data have limited application for monitoring narrow debris-covered glacier tongues. Recent databases with improved spatial resolution (50 m) (Millan et al. 2019) offer promise for monitoring of debris-covered surfaces, but this is limited by cloud cover. As for the debris thickness estimates, SAR can provide high-accuracy measurements of the direction and intensity of glacier flow in all weathers (Kumar et al. 2011) provided that corrections are applied to mitigate attitude effects and sensor distortions (Scherler et al. 2008).

Beyond debris extent and thickness, debris properties such as lithology, grain size, porosity, stratification and stability (Table 1) (Casey et al. 2012; Casey and Kääb 2012; Juen et al. 2013) are important for specific applications. Local field mapping of these properties is difficult; thus measurements are scarce. At a glacier or regional scale, many debris cover properties and features are more easily mapped using high-resolution satellite imagery than in the field. Therefore, in this section we only discuss the remote sensing mapping and monitoring of these features as summarized in Table 1, with a focus on ice cliffs, ponds and streams.

Although no formal definition for ice cliffs exists (Kneib et al. 2020) these are readily identifiable in the field as high-relief bare-ice areas interrupting the supraglacial debris layer and are often associated with a supraglacial pond. An increasing number of studies are using remote sensing techniques to identify ice cliffs from satellite data (Table 2). However, robust and transferable methods for mapping ice cliffs in a consistent manner are in their infancy. In addition, more studies are needed simply to assess the long-term changes in prevalence of ice cliffs, as well as the spatial differences in ice cliff occurrence. Of critical consideration for the above methods is the spatial resolution, and how resolved elevation models need to be to sufficiently represent ice cliffs. For example, high-resolution DEMs (c. 10 m spatial resolution) are available now at regional or global scales, some at no cost. These include the High Mountain Asia (HMA) DEM at 8 m (Shean 2017), ArcticDEM at 2 m (Noh and Howat 2015) or the TanDEM-X DEM at 12–30 m (Wessel et al. 2018) spatial resolutions. However, the spatial resolution of some of these DEMs, particularly TanDEM-X as well as other commonly available ASTER GDEMs or SRTMs (30–90 m) is not sufficient for mapping ice cliffs, which are often only a few metres wide. One possible way forward might be to validate a topographic proxy for ice cliff density and area, as nadir-view satellite imagery will have difficulty representing the total area of steep ice cliffs accurately.

These small superficial water bodies are important indicators of the debris-covered glacier's drainage system, i.e. they control the rate at which meltwater derived from the melting ice flows downstream (Irvine-Fynn et al. 2017), and they contribute to ice mass losses themselves (Sakai et al. 2000b; Miles et al. 2018b). Supraglacial ponds are considerably easier to identify in satellite imagery than ice cliffs, meaning that a number of properties can be targeted, including: (i) surface temperature (from satellite, UAV or terrestrial thermal imagery); (ii) lake volume (via sonar or topographic sink analyses to derive volume-area relationships); (iii) lake turbidity (blue index or with sub-pixel spectral analyses); (iv) changes in elevation (high-accuracy DEMs). Supraglacial ponds and their properties have received focused study over the past few years, largely with satellite data, and the optical methods to map them are well established (Table 3). Overall, few detailed field studies of supraglacial ponds exist, and more direct observations of ponds, their characteristics and their dynamics are still needed using a combination of the methods briefly outlined here. Contemporary satellite imagery can answer some of the current questions related to supraglacial ponds, including their seasonality and persistence, but efforts are needed to assess both properties and processes at the local scale, as well as their prevalence and change at the regional scale.

The inverted ablation gradient and low longitudinal gradient of debris-covered glaciers can have a strong impact on the structure and function of the glacier's entire drainage system (see the review by Miles et al. 2020). Supraglacial hydrology is directly observable in the field and with satellite data. Areas of thicker debris and lower debris are typically characterized by small catchments and discontinuous, low-efficiency drainage systems conducive to formation of supraglacial ponds (Miles et al. 2017b; Fyffe et al. 2019b) whereas areas of thinner debris and higher surface gradient can support larger catchments and efficient supraglacial stream systems (Gulley et al. 2009a; Miles et al. 2019). The relative extent of these domains is indicative of the glacier's decay and progression to stagnation (Benn et al. 2017; King et al. 2020a), but also important for understanding the diurnal and seasonal evolution of glacial discharge (Fyffe et al. 2019a). As supraglacial streams exist only where stream incision exceeds the background ablation rate (Marston 1983), they by definition directly contribute to melt; they also contribute indirectly to melt by promoting ice cliff development (Mölg et al. 2020; Kneib et al. 2021). Streams can be mapped using hydrologic analysis tools on high-precision, high-resolution topographic datasets derived from satellite stereo or UAV images (Benn et al. 2017; Miles et al. 2017b; Fyffe et al. 2020). Other efforts have mapped streams manually from satellite images or by walking their length in the field (Miles et al. 2019). Mapping supraglacial streams with DEM drainage analysis and optical imagery is similarly challenged by apparent stream discontinuities due to ice arches and flow through debris. Despite their importance for characterizing glacier drainage systems and debris-covered glacier stagnation, supraglacial streams have been addressed in fewer detailed studies than englacial conduits. Newly available high-resolution satellite images and DEMs offer the potential to better characterize supraglacial streams, but additional field investigations are needed to produce a generalized quantitative model of debris-covered glacier drainage efficiency.

Proglacial lakes have been studied through long-term and regional-scale monitoring efforts (e.g. Zhang et al. 2015; Nie et al. 2017; Shugar et al. 2020) and have been mapped in a systematic manner using established methods based on historical multispectral (optical) imagery (Fig. 3) using various water indices (Zhang et al. 2018; Zhao et al. 2018) or manual digitization (Wilson et al. 2018). In general, proglacial lakes are easier to map from remote sensing than supraglacial lakes due to their larger size. A number of these lakes have bathymetric field surveys undertaken to assess to their hazard potential (Worni et al. 2013; Haritashya et al. 2018). Here we note a few specific aspects of proglacial lake mapping that are important to consider for future studies.

• Shadows cast across the lake surface have been historically problematic for automatic lake mapping efforts, but indices have emerged recently that show improved performance (e.g. Chen et al. 2013). More problematic are the shadows that are cast across non-water surfaces, which sometimes alias as water (this is also a challenge for supraglacial pond mapping) (see Gardelle et al. 2011; Miles et al. 2017a). Some strategies to mitigate this include topographic shadow casting for correction or the use of multi-temporal data to filter out shadows not associated with water.

• Water turbidity can cause a varying spectral signal across the surface of a glacial lake, and is also useful to observe as indicative of water circulation patterns, discharge plumes and bulk suspended sediment (Wessels et al. 2002; Kraaijenbrink et al. 2016). However, the combination of surface ice and varying turbidity can cause problems for automated algorithms. Nonetheless, automated determination of surface water turbidity from satellite imagery has been accomplished in other regions (Matta et al. 2017), and could be transferred to glacial lake assessments.

• Surface ice (lake ice and icebergs) are extremely useful indicators of environmental conditions and processes, but are usually confounding factors for automated methods. Calving rates in particular have been a key target of study in marine environments and for emergent lakes, but have received relatively little attention compared with debris-covered glacier or high-mountain glacier lakes.

• Surface temperature can be useful for the delineation of large lakes, as resolution matters less, and is also a useful property itself as a controlling factor for lake water stratification and circulation (along with turbidity) and to understand the energy balance of the ice–lake–stream domain.

Of the methods mentioned in Table 1, we note here that SAR intensity mapping is extremely promising for glacial lake monitoring efforts, especially in cases where lakes are undergoing rapid change. Unlike optical data, SAR intensity mapping is insensitive to clouds and shadows, and less sensitive to turbidity and surface ice factors (e.g. Strozzi et al. 2012; Wangchuk and Bolch 2020). Furthermore, due to the size of most glacial lakes, the SAR intensity method is less affected by topographic and resolution issues than for supraglacial ponds and ice cliffs.

While we recognize that supraglacial features (ice cliffs and supraglacial ponds) and their dynamics are important to understand, it is useful to consider each in terms of their causal or controlling processes. Some features are associated with debris dynamics (e.g. differential ablation, debris deposition or emergence), while others are hydrologically/fluvially associated. We consider that mapping strategies should be driven by a specific research question. For example, ice cliffs may be studied as an indicator of debris redistribution or they can be regarded as exposed ice within the debris-covered domain in order to more accurately represent ablation rates. If areas of high meltwater production are of interest, identifying all the exposed ice on the debris cover might be more important than specifically delineating only ‘ice cliffs’. This would include other ice surface features such as ice sails (Evatt et al. 2017).

For GLOF hazard assessments, there is clearly a need for widespread screening of proglacial lakes and supraglacial ponds on a regional scale using at least semi-automated if not fully automated techniques. At local scales, proglacial lakes can be monitored in the field using bathymetric surveys (Cook and Quincey 2015; Watson et al. 2018b). So far, relatively few lakes have been the subject of such surveys and this needs further addressing. Multiple recent, current and future satellite altimeters (IceSat, Cryosat, IceSat2, SWOT) are promising avenues of operational workflow development. IceSat2 penetrates within many shallow-water bodies and might be suitable for bathymetric mapping of some proglacial lakes, although water turbidity of proglacial lakes limits the optical penetration depth; this needs direct analysis to test its feasibility.

The global pattern of glacier recession (e.g. Kargel et al. 2014) and the global nature of climate warming indicate a clear attribution to climate change (e.g. Marzeion et al. 2014; Zemp et al. 2015; Roe et al. 2017). Clearly, the mass balance of a glacier is causally linked to changes in temperature and precipitation, with accelerated negative trends of mass loss in the 21st century (Zemp et al. 2015; Solomina et al. 2016; Hugonnet et al. 2021). However, at regional scales, glaciers exhibit contrasting patterns in their response to climate changes (Sakai and Fujita 2017) due to differences in local topo-climatic factors (Salerno et al. 2017; Brun et al. 2019). Furthermore, local meteorology is usually not precisely known due to scarce measurements, leading to simplified modelling optimization schemes (Hock 2003). Complicating variables for mass accumulation include the addition of snow avalanches to mass balance and the importance of wind-blown snow from surrounding catchments.

With regard to ice ablation, site-specific losses occur via dynamic processes such as calving related to ice flow and glacier surface characteristics. In addition, surface energy balance and related melt and sublimation losses are driven by spatiotemporally varying fields of potential insolation, temperature, cloudiness, relative humidity and wind, all of which can manifest very differently depending on glacier settings and surface conditions and which are not easily characterized (Huo et al. 2021). However, these processes apply primarily in cases where the glacier surface is predominantly composed of exposed bare ice. It is observed in many mountain ranges that thickly debris-covered glacier termini persist at lower elevations than clean-ice glaciers. This indicates that the behaviour of a glacier terminus position in response to any given set of climate conditions is markedly different when the glacier has a surface debris cover compared to a clean-ice surface (Anderson and Anderson 2016). The geomorphological sensitivity of debris-covered glaciers is therefore an important and relevant concept (see Harrison 2009). While the geomorphological sensitivity of a clean-ice glacier could be established as its mass balance change over time and related to local climate change, it is not clear how we might assess the sensitivity of a debris-covered glacier, nor which metrics might be important.

The full response of debris-covered glaciers to climate forcing remains poorly understood in relation to that of clean-ice glaciers. One possible response of some mountain glaciers to climate change will be a transition from clean glaciers to debris-covered glaciers, and a further transition to rock glaciers in response to paraglacial processes increasing debris fluxes to glacier surfaces (see Monnier and Kinnard 2017; Jones et al. 2019). The long-term consequences of this transition are still largely unknown. In general, debris-covered glaciers pose a complicated case, where their behaviour and evolution are additionally related to non-climatic processes such as changes in debris flux from surrounding mountain sides or the presence of surface features such as ponds and ice cliffs. As a result, the system controlling the evolution of the debris-covered glacier system is not solely climatic in origin, but one in which paraglacial processes play an important role (Ballantyne 2002; Knight and Harrison 2014). Numerical models of glacier response to climate forcing under negative mass balance conditions suggest that debris-covered glaciers initially respond slower than clean-ice glaciers. However, the ultimate response time of debris-covered glaciers might be greater, as eventually the stagnant remnant of the glacier tongue detaches and decays in situ (Banerjee and Shankar 2013). The climate response of debris-covered glaciers is thought to be markedly asymmetrical between negative and positive mass balance conditions, with glacier adjustment rates to positive conditions matching those of clean-ice glaciers (Banerjee and Shankar 2013), such that the glacier length preferentially remembers positive mass balance phases over negative ones (Ferguson and Vieli 2020). There are very few observations of debris-covered glacier response to positive mass balance conditions (e.g. Deline 2005; Mölg et al. 2019), so the understanding gleaned from these modelling studies is unverified. However, it has been observed that substantial glacier advances can be triggered by extensive rockfall onto the glacier ablation zone. For example, following a large rock avalanche in 1920, the Brenva Glacier advanced 490 m between 1920 and 1941, whereas neighbouring glaciers in the Mont Blanc massif receded from the mid-1920s (Deline et al. 2015). This further highlights that the length of a debris-covered glacier is not a simple proxy for climate conditions alone. A better understanding of such processes is needed for long-term regional and global projections of glacier behaviour that form the basis of understanding trajectories of future meltwater availability and sea-level contribution from mountain glaciers (e.g. Kraaijenbrink et al. 2017; Rowan et al. 2017; Shannon et al. 2019).

Our understanding of how debris-covered glaciers and related landforms will evolve in the future remains limited. This means that the impact of climate change on these ice-debris systems will vary as the systems change. Viewed from the landsystem perspective, a debris-covered glacier landsystem incorporates numerous processes that respond to climate in different ways over time. This process transience of the system components presents a key challenge in simulating coupled glacier–climate behaviour (Nicholson et al. in press). For example, warming might be expected to cause a monotonic shift in precipitation phase from solid to liquid (i.e. more precipitation falls as rain rather than snow), starving the glaciers of snow accumulation while simultaneously enhancing ablation by rainfall. However, debris supply rates may show a complex non-linear response to the same warming over time. For a debris-covered glacier, the debris cover characteristics change in time as a function of supply, transfer, melt-out, thickness distribution and removal. These processes all co-evolve over time in a manner that is dependent on how the glacier geometry and ice flow dynamics adjust to the debris-modified spatial pattern of ablation. As a result, inter-relationships between these system components observed thus far might not hold into the future, and this non-stationarity means that such relationships are subject to both lags in response as well as gradual and thresholding process change, which are challenging to incorporate into a model system capable of reproducing system development over time.

Surface debris supply rate on debris-covered glaciers can be enhanced by debuttressing of rockwalls exposed by glacier recession, which can cause weakening of the valley walls and slopes. The timescale and duration of this effect is difficult to constrain and contingent on many structural, lithological and geomorphological conditions (Knight and Harrison 2018; Mancini and Lane 2020). In the longer term, debris supply may be more controlled by the rockwall area that lies within the freeze–thaw zone (Nagai et al. 2013; Banerjee and Wani 2018) and can also be influenced by heatwaves and heavy rainfall events. Secondary debris supply from debuttressed lateral moraines is an additional non-stationarity that is interesting to grapple with (van Woerkom et al. 2019). The system debris content is also affected by debris evacuation rates, which is primarily governed by the nature of the terminal deposition environment. Debris-covered glacier termini ending in outwash plains (e.g. Mayer et al. 2006) can export sediment to the foreland, while those with large, impounding latero-terminal moraine complexes (Benn and Owen 2002; Hambrey et al. 2008) cannot readily do so. Changing debris load over time will influence, together with changing ice inputs and losses, how and when debris-covered glaciers can form, and when they might transition to ice-cored rock glaciers, for example, due to increasingly inefficient supraglacial sediment evacuation (Monnier and Kinnard 2017; Jones et al. 2019; Knight et al. 2019).

The characteristic downglacier increase in debris thickness (Anderson and Anderson 2018), and the associated ablation gradient inversion toward the glacier terminus implies that maximum ablation occurs at the upper part of the debris cover and is reduced downglacier (Benn and Lehmkuhl 2000; Bisset et al. 2020). This favours glacier mass adjustment to negative mass balance conditions, by thinning instead of terminus retreat. As a result, the surface area change and terminus moraine position are poor indicators of glacier change for a debris-covered glacier. For example, in the Mont Blanc massif, the mostly clean-ice Mer de Glace retreated 2400 m since the 1820s LIA maximum while, over the same period, the debris-covered Miage Glacier retreated only 300 m (Deline 2005). Furthermore, this pattern of surface lowering ultimately causes a reduction in the downglacier surface slope, which reduces the driving stress. This causes progressive stagnation (Bolch et al. 2008b; Quincey et al. 2009) unless increased water availability induces widespread basal sliding of the glacier tongue (Pieczonka et al. 2018). Low-angled and stagnating glacier tongues featuring hummocky relief with large terminal moraines means that glacier meltwater cannot be efficiently evacuated through or from the glacier system. The glacier's hydrological network thus transitions from a moderately efficient, linked system to a discontinuous and inefficient network (Benn et al. 2017), with consequences for glacier lake formation and associated hazards (Benn et al. 2012).

Many debris-covered glaciers have developed proglacial lakes over the past several decades (Fig. 3) (Basnett et al. 2013; Racoviteanu et al. 2015; Nie et al. 2017; Shukla et al. 2018). Patterns of thinning and stagnation associated with many debris-covered glaciers suggest the further development of numerous additional glacial lakes is likely to continue over the next few decades (Quincey et al. 2007). Large lake formation and increased hazard potential is commonly associated with climate change and glacier recession (e.g. Zheng et al. 2021), but other analyses suggest no clear link to climate change (Harrison et al. 2018), so the subject remains controversial. A first step towards estimating the consequences of lake growth for hazards is to understand where new lakes will emerge, how large they will become, and how lake levels will change in relation to the surrounding landscape elements. It is important to recognize that lake expansion by itself is not the main criterion that renders a lake hazardous, and that lake elevation changes may be even more important than areal changes.

The current fundamental theory of debris cover evolution and lake development is heavily based on a few well-documented examples, notably in the Everest area of Nepal (Benn et al. 2012). During a period of sustained negative mass balance, a debris-covered glacier with large impounding moraines in this region (Fig. 4, ELA1) is expected to first undergo an upward expansion of the debris cover in response to a rise in the equilibrium line altitude (ELA), following which it will undergo a period of downwasting, stagnation and supraglacial pond formation (Fig. 4, ELA2) before fully stagnating at the terminus and forming a terminal lake into which the glacier terminus calves (Fig. 4, ELA3). Each of these stages is linked to several factors, notably specific mass balance and hydrological conditions. In the first stage, the rate of debris cover expansion is not solely related to the rise in ELA, but also conditioned by pre-existing debris content and changing supply rates. In the second stage, surface downwasting of the hummocky surface, coupled with inefficient meltwater evacuation leads to storage of water in perched lakes. It has been established by field studies and remote sensing techniques that the formation of supraglacial lakes is coupled to sustained negative glacier mass balance, substantial historical surface lowering and glacier stagnation towards the snout. Supraglacial lakes tend to occur primarily where the surface slope is less than 2° (Reynolds 2000; Quincey et al. 2007; Sakai and Fujita 2010; Linsbauer et al. 2016; Pandit and Ramsankaran 2020). The third stage (Fig. 4, ELA3) is marked by the coalescence of supraglacial lakes to form a large ice-contact lake at the englacial water table; during this regime, glacier mass losses at the terminus are strongly governed by calving and water-driven ablation processes within this ice-contact lake. This process has been documented in detail for two sites in the Nepal Himalaya: Imja Tsho (Watanabe et al. 1995, 2009) and Tsho Rolpa (Reynolds 1999; Sakai et al. 2000a). Numerous proto-lake systems have been identified at other glacier termini using remote sensing, e.g. Ngozumpa Glacier (Thompson et al. 2012). While the processes controlling pond expansion are well studied (e.g. Mertes et al. 2017), the rates of surface pond expansion and coalescence are not well understood and can change over time (e.g. Thompson et al. 2016). Further studies to determine if existing lakes contain buried subaqueous ice would be helpful in constraining lake deepening and basin volume over time. Such studies can be based on comparison of contemporary glacier lake bathymetry with historical ice thickness, in conjunction with studies of the sedimentation rates within lakes.

A further understanding of the glacier overdeepenings and slope conditions that may favour the formation of lakes and therefore may impose some controls on maximum lake volume requires accurate knowledge of debris-covered ice thickness. Consensus estimates of global ice thickness (e.g. Farinotti et al. 2017, 2021) developed within the framework of the Working Group on Glacier Ice Thickness Estimation of the International Association of Cryospheric Sciences (ITMIX) (http://cryosphericsciences.org/) are a valuable addition, but their appropriateness for debris-covered glaciers is unclear.

Having described the debris-covered glacier landsystem and its key components, in this section we turn our attention to the concept of hazard associated with these glaciers. We focus on glaciers that have receded from their terminal moraines and where moraine-dammed glacio-terminal lakes are created as they recede, because this is the projected end-member state for many debris-covered glaciers under future climate warming. If these moraine-dammed lakes drain rapidly because of dam failure or over-topping, a GLOF can occur, with potentially damaging consequences for downstream populations and infrastructure. Developing a robust framework for describing and assessing the potential glacial hazards associated with moraine-dammed lakes is therefore an important and societally relevant issue.

First of all, any discussion of hazard assessment requires a clear understanding of the terminology used. A common area of confusion is over usage of terms such as ‘hazard’, ‘risk’ and ‘vulnerability’, or ‘hazard assessment’ and ‘risk management’. Hazards are defined as potentially damaging physical events or phenomena which may cause the loss of life or injury, property damage, social and economic disruption, or environmental degradation. Vulnerability refers to a set of conditions and processes resulting from physical, social, economic and environmental factors, which increase the susceptibility of a community to the impact of hazard. Risk implies the probability of harmful consequences or expected loss (of lives, people injured, property, livelihoods, economic activity disrupted or environment damaged) resulting from actions between natural or human-induced hazards and vulnerable/capable conditions. ‘Resilience’ is defined as the capacity of a system, community or society to resist or to change in order that it may obtain an acceptable level in functioning and structure, and ‘capacity’ as the way people and organizations use existing resources to achieve various beneficial ends during unusual, abnormal and adverse conditions of a disaster event or process. Conventionally, ‘risk’ is expressed as risk = hazards × vulnerability/capacity’ (United Nations 2002); thus, ‘risk management’ implies ways in which the hazard might be mitigated as well as increasing the resilience of an affected community. On the other hand, ‘hazard assessment’ focuses on the initial physical processes involved with the hazardous situation, such as the triggering and development of a breach in a moraine dam. In the following section we specifically focus on ways to assess the hazard potential of a glacier lake system.

When assessing the hazard of a glacial lake, besides a solid understanding of the glacial landsystem as detailed earlier, it is fundamental to understand the components of a glacial lake system at the catchment scale and how each of those components behaves in response to the triggering of a GLOF. To fully assess the potential of a glacial lake hazard, it is there important to evaluate lake-specific factors, processes and dynamics of lakes at different stages of glacier lake development in relation to potential triggers in the surrounding landscape. As noted earlier, just because a glacial lake may contain a large volume of water, this does not necessarily make it inherently hazardous. Other factors within the glacier ‘system’ such as the range of landforms, possible mass movement processes, and other influencing factors from within the glacier system environment and the surrounding mountain flanks from the top of the headwall to the lowest terminal moraine dam need to be thoroughly evaluated. This requires a holistic overview of the glacial system in order to identify key components that, if present, may trigger one or more processes that might lead to the formation of a GLOF. The goal is to identify key components that may trigger one or more potentially cascading processes that might lead to a GLOF event.

When assessing the glacier lake hazard potential, two important issues exist: (a) how to assess the hazard across a region in a consistent and meaningful way and (b) how to rank them in terms of the severity of the hazard (Reynolds 2014) in a systematic, quantitative manner. In quantifying GLOF hazard, remote sensing techniques have been used to develop Tier 1 (first-pass) assessments over large areas (tens to hundreds of square kilometres) (Kääb et al. 2005; Quincey et al. 2005). Such first-pass automated assessment schemes have been developed for the Tibetan Plateau (Allen et al. 2019), the Indian Himalayas (Dubey and Goyal 2020), the Andes (Frey et al. 2018; Kougkoulos et al. 2018) and the European Alps (Huggel et al. 2004). For Tier 2 local assessments of specific glacial lake systems, very high-resolution imagery (< 1 m spatial resolution) and associated DEMs have been used for small areas (e.g. 25–100 km² or more); drones have been used to produce very high-resolution imagery and photogrammetry for this purpose (Westoby et al. 2012; Fugazza et al. 2018; Wilson et al. 2019). The UAV and terrestrial structure from motion (SfM) photogrammetry techniques bridge the gap between the difficult field campaigns and the coarse satellite data, and emerged in the last decade as a promising opportunity for estimating hazard potential and hazard management strategies. The results from such analyses can be used to complement or support field campaigns that include, for example, detailed geomorphological, geophysical, topographical and engineering geological surveying and mapping (Hambrey et al. 2008). However, a better integration of Tier 1 and Tier 2 assessments is currently needed to assess hazard potentials at multi-scales.

Despite technical guidelines on the assessment of glacier and permafrost hazards in mountain regions published by the International Association of Cryospheric Sciences and International Permafrost Association Standing Group on Glacier and Permafrost Hazards (Allen et al. 2017), a standard lake hazard assessment scheme does not exist. The existing glacial lake ranking schemes (e.g. Quincey et al. 2007; Bolch et al. 2008a; Wang et al. 2011; Worni et al. 2013; Iribarren Anacona et al. 2014; Rounce et al. 2016; Aggarwal et al. 2017; Kougkoulos et al. 2018; Dubey and Goyal 2020; Pandit and Ramsankaran 2020) all differ based on the parameters used, the weight assigned to each parameter and the source of the data used (field/remote sensing/a priori knowledge) (Emmer and Vilímek 2013; Rounce et al. 2016). Furthermore, existing schemes do not always parameterize key GLOF processes.

Consequently, there is interest in developing a standard, objective unified ranking scheme on the basis of new remote sensing data. Such a scheme would ideally be decision-based, constructed on multi-criteria and using state-of-the-art techniques such as machine learning. Furthermore, such a scheme needs to quantify both observable conditioning and triggering factors related to GLOF formation rather than on subjective criteria or derived parameters, such as lake volume. There are four threshold factors that can be used to categorize any given glacial lake system, but which on their own do not designate the existence of any hazard (Table 4). These have been designed to be used especially as a Tier 1 preliminary screening/hazard ranking tool. However, for a hazard to exist, there must be potential for a trigger event to occur that can lead to a possible GLOF. The key factors affecting the likelihood of a GLOF include: (a) minimal moraine freeboard above the lake level with narrow dam width, rendering the dam vulnerable to overtopping; (b) evidence of avalanches from valley sides and backwall, and/or from hanging glaciers directly into the lake that might induce either a seiche or avalanche push wave; (c) evidence of seepage and/or piping through the moraine dam; (d) evidence of degradation of an ice core within the terminal moraine dam that might cause progressive collapse (RGSL 2003).

For both threshold and trigger parameters, scale factors can be used to weigh how important or significant any factor is. In general, to derive a hazard score, each threshold parameter (Table 4) is ranked using one value from each of the weighting columns and summed; similarly, a trigger parameter score (Table 4) is similarly derived. This enables a hazard score to be derived using the weightings for both threshold and trigger parameters. For example, to account for relative differences in lake volume, which is often a derived value based on lake area, measured areas are used. The two scores for the threshold and trigger parameters are used as (x, y) coordinates to plot on a hazard ranking graph (see Fig. 5).

In addition to glacier and lake processes occurring upstream, a full hazard assessment scheme needs to include the impact on the populations downstream, and a full socio-economic assessment (Carey et al. 2012). In recent years there have been advances in both GLOF modelling and integration of models with robust assessments of glacial hazards and their societal impacts. For example, losses incurred from hydropower schemes following GLOFs has led the international hydropower sector to build greater resilience to climate change impacts (RGSL 2015; Reynolds 2018). The complexity of such damaging events in triggering mechanisms and in the changing processes as they propagate downstream calls for catchment-wide assessments of such geohazards. The challenged is that modelling the GLOF impact downstream requires sophisticated flow modelling which implies a number of assumptions about the characteristics of the material, lake volume, peak discharge, sediment load, channel roughness, etc. (Fig. 6) which are difficult to measure (Iribarren Anacona et al. 2015). In the last decade, multiple studies have tested and deployed a variety of modelling tools to perform numerical simulations of GLOFs downstream and to simulate different type of flows (Westoby et al. 2015). Numerical modelling approaches that coupled glacial lake impact, dam breach and flood processes are reviewed in Worni et al. (2014). One of the shortcomings of current models is that flow characteristics are complex and commonly develop as a cascade of physical processes as the flow propagates downstream. This poses the need for modelling multi-phase GLOF process cascades (e.g. Schneider et al. 2014; Worni et al. 2014; Mergili et al. 2020). Furthermore, the extreme flows are difficult to measure for calibration purposes, which entails a large degree of uncertainty (Worni et al. 2014).

Given the large uncertainties associated with the GLOF process chain in terms of timing, location and intensity of triggers (Schneider et al. 2014; Allen et al. 2017), one of the key remaining challenges is how best to communicate the changing nature of hazards (and implications for GLOF model uncertainty) to communities/stakeholders. Finally, one aspect of hazard and risk assessment that is now well established in the private sector but less so in the academic world is exposure to legal responsibility and the consequential liability arising from making statements about risk that could have outcomes affecting asset values.

Even with the substantial progress on mapping of the debris cover and associated lakes, there remain significant challenges to be addressed in terms of approaching the landsystem using a holistic approach in view of developing a standard hazard raking scheme. Here we summarize the remaining limitations and gaps in our knowledge of the system.

• Mapping of debris-covered glaciers often relies on expert knowledge, which is subjective and often subject to disagreement, especially when independent, ground truth data are not available. There is no standardized mapping for debris-covered glaciers, and existing methods are generally ‘semi-automated’ because they involve some manual correction (Racoviteanu et al. 2009; Herreid and Pellicciotti 2020). While providing important information, available global or automated methods are only suited to mapping debris within assumed glacier extents (Scherler et al. 2018). Some of the new methods for debris cover mapping have the potential to automate the mapping process of debris-covered glacier tongues, but these need further and robust testing. There remain significant gaps in high-quality debris cover outlines in some glacierized regions, and retrieval of most key debris properties from remote sensing is at a very early stage;

• Within the debris-covered glacier landsystem, debris sourcing and evacuation is a key gap in knowledge; understanding of erosion rates varies regionally, but most erosion rates are millennial-scale values. Thus, further advances need to be made to assess contemporary and recent erosion and debris supply rates within debris-covered glacier landsystems. This might include the use of fine-resolution imagery, derivation of debris supply from avalanche cones, and other creative analysis or numerical modelling approaches (e.g. Banerjee and Wani 2018; Scherler and Egholm 2020). It is key to study both singular, large debris supply events (e.g. Berthier and Brun 2019), and smaller events of debris supply using holistic efforts. Furthermore, debris flux out of the system is so far only crudely represented despite being a key property governing glacier development over time, and it may be valuable to identify the determinants of whether or not a glacier forms a large latero-terminal moraine;

• The issue of scale in remote sensing remains an important challenge. While significant progress has been made in monitoring debris cover surface properties using both field methods and remote sensing, these are often applied at different spatial scales (local to regional), making it difficult to transfer the observations from one scale to the other. Detailed field studies offer insight into specific processes (e.g. ice ablation), but they are often site-specific; remote sensing studies, on the other hand, can be applied at multi-scales, but face limitations due to spatial and temporal resolution, i.e. lack of high-resolution thermal data or surface velocities.

• The implications of increased debris supply remain unclear, for example:

  • o How do glacier thermal and dynamic regimes respond to the increased debris resulting from glacier thinning and upwards migration of debris? What are the implications of the increased debris for basal sliding, glacier thinning and stagnation, ice thickness and deformation?

  • o What are the typical glacial structures associated with increased debris supply, and what are the consequences of these structures downglacier? How will these glaciers respond to changes in terms of hydrology, ice deformation and surface debris?

  • o What are the consequences of different levels of debris sequestration on glacial landscape evolution and geomorphology? How will subglacial erosion, moraine building, lake development and sedimentation change through time?

  • o How do permafrost and debris-covered/rock glaciers interact with the glacier(s) and debris/mass fluxes through the system?

  • o How does the increased debris supply influence the formation of proglacial and supraglacial lakes impounded by lateral or terminal moraines and by supra-glacial debris deposits? How does this affect the probability of glacier-related hazards (lake outbursts floods and debris slides) under a transient climate?

• The rates at which the system transitions to different states are unclear, and proxies to project them forward in time are needed to accommodate them in glacier model projections. For example, more work is need in order to make projections of rated of debris cover expansion/thickening from an initial state of unknown debris load in the system, and with uncertain debris inputs/outputs. There is a need for improved models (i) to reasonably account for the meltwater production role of transient features such as ice cliffs, (ii) to project glacier downwasting and surface slope evolution to the threshold of supraglacial pond development and (iii) to parameterize rates of supraglacial pond expansion to allow the likely timing of pond coalescence to be estimated.

• The development and testing of a standardized, integrated lake hazard ranking scheme remains a challenge. This requires better parameterization of key GLOF processes in the glacier lake system, and the ability to capture the multi-phase GLOF process cascade.

In order to address these questions and to complete the picture of processes associated with the debris-cover glacier system, improved datasets are still needed, i.e. meteorological data from weather stations installed at high altitudes (e.g. Matthews et al. 2020), monitoring of the rates and controls on rockwall debris supply, gauging of water and sediment amounts discharged within turbid glacier streams (e.g. Heckmann et al. 2016), spatially distributed measurements of ice thickness from new technologies such as airborne ice-sounding radar suited to debris-covered glaciers (e.g. Pritchard et al. 2020) and debris thickness distributions (e.g. Nicholson et al. 2018) with which to optimized models of debris thickness, studies of the sediment discharge to the glacier foreland high-resolution regular-repeat imagery for selected debris-covered glacier landsystems with differing characteristics. There is a need for more studies that link climate, mass balance, ice dynamics, topographic evolution and hydrology to quantify how hazard potential emerges from interactions between these processes. Finally, there is a need to bridge spatial scales both in terms of connecting processes and resolving them in remotely sensed data. For example, some satellite imagery cannot resolve metre-scale features, even though features such as ice cliffs at this scale may collectively be significant to runoff generation, and the melt processes operating locally at ice cliffs need to be integrated into a glacier scale representation of the ablation regime (Ferguson and Vieli 2021).

This paper stems from a workshop supported by the Geological Society of London in 2019, that brought together researchers with a shared interest in debris-covered glaciers and related hazards with a broad range of experience and activities in approaching these landscape systems. As such, this perspectives paper draws together key insights, state-of-the-art and consensus research priorities from the exchanges fostered by the workshop. While the key state of the knowledge has been described in the preceding sections, to conclude we wish to draw out a small number of key messages.

When considering debris-covered glaciers, we argue that it is vital to adopt a landsystems approach that includes the flux of both solid and liquid water and sediment within catchments as well as estimates of how these processes influence and are influenced by glacier behaviour over time. Despite key developments and advances in the use of satellite remote sensing to estimate these processes, there remain gaps in the validation of these tools using field-based measurements as these remain scarce. There remains much work to be done to develop robust tools to upscale the knowledge gained from small process studies to a landsystems scale so that it can be integrated in satellite monitoring and numerical models of larger spatio-temporal scale glacier and landscape development.

Debris-covered glaciers are projected to increase in number proportionately as mountain glaciers diminish, but the specific trajectories of glacier development are elusive due to the complex coupling of non-stationary processes and feedbacks within the debris-covered glacier system. Critically, some glaciers form large impounding latero-terminal moraines that drive local hydrological processes and which have implications for glacier hazards, while others do not, and we lack a clear method of discriminating which pathway a given glacier or glacierized region will follow.

Finally, we suggest that consideration of cascading hazards within the wider landsystem is critical for developing meaningful glacial lake hazard assessment. There is a need to address this issue due to communication failures in the past, so a better interaction between the debris-covered glacier community and the geomorphological and climate science communities is needed for this perspective framework to be successful.

We thank the Geological Society for their support in hosting the ‘Debris-covered Glaciers and Related Lakes: Understanding the Challenges’ workshop held in London in September 2019. We thank all the speakers in the workshop: T. Bolch, J. Carrivick, E. Miles, D. Quincey, J. Reynolds, M. Westoby and S. Watson for their insights during the workshop discussions, and all workshop participants for their efforts and enthusiasm.

AER: conceptualization (lead), investigation (equal), writing – original draft (lead), writing – review & editing (equal); LN: conceptualization (equal), writing – original draft (equal), writing – review & editing (supporting); NFG: project administration (lead), supervision (lead), writing – review & editing (supporting); EM: conceptualization (supporting), methodology (supporting), writing – original draft (equal), writing – review & editing (equal); SH: writing – review & editing (supporting); JMR: writing – review & editing (supporting)

AR's research was supported by the European Union's Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No 663830 and SCoRE Cymru. JR was supported by Reynolds International Ltd. LN was supported by Austrian Science Fund (FWF) Grant P28521.

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Data sharing is not applicable to this article as no datasets were generated or analysed during the current study.

Scientific editing by Philip Hughes