Cooling glaciers in a warming climate since the Little Ice Age at Qaanaaq, northwest Kalaallit Nunaat (Greenland)

The centennial response of land‐terminating glaciers in Greenland to climate change is largely unknown. Yet, such information is important to understand ongoing changes and for projecting the future evolution of Arctic subpolar glaciers, meltwater runoff, and sediment fluxes. This paper analyses the topography, geomorphology, and sedimentology of prominent moraine ridges and the proglacial areas of ice cap outlet glaciers on the Qaanaaq peninsula (Piulip Nunaa). We determine geometric changes of glaciers since the neoglacial maximum; the Little Ice Age (LIA), and we compare glacier behaviour during the LIA with that of the present day. There has been very little change in the rate of volume loss of each outlet glacier since the LIA compared with the rate between 2000 and 2019. However, the percentage of each glacier that is likely composed of cold‐based ice has increased since the LIA, typically by 20%. The LIA moraines comprise subrounded, striated, and faceted clasts that evidence subglacial transport, and outwash plains, flutes, kames, and eskers that evidence subglacial motion and meltwater within temperate ice. Contrastingly, contemporary ice margins and their convex ice surfaces comprise pronounced primary foliation, ephemeral supraglacial drainage, sediment drapes from thrust plane fractures, and an absence of open crevasses and moulins. These calculations and observations together lead us to interpret that these outlet glaciers have transitioned towards an increasingly cold‐based thermal regime despite a warming regional climate. Thermal regime transitions control glacier dynamics and therefore should be incorporated into glacier evolution models, especially where polythermal glaciers prevail and where climate is changing rapidly.


| INTRODUCTION
Mountain glaciers and ice caps account for between 30% and 60% of contributions to global sea level rise (e.g., Bamber et al., 2018;Gregory et al., 2013;Meier et al., 2007;Radi c & Hock, 2011).They will contribute between $80 mm (RCP2.6)and $160 mm (RCP8.5) to global sea level rise by the year 2100 (Marzeion et al., 2020).Because the Arctic is the one of the fastest-warming parts of our world (IPCC, 2021), there is an urgent need to understand the evolution of Arctic mountain glaciers and ice caps, which account for 39% of the total number and 62% of the total area of glaciers and ice caps worldwide (RGI Consortium, 2017;chp. 6).However, Arctic glacier evolution is complex not least because (i) most Arctic glaciers have a polythermal regime and (ii) the distribution of cold and warm parts of a polythermal glacier evolves seasonally and over centennial timescales.Polythermal refers to glaciers composed of parts with permanently frozen beds and parts with 'warm' or 'temperate' ice.
Temperate ice has basal meltwater present either due to pressure melting or to hydrological pathways at the glacier bed.
Glacier thermal regime is not represented in global glacier evolution models.That omission is a problem because thermal regime controls subglacial hydrology and so is fundamental in influencing spatiotemporal ice dynamics at both short (daily and seasonal) and long temporal (decadal to centennial) scales (e.g., Naegeli et al., 2014;Rabus & Echelmeyer, 1998;Rippin et al., 2003).It is imperative to understand the changing distribution of cold and temperate ice to understand hydrological and ice dynamical processes and long-term changes, which are especially poorly constrained (Irvine-Fynn et al., 2011).Specifically, temperate glaciers have liquid water at the ice-bed interface, which enables basal sliding, high ice velocities (hundreds to thousands of mÁyear À1 ), and high rates of geomorphological work.Cold-based glaciers are frozen to the bed and therefore move slowly (tens of mÁyear À1 ) and have low sediment transport rates.
Despite this imperative, studies of the thermal structure and evolution of subpolar or polythermal type mountain glaciers are few and far between.Almost all are either on a single glacier (Glasser & Hambrey, 2001, excepted for example), and most are from Svalbard either focussing on interpreting ground penetrating radar (GPR) data (e.g., Baelum & Benn, 2011;Hodgkins et al., 1999;Irvine-Fynn et al., 2006;Karušs et al., 2022;Moore et al., 1999;Stuart et al., 2003;Temminghoff et al., 2019) or on combining glaciological, topographical, and geomorphological-sedimentological analyses (see tab. 1 and citations within Evans et al., 2022).Radio-echo sounding has been used to interpret the thermal regime of the largest glaciers in Iceland (e.g., Bjornsson et al., 1996), and GPR has been used on single glaciers in Iceland and Arctic Scandinavia (e.g., Gusmeroli et al., 2010;Pettersson et al., 2003;Reinardy et al., 2019), Arctic North America (e.g., Bingham et al., 2005;Brown et al., 2009), and in the Southern Hemisphere on King George Island (Blindow et al., 2010).Borehole temperature measurements and other in situ hydrological observations (e.g., Skidmore & Sharp, 1999;Sobota, 2009) used to suggest thermal regime are particularly rare.A multidisciplinary approach including geomorphological-sedimentological analyses has been conducted on James Ross Island off the Antarctic Peninsula (Carrivick et al., 2012).There remains (since Hambrey & Glasser, 2012, see their fig.1) a dearth of information on the thermal regime of glaciers in Greenland.Greenland glaciers and ice caps occupy 20% of all Arctic glacier area, and they represent 9% of the number and 13% of the area of all glaciers and ice caps worldwide, second only to the Antarctic Peninsula (RGI Consortium, 2017;chp. 6).
The aim of this study is to synthesise the character and behaviour of a sample of Arctic glaciers in northwest Greenland during the Little Ice Age (LIA) and at present.We achieve this aim by using (i) spatial analyses to quantify the geometric and physical evolution of glaciers and (ii) high-resolution topography, geomorphological, and sedimentological observations to infer process-form relationships and hence thermal regime relating to both past and present glacier extents.

| STUDY SITE
Qaanaaq ice cap covers an area of $260 km 2 and is located on a peninsula (Piulip Nunaa) in the Thule region of northwest Kalaallit Nunaat (Greenland) (Figure 1).Geologically, the Piulip Nunaa is composed of Mesoproterozoic sedimentary rocks of several subtypes with grainsizes from siltstones to conglomerates.Following the Last Glacial Maximum, ice retreat initiated at $12 ka in the fjords, whereas deglaciation of the Taserssuit Valley (Figure 1) was $4000 years later (Søndergaard et al., 2019).The Qaanaaq ice cap outlet glaciers can be assumed to have advanced in the Neoglacial, if it is supposed that they behave in accordance with glaciers across Greenland (Kjaer et al., 2022).Qaanaaq ice cap was smaller than at present until 900 years ago after which outlet glaciers expanded towards their LIA maximum (Søndergaard et al., 2019).On average and across Greenland, the LIA terminated in year 1900 (Kjaer et al., 2022).
The summit of the present-day Qaanaaq ice cap is $1110 m asl., and 17 outlet glaciers extend down to <100 m asl., the largest of which include Fan Glacier and Scarlet Heart Glacier (Figure 1).Direct glaciological investigations have been undertaken at Qaanaaq glacier (Figure 1) over 10 days in July 2012, in which ablation stake surveys quantified (summer) melt rates of $40 mmÁday À1 and surface velocities of $50 mmÁday À1 but with significant (50% to 150% of the mean) diurnal variability, and GPR data revealed ice thickness typically of $120 m (Sugiyama et al., 2014).

| LIA reconstruction
Prominent moraine ridges at all major land-terminating outlet glaciers of the Qaanaaq ice cap were identified on a hillshaded image of the 2-m resolution ArcticDEM (Porter et al., 2018).This is a mosaic Digital Elevation Model (DEM), and images used to construct it in the Qaanaaq region mostly date from 2015.The moraine ridges have not been geochronologically dated but are ascribed to the LIA based on the following: (i) their position within/below a neoglacial trimline (Figure S1); (ii) their prominence, which is larger and distinct than icecored moraines abutting the contemporary ice margin (Figure S1); and (iii) their composition, which contrasts with vegetated surfaces beyond the proglacial area and till plains and outwash plains that are situated between them and the contemporary ice margin (Kjaer et al., 2022).Following the workflow of Carrivick, Boston, et al. (2019), Carrivick et al. (2020Carrivick et al. ( , 2022)), and Lee et al. (2021), these moraine crests were manually digitised by reshaping the RGI_v6.0outlines and extended along any trimlines.In this region, the RGI_v6.0outlines, which are available from Global Land Ice Measurements from Space (GLIMS) (Raup et al., 2007), pertain to year 2000.We analysed our LIA outlines and the RGI_v6.0glacier outlines to derive an equilibrium line altitude and hence an ablation area, using an accumulation area balance ratio of 2.24, as generally applicable to Arctic glaciers (Rea, 2009).Our LIA ablation area outlines were converted to points to extract elevations from the ArcticDEM, and a surface interpolated between those points using a natural neighbour algorithm produced a simple three-dimensional reconstruction of the ice surface during the LIA.Differencing the reconstructed LIA ice surface with the ArcticDEM produced an elevation change map and hence an estimate of glacier volume change.Our change estimates are a minimum due to excluding higher parts of the outlet glaciers and ice cap in the accumulation areas.

| Ice thickness, bedrock topography, and thermal regime
Previous studies on the thermal regime of Arctic glaciers have relied on radio-echo sounding or GPR data for ice thickness, and a single longitudinal profile of GPR-derived ice thickness is available from Sugiyama et al. (2014).However, we herein exploit the ice thickness model dataset of Millan et al. (2022) because their velocity inversion model is especially suitable for application to ice caps and low-angle glaciers, such as are prevalent at Qaanaaq, and because it covers the entire ice cap.We use the ice thickness to derive subglacial bedrock topography by subtracting it from the surface ArcticDEM.Ice thickness during the LIA maximum, which we take as year 1900 following Kjaer et al. (2022), was estimated as the difference between our LIA ice surface and the bedrock topography.
Earlier research has suggested that a threshold ice thickness of $80 to 100 m is required for temperate ice to persist year-round in Arctic or subpolar glaciers (Murray et al., 2000).Although some numerical modelling has shown that subglacial meltwater pathways and a lag effect of past thermal conditions can create discrepancy between theoretical and actual basal conditions (i.e., pressure melting point being reached or not; Wohlleben et al., 2009), field data on ice thickness and glacier thermal regime have shown a suggested ice thickness threshold of $90 m to be a suitable simple rule for glacierwide patterns of thermal regime (e.g., Karušs et al., 2022).Therefore, whilst mindful of its simplicity and regarding its utility for quickly identifying spatio-temporal patterns, we herein apply that threshold to our reconstructed ice thickness from the LIA to the present to infer the spatial distribution of ice is likely to be cold-based.
F I G U R E 1 Location of Qaanaaq ice cap and names of outlet glaciers focussed on in this study.Glaciers (year 2000 from the RGI_v6.0)are from GLIMS (Raup et al., 2007), ice thickness is from Millan et al. (2022), and the background is a hillshaded image of the ArcticDEM mosaic (Porter et al., 2018).[Color figure can be viewed at wileyonlinelibrary.com]

| UAV-based high-resolution imagery and topography
The ArcticDEM (Porter et al., 2018) is exceptional in its spatiotemporal coverage, and the 2-m horizontal resolution is appropriate for many terrain analyses.However, it is too coarse to resolve microtopography, surface texture, and surface composition details.Furthermore, coverage of high-resolution satellite imagery, such as PlanetScope (https://www.planet.com),does not extend this far north.Therefore, in August 2022, we conducted our own UAV-based Structure-from-Motion photogrammetry field surveys at four outlet glaciers of the Qaanaaq ice cap, selected for their distribution around the ice cap.
The purpose of these UAV surveys was to provide higher resolution imagery and topographic data products to assist in the identification of surface features and processes.Local ground control was unavailable during the field campaign; thus, the surveys exhibit uncertainties in the absolute point locations and would be inappropriate data sources on which to base extensive quantitative topographic analyses.Nevertheless, the relative point positions and surface shape and texture within the survey area are consistent and can be used to supplement field mapping, albeit with greater likelihood of systematic errors arising from 'weak' georeferencing (i.e., large uncertainties in model scale, translation, and rotation) (James et al., 2017).
We derived DEMs of 0.2-m horizontal resolution alongside orthomosaic images (0.05-to 0.2-m resolution) for each outlet glacier survey, using imagery acquired from a DJI Mavic 2 Pro Uncrewed Aerial Vehicle (UAV).Details of the UAV and camera sensor are provided in Table S1.The UAV was flown in an approximate grid pattern above the survey area, ensuring sufficient overlap and sidelap between images for photogrammetric reconstruction.The surveys were 'directly georeferenced' (in WGS84) using the relatively lowquality (metre-scale accuracy) on-board UAV Global Positioning System (GPS) and processed using Agisoft Photoscan Professional 1.4.0.
Camera alignment errors were typically between 4 and 8 m.The resulting orthomosaics and DEMs were manually aligned with ArcticDEM data to correct for the poorly constrained georeferencing.Details of each survey and errors are provided in Table S2.Our Structure-from-Motion-derived DEMs and orthomosaics are mapped at their full extent in Supporting Information S1.
These datasets were sufficient to enable spatial analyses, to derive microtopography, which we obtained as elevation deviations from a local 10-m mean, and to interpret geomorphology, surface texture, and glacier surface structure at each outlet glacier.

| Geometry and morphology of outlet glaciers
Fragmentation of former tributaries to individual glaciers (e.g., at Sydgletscher) means there were 17 glaciers during the LIA compared with the 19 glaciers present now (Figure 1).The area change from 322 km 2 during the LIA to 301 km 2 in $2015 equates to a 6.5% area loss.Total surface lowering since the LIA is typically <25 m, but a few parts of some glaciers have been up to 120 m, equating to a rate of typically <1.0 mÁyear À1 (Figure 2a).There has therefore been no change in the rate of surface lowering when comparing between the LIA to 2015 (Figure 2a) and 2000 to 2019 (Figure 2b).
Volume loss across the entire Qaanaaq ice cap between the LIA maximum and 2015 has been at least 1.8 km 3 (±20%) at a rate of at least 0.015 km 3 year À1 , and between 2000 and 2019 was 0.185 km 3 (±50%) at a rate of 0.009 km 3 year À1 .Whilst mindful of uncertainty, comparison of these rates between the time periods is suggestive of a decrease in the rate of volume loss.The mean equilibrium line altitude across Qaanaaq has risen by 57 m from 682 m asl.during the LIA to 739 m asl.during approximately year 2015.Comparing modelled ice thickness within glacier ablation areas during the LIA (Figure 2c) and between 2000 and 2019 (Figure 2d) reveals the expansion of probable cold-based parts of glaciers, typically by $20% more coverage for all but the largest glaciers (e.g., Fan Glacier only 20% to 28%).

| Proglacial geomorphology and sediments
Sydgletscher, the unnamed glacier, Fan Glacier, and Scarlet Heart Glacier (Figure 1) all have two sets of major moraine ridges: an outermost single-crested sharp ridge and an inner ridge that in contrast has a more subdued relief but is larger/more voluminous overall (Figures 3a,b,4a,b,5a,b,and 6a,b).The subdued relief of the inner moraine could reflect having been overridden by glacier ice (during the LIA) making it an older neoglacial moraine (cf.nearby on east Baffin Island; Young et al., 2015).Sydgletscher proglacial area is dominated by ice-cored moraines, stacked kame-moraine ridges, and braided outwash plains (sandur) (Figure 3a).Most uniquely, eskers and kames are present in the Sydgletscher proglacial area (Figures 3b, S5, and S6).In contrast, the distal part of the proglacial area of the unnamed glacier is narrow and steep and is characterised by the imprint of both erosional processes in the form of scoured bedrock and minor meltwater channels on the east-facing hillslope and by depositional processes in the form of several single-crested moraine ridges especially on the west-facing hillslope (Figure 4a).The proglacial area at Fan Glacier has an outwash plain (sandur) beyond the LIA limit (Figure 5a), a till plain with crude fluting within that LIA moraine (Figure S7) limit (Figure 5b), and, with the exception of two large gorges within the LIA lateral moraine, progressively less (coverage of) evidence of meltwater with proximity to the contemporary ice margin (Figure 5c).Similarly, the proglacial area of Scarlet Heart Glacier is also dominated by a till plain (Figure 6a).

| LIA moraine ridges and sediments
Moraine ridges around the outlet glaciers of Qaanaaq ice cap vary considerably in morphology and composition, both intersite and intrasite.Fan Glacier has (LIA) lateral moraine complexes that are typically >250 m wide and formed of multiple hummocky pseudo-parallel ridges each $5 m high (Figure 7a) and composed of matrix-supported diamict containing subangular cobbles and small boulders (Figure S7).
Away from the hillslopes, the former terminus of Fan Glacier is marked by a single low relief ($5 m) ridge of clast-supported diamict containing subangular boulders (Figure 7e).Fan glacier is also unusual for two relatively deep and narrow gorges running parallel to the western lateral moraine ridges (Figure 5c).Given that these gorges cut through the moraine, they must have been formed more recently than the moraine building phase.The gorges are incised into siltstones and are strewn with well-rounded imbricated sandstone boulders, that is, erratic boulders that have been fluvially transported (Figure S2).
Scarlet Heart Glacier has two sets of arcuate moraine ridges, which is not unusual in N. Greenland and which most likely reflects a neoglacial ice margin advance and a subsequent LIA advance (Kjaer et al., 2022).Both ridges are small in amplitude ($5 m) but nevertheless distinct from surrounding till surfaces by a near-unbroken arcuate morphology as well as in clast size (>1-m diameter boulders) and clast shape (subangular).This almost entirely intact (i.e., not dissected or reworked by outwash) till plain composed of faceted and striated clasts is essentially similar to that at Fan Glacier but contrasts with the extensive glacifluvial outwash plains at Sydgletscher and the steeper outwash fans at Qaanaaq Glacier, for example.The steepest and highest relief (50 m) moraine ridges can be found in front of a small unnamed glacier immediately to the west of Scarlet Heart Glacier (Figures 6c and   T A B L E 1 Landforms and processes that can be ascribed to distinctive parts of a polythermal glacier.4c and 7c).They are ice-cored (Figure S3), composed of subangular cobble and boulder clasts, and are experiencing large-scale degradation/collapse (Figure 7c).
LIA lateral moraine ridges at Sydgletscher are distinct in colour from the hillslopes (Figure 7d) but in their lower parts are the most degraded of any we observed; active slumps, slides, and flows dominate their present-day morphology (Figure 3c).Natural exposures due to slumps of surficial till reveal that these lateral moraines are icecored.Moraine degradation combined with the exceptionally soft and distinctively red-coloured sandstone-derived sediment means that the sediment flux from this catchment is likely very high, as we witnessed exceptionally turbid rivers and noted a large active delta and sediment plume into the fjord (Figure S4).

| Contemporary ice margins
Ice margins of Qaanaaq outlet glaciers presently have long-profile slopes that are gentle (typically <5 ) and merge without a change of slope to proglacial sediments.In contrast, lateral or cross-sectional profiles are convex, and the lateral ice margins often have an abrupt near-vertical wall up to $10 m high (e.g., Figure 7a).Scarlet Heart Glacier has an ice margin with only a few small drainage channels emanating from it; the largest is associated with the medial moraine.
However, there are abandoned subglacial channels (Figure 8a) and associated abandoned proglacial streams (Figure 6c).The abandoned channel reveals basal ice with high concentrations of matrixsupported diamict with subrounded cobbles (Figure 8a).Sydgletscher has developed (since $year 2000) an embayment behind a near-  Surface sediment types range from widely distributed aeolian dust and cryoconite deposits (cf.Matoba et al., 2020;Takeuchi et al., 2014;Uetake et al., 2010) to discrete deposits of coarser sediment, namely, geometric ridge networks (Figure 9a) and dirt cones (Figure 9b,c).The geometric ridge networks are composed of siltysand ridges, and we noticed that one of those on Sydgletscher has an odorous biogenic component.These geometric ridge sediments emanate from the glacier beds via structural cracks in the lowermost parts of glaciers, which are most likely thrust planes.The ridges are up to 2 m high, and they insulate the ice below from melting (Figure 9a).
They probably do not have a long preservation potential, perhaps only months or a few years.We also observed dirt cones in the lowermost parts of glacier ablation areas, and some cones were impressively high (14 m) (Figure 9b), and so, we consider that they are quite mature features because they also evidence the insulation of ice and the contemporary ice surface lowering is slow ($ 0.1 mÁyear À1 ) (Figure 2b).
Those dirt cones in Figure 9b contained sediment that had been transported by a supraglacial stream sourced from a medial moraine.Sediment emerging from the glacier bed onto the ice surface via thrust planes also forms sediment drapes on the ice surface at Fan Glacier (Figure 9d).
Contemporary hydrology and drainage of the Qaanaaq glaciers is mostly via surface runoff in minor (<0.5 m deep) rivulets formed from diffuse flow within firn and superimposed ice on glacier ablation areas.
However, as well as the aforementioned abandoned subglacial conduit(s) emanating from glacier termini (e.g., Figure 4b), we also observed a few abandoned englacial conduits (Figure 9e).These  et al., 2009) but probably due to a migrating piezometric surface, both seasonally and longer term (cf.Irvine-Fynn et al., 2006).The few active supraglacial streams that we found were deeply incised and extremely sinuous, and most bends were greatly undercut (Figure 9f), which we interpret to illustrate progressive thermo-mechanical fluvial erosion and thus persistence and longevity of the channel pathway through many melt seasons.

| DISCUSSION
During the LIA, the land-terminating outlet glaciers of Qaanaaq ice cap had greater proportions of their beds with temperate ice than at present.Indeed, the probable coverage of cold ice at the bed has increased by >30% for some glaciers when comparing the LIA state (Figure 2c) to between 2000 and 2019 (Figure 2d).Additionally, we observe that the composition of LIA moraine ridges of poorly sorted matrix-supported diamict containing heterolithic boulder-sized, faceted, and striated subrounded clasts evidence active subglacial transport and deposition with longitudinal compression and probable glaciotectonic thrusting (Cofaigh et al., 1999;Hambrey & Glasser, 2002) and therefore 'warm' or temperate ice during the LIA maximum at the outlet glacier termini and ice margin(s).The gorges incised through moraine ridges at Fan Glacier (Figure S2) evidence high volumes/energies of (proglacial) water, which could be ablationfed meltwater (cf. the ice-marginal channels in the Dry Valleys; Atkins & Dickinson, 2007).However, the imbricated well-rounded boulders suggest considerable hydraulic power and sustained fluvial abrasion, so we interpret these gorges to be flood-related.Additional  Sydgletscher, and Fan Glacier (Figure 2d) corresponds to the locations of where we observed geometric ridges of fine-grained basal sediment (Figure 7a), dirt cones (Figure 7b), abandoned englacial conduits (Figure 7e), and supraglacial streams (Figure 7f).Therefore, it seems likely that these landforms and sediments are the product of a surface cold ice layer (cf.Blatter & Hutter, 1991;Irvine-Fynn et al., 2011)   The pattern and spatial extent of surface lowering is also important; during the $120 years since the LIA, it has been restricted to the lowermost parts of the outlet glaciers (Figure 2a) but between 2000 and 2019 has been approximately uniform across all parts of all ablation areas (Figure 2b).The surface lowering is important because where ice becomes thin, the winter cold (air) wave can penetrate to the bed and retard pressure melting and make basal sliding unlikely.
However, the distribution of cold and temperate ice within a glacier is a result not only of the present-day conditions but also of the past.For example, temperate ice found within the interior of a polythermal glacier could have been advected from the accumulation area, which might include a lag in the thermal evolution (Rippin et al., 2011), or it could have been generated in situ (Irvine-Fynn et al., 2006, 2011).
Given these centennial and decadal scale changes, the future evolution of the Qaanaaq ice cap is therefore uncertain and that is not least due to the polythermal state of these glaciers.We realise that 'cooling glaciers' and unchanged (no significant increase) ice mass loss rates are unusual in a global context given the general acceleration of ice mass loss in temperate world regions since the LIA (e.g., Carrivick et al., 2020Carrivick et al., , 2022;;Davies & Glasser, 2012;Lee et al., 2021) and for NE Greenland (Carrivick, Heckmann, et al., 2019).However, although the Arctic is warming extremely rapidly, glacier responses across Greenland exhibit great spatio-temporal variability (e.g., Hugonnet et al., 2021;Khan et al., 2022), for example, when comparing inlandcontinental versus coastal-maritime environments (e.g., Osman et al., 2021) and with a suite of glaciological controls such as hypsometry and terminus environment (e.g., see citations within Kjaer et al., 2022).In a recent study from the high Arctic, Kochtitzky et al.
(2022) also found no change in rate of mass loss of a single glacier on Ellesmere Island and attributed this to its cold-based condition.Unfortunately, a lack of volumetric analysis of glacier changes since the LIA hinders other comparisons.

| SUMMARY AND CONCLUSIONS
We have made a reconstruction of LIA glacier extent and LIA ice thickness, and we have interrogated the geomorphology of LIA moraine ridges, proglacial areas, and contemporary ice margins.Overall, we interpret that these multiple lines of evidence show that the Qaanaaq ice cap land-terminating outlet glaciers had more extensive temperate basal conditions, subglacial water movement, and produced landforms by both active ice and glacifluvial processes during the LIA than at present.Indeed, the glaciers at present contain abundant evidence of widespread cold-based ice.Our suggestion of 'cooling glaciers' is not new and has been suggested conceptually by several studies (cf.Carrivick et al., 2012;Glasser & Hambrey, 2001;Hodgkins et al., 1999;Karušs et al., 2022;Sevestre et al., 2015;Van Pelt et al., 2016).The concept is that thicker glaciers experience larger driving stresses and hence greater velocity thereby more likely generating pressure melting at the bed.In contrast, as some Arctic, subpolar, Our findings therefore improve our understanding of past glacier changes and should assist with projecting future glacier changes.We contend that glacier thermal regime (and transitions) should be included in glacier evolution models as a key factor controlling glacier dynamics.This inclusion will have importance in other world regions with polythermal glaciers, across the Arctic and the Antarctic Peninsula, for example, which are also the world regions where climate is changing especially fast.
impressive dimensions of this ridge most likely reflect the exceptionally steep and proximal cliffs narrowly bounding and funnelling this glacier and providing abundant angular boulders.Moraine ridges attributed to the LIA at Qaanaaq glacier have degraded to form distinct cones or other remnant parts of terraces (Figures

F
I G U R E 2 Rates of surface elevation change within glacier ablation areas between the LIA and approximately year 2015 (a) and between 2000 and 2019 from Hugonnet et al. (2021) (b).The coverage (with percentage of each outlet glacier ablation area labelled) of glacier beds that are likely cold-based during the LIA (c) and between 2000 and 2019 (d) as based on modelled ice thickness.[Color figure can be viewed at wileyonlinelibrary.com]

F
I G U R E 3 (a) Microtopography, geomorphology, land surface texture, and glacier surface structure at Sydgletscher.We highlight parts of the proglacial area (b) and the contemporary ice margin (c) using a hillshaded image of a high-resolution (0.2-m XY grid) digital elevation model, a mosaic high-resolution (0.02 m) orthophotograph, and as elevation deviations from a 10-m spatial mean.The full extent of these datasets is depicted in Supporting Information S1. [Color figure can be viewed at wileyonlinelibrary.com] abandoned glacier ice ramp, and that embayment is walled by very steep ice slopes and floored with silty-sand sediment (Figure3c).All glacier ice margins are to a greater or lesser degree being submerged in debris, as basal sediment (Figure8b,c) emerges onto and drapes the ice surface, mostly via thrust planes (Figure8d,e), producing till plains (Figure9d).

F
I G U R E 4 (a) Microtopography, geomorphology, land surface texture, and glacier surface structure at an unnamed glacier west of Qaanaaq Glacier.We highlight parts of the proglacial area (b) and the contemporary ice margin (c) using a mosaic high-resolution (0.02 m) orthophotograph, a hillshaded image of a high-resolution (0.2-m XY grid) digital elevation model, and as elevation deviations from a 10-m spatial mean.Arrows in (c) denote major abandoned water pathways.The full extent of these datasets is depicted in Supporting Information S1. [Color figure can be viewed at wileyonlinelibrary.com]F I G U R E 5 (a) Microtopography, geomorphology, land surface texture, and glacier surface structure at Fan Glacier.We highlight parts of the proglacial area (b) and the contemporary ice margin (c) using a mosaic high-resolution (0.02 m) orthophotograph, a hillshaded image of a highresolution (0.2-m XY grid) digital elevation model, and as elevation deviations from a 10-m spatial mean.Arrows in (c) denote dry gorges and major abandoned water pathways.The full extent of these datasets is depicted in Supporting Information S1. [Color figure can be viewed at wileyonlinelibrary.com]4.5 | Contemporary ice surfacesAlthough from a distance, the Qaanaaq glaciers appear to have extremely clean ice surfaces, closer inspection reveals a variety of structural, sedimentary, and hydrological features.All the glacier surfaces have relict crevasse traces or thrust planes (FigureS5) and longitudinal foliation, which is a layered structure that results from deformation(Jennings & Hambrey, 2021) (Figure9f).Together with a complete absence of open crevasses, these ice surface features indicate low shear stresses due to slow ice velocity (generally <50 mÁyear À1 fromMillan et al., 2022, dataset).

F
I G U R E 6 (a) Microtopography, geomorphology, land surface texture, and glacier surface structure at Scarlet Heart Glacier.We highlight parts of the contemporary ice margin (b) and the proglacial area (c) using a mosaic highresolution (0.02 m) orthophotograph, a hillshaded image of a high-resolution (0.2-m XY grid) digital elevation model, and as elevation deviations from a 10-m spatial mean.Arrows in (c) denote major abandoned water pathways.The full extent of these datasets is depicted in Supporting Information S1. [Color figure can be viewed at wileyonlinelibrary.com] englacial conduits are stacked vertically over $10 m, and their elliptical morphology and arrangement suggests successive down-cutting and abandonment, akin to the cut-and-closure model (Gulley support for this interpretation comes from the observation of contemporary floods from Qaanaaq glaciers caused by high air temperatures F I G U R E 7 Examples of the morphology and composition of prominent moraine ridges attributed to the LIA around Qaanaaq ice cap.Complex of subparallel lateral moraine ridges at Fan Glacier (a), relatively steep and high (50 m) relief ridges at an unnamed glacier immediately to the west of Scarlet Heart Glacier (b), degraded and collapsing ridge/terrace of subrounded boulders at an unnamed glacier immediately to the east of Qaanaaq town (c), Sydgletscher with obvious lateral moraine ridges and its ice-marginal embayment floored with silty-sand sediment (d), and subrounded boulders and abundant matrix within the single ridge that forms the most westerly and $5 m high part of that moraine ridge at Fan Glacier (e).[Color figure can be viewed at wileyonlinelibrary.com]and intense rainfall(Kondo et al., 2021).The proglacial till plains at Scarlet Heart Glacier and at Fan Glacier evidence prolonged uninterrupted 'active' ice-margin retreat since the LIA maximum as well as melt-out of basal debris.Correspondingly, these two glaciers are those with the lowest proportions of likely cold-based ice (Figure2c,d).

F
I G U R E 8 Examples of contemporary ice margin composition (August 2022).Subglacial conduit at terminus of Scarlet Heart Glacier with silt/ sand-sized sediment within glacier ice planes (in conduit roof) and 2-m-thick frozen basal sediment composed of subrounded heterolithic clasts (a), englacial sediment with gravel and cobble-sized subangular clasts emerging; projecting upwards out of image, at terminus of an unnamed glacier immediately east of Qaanaaq town (b), Sydgletscher ice margin with debris bands and basal sediment emerging; planes projecting upwards right to left in image (c), Fan Glacier ice margin with boulder-strewn surface as basal sediment emerges over a $50-m-wide arc.Note low concentration of clasts in lower 1.6 m of individual layers (d) and 8-m-high exposure of emerging (thrust upwards) basal glacier ice at Qaanaaq Glacier (e).[Color figure can be viewed at wileyonlinelibrary.com]F I G U R E 9 Examples of contemporary (August 2022) ice surface morphology and composition.Fine-grained (silty sand) sediments (with acute anaerobic biogenic odour) from a thrust fault plane at Sydgletscher (a), exceptionally tall dirt cone (14 m) composed of silty-sand sediment at Sydgletscher (b), thrust fault planes and associated sediment and medial moraine (beyond the person) at Sydgletscher (c), basal sediment emerging from thrust fault planes onto the ice surface across a 50-m-wide arc and causing burial of the terminus of Fan Glacier (d), multiple abandoned englacial conduits at Fan Glacier (e), and longitudinal foliation and supraglacial channel on Qaanaaq Glacier (f).[Color figure can be viewed at wileyonlinelibrary.com]In contrast, present-day Qaanaaq outlet glacier termini and ice margins and ice surfaces display abundant evidence of being coldbased and composed of cold ice, respectively.The low-angle ice surfaces and an absence of crevasses (e.g., Figure 7a) both indicate low shear stresses due to ice velocities that are sufficiently low (cf.dataset of Millan et al., 2022) enough to be associated with an absence of basal sliding and hence with cold-based ice.The boundary between likely cold-based and temperate ice that is 'modelled' to be about 500 m from glacier termini at Scarlet Heart Glacier, glaciers and that coverage has not changed significantly between the LIA (48% cold-based) and the present (45% cold-based) (Figure2c,d).It cannot be said whether the negligible changes in ice thinning between the LIA and 2015 and between 2000 and 2019 (Figure2a,b) are a result of the expansion of parts of outlet glaciers that are likely cold-based or whether the increased coverage of cold-based parts has slowed glacier motion and thereby also mass loss.Indeed, this feedback between ice thickness and subglacial condition/processes could be initiated in either direction.However, irrespective of which processes triggered the other, if we assume that the climate change across the Qaanaaq peninsula has been uniform since the LIA, then the interglacier variability in glacier elevation changes (Figure2) and dynamics (from the geomorphology and sediments) indicate a strong glaciological control on basal thermal regime.The surface mass balance of the ice cap has been measured from 2012 to present and shows significantly large annual variations(Tsutaki et al., 2017).Thinning of the ice cap from 2006 to 2010 at a rate of À1.1 mÁyear À1(Saito et al., 2016) accords with the mean rate suggested by the Hugonnet et al. (2021) dataset 2000 to 2019 (Figure 2b), and Bolch et al. (2013; their tab.2) reported a mean surface lowering of 0.6 mÁyear À1 for northwest Greenland landterminating glaciers between 2003 and 2008.
polythermal glaciers have thinned, they could have become cold(er) due to reduced driving stress and more likely penetration of cold air temperatures to the bed.However, our study is novel for drawing multiple lines of evidence (an LIA ice surface reconstruction, ice thickness changes, LIA geomorphology, and contemporary glaciological observations) together to consider this conceptual model.Testing this conceptual model is important because (i) glacier thinning as a result of rising air temperatures might result in a reduction or an increase in the spatial extent of cold-based glacier ice (see fig. 2 and citations in Irvine-Fynn et al., 2011); (ii) glacier thermal regime evolution affects glacier dynamics and hence mass loss rates; and (iii) thermal regime feedbacks are not accounted for in global glacier evolution models.