The Greenland Ice-Marginal Lake Inventory Series from 2016 to 2023

Ice-marginal lakes form at the edge of the Greenland Ice Sheet and its surrounding peripheral glaciers and ice caps (PGIC), where outflowing glacial meltwater is trapped by a moraine, or by the ice itself, and create a reservoir that is in contact with the adjacent ice. While glacial meltwater is typically assumed to flow directly into the ocean, ice-marginal lakes temporarily store a portion of this runoff, influencing glacier dynamics and ablation, ecosystems, and downstream hydrology. Their presence, and change in abundance and size, remain under-represented in projections of sea level change and glacier mass loss. Here, we present an eight-year (2016–2023) inventory of 2918 automatically classified ice-marginal lakes (≧0.05 km²) across Greenland, tracking changes in lake abundance, surface extent, and summer surface temperature over time. Fluctuations in lake abundance were most pronounced at the north (22 %) and northeast (14 %) PGIC margins and the southwest Ice Sheet margin (8 %). Over the study period, an increase in surface lake area was evident at 283 lakes, a decreasing trend was evident at 240 lakes, and 1373 remained stable (±0.05 km²). The northeast region contained the largest lakes, with a median size of 0.40 km² at the ice sheet margin and 0.24 km² at PGIC margins. Average summer surface temperatures fluctuated between 3.8 °C (2018) and 5.3 °C (2023), with spatial and temporal trends identified with possible links to lake setting and size. Validation against manually identified lakes showed 64 % agreement, yielding an error estimate of −809 lakes (36 %), while lake area uncertainty was ±5 %. Surface temperature estimates showed strong agreement with in situ measurements (r²=0.87, RMSE =1.68 °C, error ±1.2 °C). This dataset provides a crucial foundation for quantifying meltwater storage at ice margins and refining sea level contribution projections while supporting research on glacier-lake interactions, Arctic ecology, and environmental management. The inventory series is openly accessible on the GEUS Dataverse (https://doi.org/10.22008/FK2/MBKW9N, How et al., 2025) with full metadata and documentation, and a reproducible processing workflow.

1 Introduction

The Greenland Ice Sheet and its peripheral glaciers and ice caps (PGICs) are forecast to be the largest cryospheric contributor to sea level rise over the coming century (Arctic Monitoring and Assessment Programme, 2021). At present, these projections assume that meltwater from the Greenland Ice Sheet flows directly into the ocean, yet a portion of this is known to be stored temporarily in ice-marginal lakes, along the ice edge of the Greenland Ice Sheet and in front of, and beside, surrounding PGICs. An ice-marginal lake is a reservoir of meltwater that is dammed by a moraine or by the ice itself, therefore the trapped meltwater is in contact with the ice margin. The delayed release of meltwater at the ice margin is a significant, dynamic component of terrestrial storage, as well as a substantial CO₂ sink and part of the hydrological system (St. Pierre et al., 2019).

Many ice-marginal lakes can be persistent and stable (Carrivick and Quincey, 2014), but an increasing number are recognised to be highly dynamic systems (Dømgaard et al., 2024). For example, many ice-marginal lakes are prone to sudden and short-lived drainage events, thereby producing GLOFs (Glacial Lake Outburst Flood events) (e.g., Taylor et al., 2023) which are also referred to as jökulhlaups (Icelandic) (e.g., Eibl et al., 2023) or sermimit supinerit (direct translation into Kalaallisut, West Greenlandic; in singular sermimit supineq) (Oqaasileriffik, personal communication, November 2024). GLOFs in Greenland have caused glacier speed-up events (e.g., Kjeldsen et al., 2017), influenced downstream erosion and sedimentation rates (e.g., Carrivick et al., 2013; Tomczyk et al., 2020), and water salinity (e.g., Grinsted et al., 2017; Kjeldsen et al., 2014).

The presence of an ice-marginal lake introduces a suite of thermo-mechanical processes, including lacustrine submarine melting and calving, that can dictate glacier margin morphology, dynamics and exacerbate ice mass loss (Warren and Kirkbride, 2003; Röhl, 2006; Carrivick and Tweed, 2019; Sutherland et al., 2020; Zhang et al., 2023b). With continued retreat of the Greenland Ice Sheet under a warming climate, ice-marginal lakes are expected to become more abundant, larger and warmer; and will likely amplify lacustrine-driven proglacial melt rates, GLOF events, and terrestrial meltwater retention (Carrivick and Tweed, 2016; Grinsted et al., 2017; Shugar et al., 2020; Carrivick et al., 2022; Dye et al., 2021; Dømgaard et al., 2023; Lützow et al., 2023; Rick et al., 2023; Veh et al., 2023; Holt et al., 2024; Zhang et al., 2024). However, ice-marginal lakes and their associated processes are largely absent from sea level change projections, which assume an immediate meltwater contribution to the ocean. This assumption overlooks the role of ice-marginal lakes as intermediary storage, and changes in lacustrine and hydrological conditions, caused for instance when glaciers retreat onto land.

Mapping ice-marginal lakes is challenging due to the variability in lake characteristics. Remote sensing has been a viable approach for mapping the presence and surface extent of ice-marginal lakes, as demonstrated by inventories in Greenland (How et al., 2021), Alaska (Rick et al., 2022), Norway (Andreassen et al., 2022), Svalbard (Wieczorek et al., 2023) and High Mountain Asia (Chen et al., 2021). In general, classification approaches have been established to identify water bodies from SAR (Synthetic Aperture Radar) and multi-spectral (i.e. red, green, blue, near-infrared, shortwave) imagery, along with water potential identification using sink analysis from Digital Elevation Models (DEMs). As Greenland covers a large latitudinal range, ice-marginal lakes have very varying conditions which make them difficult to classify through one adopted method (How et al., 2021). For instance, surface sediment load, ice, and snow cover can vary significantly, with perennial ice cover in some cases at high latitudes and elevations (e.g., Mallalieu et al., 2021). Accordingly, multi-method classification approaches are necessary to capture this diversity (How et al., 2021).

Existing ice-marginal lake inventories are often static and therefore do not capture the dynamic nature of these lakes, nor capture new lakes and retire detached lakes once the margin has retreated. Given that ice-marginal lakes are projected to increase in size and abundance over time (Shugar et al., 2020; Zhang et al., 2024), it is of utmost importance to generate time-series that adequately capture ice-marginal lake change and could potentially contribute to future sea level assessments.

Here, we present an annual series of Greenland ice-marginal lakes from 2016 to 2023, classified using an established multi-method remote sensing approach. Each annual inventory maps the number of lakes (i.e. abundance) and lake surface area, along with attributes such as known lake name and water surface temperature estimations. These inventories reveal evolving lake conditions that support future assessments of sea level contribution, lake response to climate change, ecosystem productivity, and biological activity associated with the Greenland Ice Sheet and the PGICs.

2 Data description

2.1 Dataset overview

The Greenland ice-marginal lake annual inventory series is a follow-on from the 2017 Greenland ice-marginal lake inventory, largely adopting the same classification approach, data structure and formatting (How et al., 2021). How et al. (2021) provided the first Greenland-wide ice-marginal lake inventory as a static dataset, building upon regional multi-temporal efforts, such as the southwest inventory classified from Landsat imagery (Carrivick and Quincey, 2014).

The dataset presented consists of a series of annual inventories, mapping the presence and extent of ice-marginal lakes across Greenland (Fig. 1). Ice-marginal lakes are defined as water bodies ≥0.05 km² (based on satellite image spatial resolution), which are immediately adjacent to (and therefore in contact with) the Greenland Ice Sheet and/or the PGICs of Greenland. The annual inventory series spans the entirety of Greenland, including all terrestrial regions. Thus far, there are 8 annual inventories, covering the Sentinel satellite era from 2016 to 2023, where one inventory represents one year.

Figure 1 An overview of the abundance of lakes in ice-marginal lake inventory series, 2016–2023. Each mapped point denotes a unique lake, mapped across the Greenland Ice Sheet margin (blue) and the surrounding PGIC margins (white). The tables associated with each region present general statistics (with average surface temperature values provided in °C). Red starred points on the map correspond to the largest lake of each region, with placenames sourced from the placename database provided by Oqaasileriffik (the Language Secretariat of Greenland) and inventory identification numbers presented where a name is not given. It is noted that the name of the largest lake in the CE region (Catalina Lake) is not present in the placename database, and instead we adopt the lake name from Grinsted et al. (2017). The catchment regions are those defined by Mouginot and Rignot (2019). Base maps for plotting are from QGreenland v3.0 (Moon et al., 2023).

2.2 Data format and structure

The inventory series data are distributed as polygon vector features in an open GeoPackage format (.gpkg), with coordinates provided in the WGS NSIDC Sea Ice Polar Stereographic North (EPSG:3413) projected coordinate system. File names follow the convention defined in the original 2017 Greenland ice-marginal lake inventory (Wiesmann et al., 2021; How et al., 2021): <inventory-year>-<funder>-<project-acronym>-IML-f<version-number>.<file-extension$>.$

Each GeoPackage file contains metadata regarding the lake description, physical measurements, lake surface temperature, method/s of classification, verification and possible editing (Table 1). A key piece of metadata to highlight is the lake identification number (“lake_id”, Table 1), which are assigned to each classified ice-marginal lake, often consisting of multiple polygon features and/or classifications. These unique identifications are compatible across inventory years, therefore supporting time-series analysis and comparison across inventories.

Table 1 Summary of metadata included with each ice-marginal lake inventory in the annual series.

Information is included regarding whether the adjacent ice margin is either the ice sheet or PGIC (“margin”, Table 1). This margin information originates from the MEaSUREs (NASA Making Earth System Data Records for Use in Research Environments) GIMP (Greenland Ice Mapping Project) 15 m ice mask, where the provided GeoTiff files were merged, transformed into a vector format and then assigned as either an ice sheet or PGIC margin (Howat et al., 2014; Howat, 2017). In addition, each ice-marginal lake is assigned a region – north-west (NW), north (NO), north-east (NE), central-east (CE), south-east (SE), south-west (SW), and central-west (CW) (“region”, Table 1). These regions and their corresponding names are based on ice sheet catchment regions from Mouginot and Rignot (2019), which are used to also extend to the terrestrial periphery beyond the ice sheet. By doing so, regional trends can be identified from ice-marginal lakes with a PGIC margin as well as the ice sheet (Carrivick et al., 2022).

Lake names (“lake_name”, Table 1) are assigned in instances where a name is available, with preference to West Greenlandic (Kalaallisut) placenames followed by Old Greenlandic and alternative foreign placenames. Placenames are provided by Oqaasileriffik (the Language Secretariat of Greenland) placename database (https://nunataqqi.gl/, last access: 12 November 2025), filtering placenames to those associated with lake features. The placename database is distributed with QGreenland v3.0 (Moon et al., 2023).

A readme file is included with the dataset that outlines the data file contents and terms of use. An additional data file is included which is a point vector GeoPackage file representing all identified lakes across the inventory series (presented in Fig. 1). This includes manually identified lake locations that are not captured in the inventory series using the automated classification approaches.

3 Methodology

3.1 Lake classification

Ice-marginal lakes (≧0.05 km²) were identified using three established classification methods, from Synthetic Aperture Radar (SAR) and multi-spectral imagery, and Digital Elevation Models (DEMs). This approach is based on the production of the 2017 Greenland ice-marginal lake inventory, which is summarised in Fig. 2 (How et al., 2021). The main progression (and therefore difference) is that the processing pipeline is now unified and operates through Google Earth Engine to conduct the heavy image processing, and filtering/post-processing conducted with open-source spatial packages in Python, namely geopandas (Jordahl et al., 2020) and rasterio (The Rasterio development team, 2013–o). The Python pipeline is deployable as a package called GrIML, which is accompanied by thorough documentation and guidelines (https://griml.readthedocs.io, last access: 12 November 2025) (How, 2025b). This ensures a high level of reproducibility and transparency that adheres to the FAIR (Findability, Accessibility, Interoperability, and Reusability) principles (Wilkinson et al., 2016).

Figure 2 A visualisation of the processing workflow for the generation of the ice-marginal lake inventory series, including components performed with © Google Earth Engine (“Classify”) and the Python package GrIML (How, 2025b), which utilises Python spatial data handling packages geopandas (Jordahl et al., 2020) and rasterio (The Rasterio development team, 2013–o). The workflow is based on How et al. (2021). The annotated rectangles refer to process stages (reading from left to right), the coloured annotated circles represent intermediary outputs to the corresponding process stages in the workflow, and the grey annotated circles represent workflow inputs.

3.1.1 SAR backscatter classification

Water bodies were classified from Sentinel-1 GRD (Ground Range Detected) scenes, which are dual-polarization C-band SAR images (Table 2). Scenes were pre-processed using the Sentinel-1 Toolbox to generate calibrated, ortho-corrected data, specifically thermal noise removal, radiometric calibration and terrain correction using either the SRTM 30 or ASTER DEM. Scenes were then filtered to IW swath and HH polarisation, with image acquisitions limited to the summer months (1 July to 31 August) of each inventory year (2016 to 2023).

Table 2 Summary of satellite data sources.

^* Cloud cover percentage refers to the granule-specific cloudy pixel percentage for the individual scene footprint.

Averaged mosaics for each year were derived from all summer scenes for each year of the inventory series at a 10 m spatial resolution. These mosaics were smoothed using a focal median of 50 m. Classifications were derived from these averaged and smoothed mosaics using a static threshold trained for detecting open water bodies (How et al., 2021). The derived lake extents are direct classifications of the presence of water. Therefore, such lakes are suitable for use in analysis of lake area changes.

3.1.2 Multi-spectral indices classification

Water bodies were classified from Sentinel-2 MSI (Multispectral Imager), TOA (Top-Of-Atmosphere), Level-1C scenes acquired for the summer months (1 July to 31 August) of each inventory year (2016 to 2023) with a maximum cloud cover of 20 % (Table 2). Clouds were masked using the cloud mask provided with each scene (QA60), masking out opaque and cirrus clouds. The bands of interest were extracted, specifically blue (B2), green (B3), red (B4), near-infrared (B8), and the two shortwave infrared bands (B11, B12) (Table 3). The shortwave infrared bands were resampled from 60 to 10 m spatial resolution, and then averaged band mosaics were produced from all summer scenes for each inventory year.

Table 3 Summary of multi-spectral indices for ice-marginal lake classification from Sentinel-2 Level 1C scenes^*.

^* B2 = Blue (10 m); B3 = Green (10 m); B4 = Red (10 m); B8 = Near-Infrared (NIR) (10 m); B11 = Shortwave Infrared 1610 nm (SWIR1) (20 m); B12 = Shortwave Infrared 2190 nm (SWIR2) (20 m).

Five well-established spectral indices were used to classify open water bodies: (1) Normalised Difference Water Index (McFeeters, 1996); (2) Modified Normalised Difference Water Index (Xu, 2006); (3) Automated Water Extraction Index (with shadow) (Feyisa et al., 2014); (4) Automated Water Extraction Index (no shadow) (Feyisa et al., 2014); and (5) Snow brightness ratio (How et al., 2021). The equations, band information, thresholds and strengths of each indices is summarised in Table 3. Thresholds for the indices were chosen based on previous studies of ice-marginal lakes (How et al., 2021; Shugar et al., 2020), where a pixel was classified as water only if it met the threshold criteria for all five indices.

Classifications using the multi-spectral indices approach are suitable for detecting lake extent and area change given that it is a direct classification of water. In addition, the lake extents serve as valuable comparison to those derived using the SAR backscatter classification approach. For example, an error estimate of lake area can be calculated according to the difference in lake area derived from both approaches.

3.1.3 Sink classification

Water bodies were classified from the ArcticDEM 2-metre mosaic (version 3), which is processed from high-resolution optical stereo satellite imagery by an adapted version of the Ames Stereo Pipeline. The ArcticDEM mosaic is then compiled from the best quality ArcticDEM strip files and manually adjusted to form a static data product (Shean et al., 2016; Porter et al., 2018).

The mosaic was smoothed using a focal median of 110 m, and DEM depressions (i.e. where water pools) were filled over a 50-pixel moving window and subsequently subtracted from the original mosaic; producing the outline of a lake (Table 2).

It is noted that this is an indirect water classification method compared to the former two approaches (which directly detect water). Therefore, validation was required to confirm the presence of water in classified DEM sinks, which will be elaborated further in Sect. 3.1.4. Whilst lakes classified using this approach are unsuitable for lake extent analysis, they are highly useful for analysis of lake abundance, and the change in the number of lakes over time.

3.1.4 Inventory compilation

Identified water bodies were compiled for each inventory year and filtered via two automated strategies: (1) by location; and (2) by size (Fig. 2):

1. Lakes were filtered based on their location relative to the ice margin. Here, a 1 km buffer was derived around a vectorised version of the MEaSUREs GIMP 15 m ice mask and classified water bodies were retained if they were located within the buffer (Howat et al., 2014; Howat, 2017)

2. Lakes were filtered by size, only retaining lakes above a minimum size threshold of 0.05 km² based on the spatial resolution of the source satellite imagery

Finally, each inventory year dataset was manually curated for quality purposes. Specifically, the manual curation was targeted at:

• Removing mis-classifications

• Editing classifications where, for example, the shadowing mask did not adequately remove shadowing effects

• Removing detected water bodies that did not hold water in specific years (for example, from the DEM sink classification approach)

• Removing water bodies that were detached from the ice margin

• Assigning common lake identifications in instances where a lake is composed of several bodies/polygons. This is supported by automated intersection analysis, where identifications are initially defined as overlapping water bodies across all inventory years

This manual curation was carried out via visual inspection of Sentinel-2 TOA Level-1C true colour composites from each inventory year, using averaged composites composed of scenes collected between 1 July and 31 August. Metadata for each lake was then populated across the inventory series, as outlined in Sect. 2.2 and summarised in Table 1.

For the purpose of the lake abundance and area analysis presented subsequently, common water bodies (i.e. classified with more than one method) are dissolved based on lake identification number to form the maximum possible extent. Geometries with the same lake identification number are firstly merged into a single geometry using a union join, and then all geometry attributes are aggregated and combined (Jordahl et al., 2020). This is performed in an automated and consistent manner, where all successful classifications with the same identification number are dissolved in all instances.

3.2 Summer water surface temperature estimation

A summer water surface temperature estimate was provided for each classified lake across inventory years. Water surface temperature estimates were derived from the Landsat 8 and Landsat 9 OLI/TIRS surface temperature data, which is a Collection 2, Level-2 science product that is part of a large Landsat re-processing effort (Table 2). Surface temperature estimates were generated from descending, day-lit Landsat $\mathrm{8}/\mathrm{9}$ acquisitions with thermal infrared band information (100 m spatial resolution) and auxiliary data (i.e. Top Of Atmosphere reflectance and brightness temperature), along with ASTER datasets (global emissivity and normalised difference vegetation index) and MODIS and VIIRS atmospheric auxiliary data (geopotential height, specific humidity and air temperature) (Earth Resources Observation and Science (EROS) Center, 2020; Malakar et al., 2018; U.S. Geological Survey, 2023).

Due to the lack of in situ lake surface temperature measurements in Greenland, the scheme proposed by Dyba et al. (2022) was adopted, whereby surface temperature values (LST_land) were corrected to water surface temperature (LST_water) using the following calibration:

$$\begin{array}{}\text{(1)}& {\displaystyle }& {\displaystyle }{\text{LST}}_{\text{scaled}}={\text{LST}}_{\text{land}}\times \mathrm{0.00341802}+\mathrm{149.0}\text{(2)}& {\displaystyle }& {\displaystyle }{\text{LST}}_{\text{water}}=(\mathrm{0.806}\times {\text{LST}}_{\text{scaled}}+\mathrm{54.37})-\mathrm{273.15}\end{array}$$

Where LST_scaled is the applied scale factor for computing temperature in Kelvin (K) units, and LST_water is the calibrated water surface temperature in °C (Ermida et al., 2020; NASA Applied Remote Sensing Training (ARSET) program, 2022; Dyba et al., 2022). This calibration has previously shown strong correlations against in situ measurements (average RMSE =2.8 °C and R²=0.8) from 38 lakes in Poland, highlighting accurate estimates through a simple linear calibration (Dyba et al., 2022).

A summer average water surface temperature estimate was derived using this approach, calculating an average, maximum and minimum water surface temperature value for each lake extent over each inventory year, along with the standard deviation. Scenes were filtered by a maximum cloud cover of 20 %, with acquisitions limited to the month of August to reduce the probability of ice-covered lake conditions. Lake extents were cropped by a border pixel (i.e. 100 m) to reduce the impact of edge effects. All unrealistic estimates below freezing (i.e. <0 °C) were removed, therefore reducing the influence of icebergs and ice cover on the water temperature estimates.

4 Results

4.1 Lake abundance

In total, the dataset identifies 2918 automatically delineated ice-marginal lakes across all inventory years (2016–2023) (Fig. 1). Of these lakes, 2054 share a margin with the ice sheet whilst 864 are in contact with PGICs (Fig. 3). The SW region holds the most lakes compared to other regions, with 786 classifications (640 classified at the ice sheet margin and 146 at the PGIC margins). A high abundance of PGIC lakes is found in the NO region, with 278 lakes, compared to only 37 PGIC lakes in SE. This reflects the presence of more PGICs in northern Greenland, compared to the greater ice sheet cover in the southeast.

Figure 3 Change in the abundance, total (summed) area and median lake area (km²) of ice-marginal lakes around the ice sheet margin (a) and PGIC margins (b). Each of the coloured bars denote lake abundance per region for a given year of the inventory series (2016–2023), with annotated numbers corresponding to the number of lakes classified for each region. Each line plot indicates the total and median lake area per region for a given year of the inventory series. Lake area statistics are compiled from all lakes classified from SAR and multi-spectral imagery, as DEM classifications are not a direct detection of water bodies.

Small fluctuations in the abundance of lakes are evident, fluctuating between 1963 (inventory year 2021) and 1827 (inventory year 2022) lakes (Fig. 3). The largest variation in lake abundance at the ice sheet is evident at the SW margin, with a fluctuating range of 48 lakes (8 %) between 567 (2018, 2023) and 615 lakes (2021) (Fig. 3a). The CW and SE margin experienced the least variability, only varying by 16 lakes and 15 lakes, respectively. Changes in lake abundance at the ice sheet margin do not follow any spatial or temporal trends, with fluctuations unconnected to inventory year or margin region.

Lakes at the margins of periphery ice caps and glaciers vary between 723 (inventory year 2022) and 806 (inventory year 2019), with an average range of 11 lakes at the margins of periphery ice caps and glaciers (compared to an average range of 26 lakes at the ice sheet margin) (Fig. 3b). The largest fluctuations in lake abundance are seen in the NO and NE regions, fluctuating by 22 (8 %) and 29 lakes (14 %), respectively. This is linked to the higher number of lakes in these regions, which is supported by the smallest fluctuations evident in the regions where fewer lakes are present (i.e. NW, SE, CE and CW).

4.2 Lake surface extent

The largest lake in the inventory is Romer Sø (80°${\mathrm{59}}^{\prime }{\mathrm{54}}^{\prime \prime }$ N, 19°${\mathrm{09}}^{\prime }{\mathrm{21}}^{\prime \prime }$ W), located in northeast Greenland, with a total area of 126.86 km² (Fig. 1). The average lake size is 1.29 km², and the median lake size is 0.27 km² with 2395 lakes between 0.05–1.00 km² (82 %). Only 59 lakes in the inventory series have a total area above 10 km² (2 %). When lake surface extent is subcategorised by region then the NE region holds the largest lakes with a median lake area of 0.34 km² (mean: 1.63 km²). The largest lakes sharing a margin with the ice sheet reside in the NO and NE regions, with a median lake size of 0.43 (mean: 1.78 km²) and 0.40 (mean: 1.43 km²), respectively. The largest PGIC lakes are located in the NW region with a median lake size of 0.34 km² (mean: 0.52 km²).

The inventory series also holds information on the change in lake area over time, by comparing corresponding lake extents from one of the direct classification methods (i.e. from SAR and/or multi-spectral imagery). Lake areas are first dissolved across the SAR and multi-spectral classifications, and then statistics are derived by comparing lake area through the inventory years. The total (summed) lake area for each region is presented for the purpose of assessing changes in lake area across the inventory series. In addition, median lake area statistics are presented instead of averages in order to lessen the influence of outliers, especially those originating from large lakes (Fig. 3).

The total (summed) lake area is categorised by ice margin type (ice sheet and PGIC) in Fig. 3, with subsequent time steps revealing absolute area change across each region. No inherent trends in total lake area between 2016 and 2023 are evident across all regions, such as uniform change through time. The lack of a strong trend could reflect the high climatic variability across Greenland given its large latitudinal range and exposure to variable atmospheric and oceanic climatic forcing. It could also be related to variability in lake contact with the ice margin over time, and/or varying iceberg presence which is not considered in the SAR and multi-spectral classification approaches.

Regional inter-comparison provides valuable insights into terrestrial meltwater storage. At the ice sheet margin, the SW region makes up the largest component of total lake area, fluctuating between 327.15 km² (2018) and 733.89 km² (2017) (Fig. 3a). Total lake area across PGIC margins is considerably smaller, with regional lake area varying between 2.19 km² (CW region, 2022) and 73.45 km² (NO region, 2019) (Fig. 3b). Two exceptions are the NE and SW regions, where total lake area is notably higher and likely associated with the higher number of PGIC lakes in these regions.

Change in median lake area over the ice sheet margin is relatively consistent across the inventory series, with the median across all regions fluctuating between 0.14 km² (SE region, 2018) and 0.38 km² (CW region, 2020) (Fig. 3a). Outliers are evident at the CE and NE regions in 2022, where the median lake area reaches peaks of 0.59 and 0.44 km², respectively. As a result, the CE and NE regions demonstrate most variability in median lake area compared to others. PGIC lakes have comparative median lake areas, fluctuating similarly between 0.13 km² (SE, 2018) and 0.40 km² (CW, 2022), with a high outlier in the CE region (2017, 0.67 km²) (Fig. 3b). PGIC lakes in the SW region have markedly higher median areas (0.24–0.37 km²) compared to other PGIC lake regions and their ice sheet region counterparts.

Overall, the rate of change in lake area could be analysed across 1696 lakes in the inventory series, representing 31 % of all mapped lakes. Of these, 83 lakes showed growth between 2016 and 2023, with an increase in area of greater than 0.05 km² yr⁻¹, while 240 lakes experienced a decline, with a decrease in area of greater than 0.05 km² yr⁻¹. The remaining 1373 lakes exhibited no significant change, with area variations limited to ±0.05 km² yr⁻¹ (Fig. 4). The largest rate of area change was observed in the larger lakes. The average area of lakes that expanded or contracted was 6.18 and 7.77 km², respectively, with a total combined area of 2378.75 km². In contrast, the average area of lakes that remained stable was much smaller, at 0.41 km², contributing a total combined area of 560.39 km². Therefore, while the majority of lakes experienced minimal changes in area, suggesting stability in size, they account for a much smaller proportion of the total ice-marginal lake volume compared to the larger, more dynamic lakes. It is likely that the regional statistics on lake area change are more influenced by shifts in the largest lakes than by trends across all regions and lake sizes.

Figure 4 Rate of change in lake area (km² yr⁻¹) across the ice-marginal lake inventory series, 2016–2023. Example regions are highlighted from NO (a), NE (b), CW (c), and SW (d) regions of both the ice sheet and the PGICs. Positive rate of change (over 0.05 km² yr⁻¹; blue circles), negative rate of change (less than −0.05 km² yr⁻¹; orange squares), and minimal rate of change (between −0.05 and 0.05 km² yr⁻¹; white triangles) are mapped, with the size of the symbol denoting the amplitude of the rate of change (km² yr⁻¹). Each point denotes the rate of change in lake area, calculated as a linear regression slope across all lakes classified from SAR and multi-spectral imagery in the inventory series. Lakes with no available area data (i.e. not classified using the SAR and multi-spectral imagery methods) are marked with crosshairs. The catchment regions are those defined by Mouginot and Rignot (2019). Base maps for plotting are from QGreenland v3.0 (Moon et al., 2023).

The inventory series demonstrates changes to lake morphology (and the corresponding change in ice margin morphology), of which four example scenarios are presented in Fig. 5. A classic terminus basin retreat style is evident across many ice-marginal lake extents, as presented in Fig. 5a, where terminus retreat/lake expansion is marked in the central section of the glacier outlet, leaving a trailing terminus morphology at the lateral margins. Peripheral terminus retreat is highlighted in Fig. 5b, where terminus retreat/lake expansion is focused at the lateral margins. There are instances where the presence of a lake affects the boundary conditions of two glacier termini, as demonstrated in Fig. 5c, where two glaciers terminate into the same common ice-marginal lake. And finally, there are instances displayed in the inventory series where there is margin retreat/lake expansion focused around a discrete zone, such as in Fig. 5d where a marked embayment has formed at a particular point in the north region of the glacier terminus.

Figure 5 Examples of lake morphology change, and the corresponding evolution of ice termini morphology, from the ice-marginal lake inventory series. These examples highlight basin margin retreat (a), peripheral margin retreat (b), bilateral margin retreat (c), and focused margin retreat (d). It is noted that the example from (a) is a lake with persistent ice cover throughout the summer season. Holes in the lake delineations are attributed to consistent features, such as islands, and transient features, such as icebergs. The background satellite imagery presented is from a Sentinel-2 10 m 2022 mosaic (Styrelsen for Dataforsyning og Infrastruktur, 2024). The base layers for the inset map plotting are from QGreenland v3.0 (Moon et al., 2023).

4.3 Lake surface temperature

An average surface temperature estimate was derived for each inventory lake from all available Landsat $\mathrm{8}/\mathrm{9}$ scenes acquired in the month of August for each inventory year (see Sect. 3.2). This information is provided in the metadata of the ice-marginal lake inventory series. Examining the average lake surface temperature estimate across all lakes (i.e. the sum of all lake averages divided by the number of lakes), the average lake surface temperature fluctuates between 3.8 °C (2018) and 5.3 °C (2023) (Fig. 6). The median lake surface temperature follows a similar trend, fluctuating between 3.0 °C (2018) and 4.4 °C (2023). Fluctuations year on year vary, with instances of lake temperature being lower between annual time steps (e.g. from 4.5 to 3.8 °C from 2017 to 2018), higher (e.g. from 4.8 to 5.3 °C from 2022 to 2023), and remaining consistent (e.g. 4.8 °C for 2021 to 2022). Overall, there is no evident trend in average lake surface temperature across the period. However, there is a notable anomaly in 2018 which could reflect a true climatic event but needs to be investigated further.

Figure 6 Average lake surface temperature estimates from the month of August at each inventory lake for each inventory year (2016–2023) (grey). The average of all lakes (black) is the sum of all lake averages divided by the number of lakes, corresponding to the values reported on the plot. Boxplots (blue) indicate the median and the interquartile range (25 %–75 %), with the shaded band representing the variation across the annual interquartile ranges. Lake surface temperature is derived from Landsat 8 and Landsat 9 OLI/TIRS Collection 2 Level 2 surface temperature data product. Averages are calculated from all available scenes acquired from the month of August to limit the risk of mis-estimates due to ice-covered conditions, with all estimates below 0 °C removed.

Average surface temperature can be examined spatially across each lake in the inventory series (Fig. 7a). This reveals an apparent latitudinal trend, with lakes across the northern regions (NO, NW, and NE) being cooler, on average, than those in the southern regions (SW and SE). The northern regions have a higher abundance of lakes with an average surface temperature between 0 and 4 °C, whereas lake temperature in the southern regions tends to be between 4 and 10 °C. There are visible exceptions to this spatial trend, such as ice-marginal lakes present in nunatak areas which are generally cooler because of a greater presence of ice surrounding them and/or beneath them.

Figure 7 Spatial and temporal evolution of average lake surface temperature estimates at each inventory lake. Lake surface temperature estimates for each ice-marginal lake are represented by colour in (a), where a surface temperature estimate for each lake is the average value across all inventory years. The temporal variability in lake surface temperature estimate is shown in (b), categorised by lake size across six size groups (up to 0.1, 0.1–0.2, 0.2–0.5, 0.5–1.0, 1.0–1.5, 1.5–5.0, and above 5.0 km²). Individual lake temperature estimates (grey) are overlain by the average of all lakes within each size group (black) (i.e. the sum of all lake averages in each group divided by the number of lakes). All estimates below 0 °C are removed. The base map for plotting is from QGreenland v3.0 (Moon et al., 2023), with catchment regions as defined by Mouginot and Rignot (2019).

When assessing lake temperature change from 2016 to 2023, lake size appears to influence the average surface temperature and the rate of temperature change across each inventory year (Fig. 7b). The smallest lakes (≦0.1 km²) are warmest, on average, with an average temperature of 5.45 °C and the overall average varying between 4.9 and 5.9 °C. The largest lakes (5.0–150.0 km²) are colder, with an average temperature of 3.9 °C, varying between 3.5 and 4.5 °C. Lakes with a smaller surface extent (≦0.5 km²) remained relatively consistent temperatures across each inventory year, with the average fluctuating by a maximum of 0.2 °C between 2016 and 2023 °C (Fig. 7b). The largest lakes in the inventory series (5.0–150.0 km²) experienced the largest temperature change between 2016 and 2023, cooling by an average of 0.5 °C.

5 Data quality and validation

5.1 Data quality control

Classification information is provided with the ice-marginal lake inventory series, so that the performance of each classification method can be evaluated (Fig. 8). 14 020 of all detections in the inventory series (66 %) were classified using only one of the methods, composed largely from the DEM method (Fig. 8a). 5156 detections were classified using two methods (24 %), and 2264 detections were classified using all three methods (11 %). It is noted that the number of classification methods is not a measure of certainty but instead should be interpreted as a reflection of lake appearance and its adherence to the criteria of each classification method, as well as satellite data availability.

Figure 8 Lake classifications by method across the ice-marginal lake inventory series, with an overview of lake classifications over all inventory years (a) and classifications by region for lakes adjacent to the ice sheet (b) and the PGICs (c). SAR refers to the SAR backscatter classification method from Sentinel-1 imagery, VIS refers to the multi-spectral classification method from Sentinel-2 imagery, DEM refers to the DEM sink detection method from the ArcticDEM, and listed methods refer to instances where more than one method has been used to classify a lake (e.g. “SAR & DEM”, “SAR, VIS & DEM”). The legend and colour scheme in (a) correspond to (b) and (c). Values in brackets in (a) are the absolute number of lakes corresponding to the provided percentages. The values printed on top of the bars in (b) are the total number of classifications in the given region.

The SW region was typically where most lakes were classified with all three classification methods; across both the ice sheet margin (Fig. 8b) and the PGIC margins (Fig. 8c). This is likely because the classification methods have been extensively applied and developed in the SW region, making them particularly well-suited for use there compared to other regions (e.g., Carrivick and Quincey, 2014; Carrivick et al., 2017; Kjeldsen et al., 2017). The DEM method was heavily relied upon in the NO and NE regions where direct classification of open water was challenging as lakes were more likely to be consistently ice/snow covered, and satellite image availability from Sentinel-1 and Sentinel-2 was limited (How et al., 2021).

5.2 Lake abundance error estimation

All ice-marginal lakes across Greenland were manually verified for each inventory year, including those that were not classified using the automated methods. All manually identified lakes are included in the dataset as point locations alongside centroid locations for those automatically classified, with metadata flags for differentiating them.

Across all inventory years, 4543 ice-marginal lakes were manually identified in total, of which 2915 (64 %) are captured by the automated classification approach. This forms an abundance error estimation of −809 (36 %), which reflects the undercounting of lakes by the automated classification (i.e. a negative bias). However, manually classified lakes include those under the size threshold (i.e. <0.05 km²) adopted in the automated classification approach. The under-estimation of ice-marginal lakes within the inventory series therefore, in part, reflects smaller lakes that are removed from the dataset automatically due to the minimum area filtering. Therefore, the error estimate reported is a conservative lower bound, which reflects the underestimation of ice-marginal lakes due to false negatives (i.e. where lakes are missed by the automated classification method).

5.3 Lake size error estimation

Lakes classified with both the SAR backscatter and multi-spectral indices classification approaches were compared to assess variability in footprint size and provide an error estimation of lake extent. Classifications with the sink detection approach were excluded from this analysis, as the sink detection approach is an indirect measurement of lake extent rather than a direct classification of water. The use of independently derived classification methods has previously been adopted and offers a reasonable proxy for estimating error in lake delineations (e.g., Leeson et al., 2013; Williamson et al., 2018; Moussavi et al., 2020; Lesi et al., 2022; Tom et al., 2022).

Across all ice-marginal lake classifications within the inventory series (2016–2023), 3070 lakes were successfully classified with both the SAR backscatter and multi-spectral indices classification approaches. On average, there is a difference of 1.54 km² between the two classifications, with a median difference of 0.17 km². When all lake areas are summed and compared, the total lake area error (i.e. the total areal difference between the SAR and multi-spectral classifications) is 774.94 km², forming a central percentage error estimate of ±5 % (Table 4).

Table 4 A summary of lake size error estimation by size groupings.

Total absolute error is the difference between the total lake area from SAR and multi-spectral classifications. The total percentage error is calculated from the difference in total lake area and the total maximum lake area (i.e. the largest extent classified by either the SAR or multi-spectral method). The total percentage error is provided as a central error estimate.

The absolute error in classified lake areas propagates according to lake size, with lakes over 5.00 km² (339 lakes) having a total error of 405.98 km² compared to a total error of 1.89 km² across lakes between 0.05 and 0.10 km² (195 lakes) (Table 4). However, the percentage (i.e. relative) error varies marginally according to lake size, with the central percentage error estimate lying between ±6 % (0.05–0.10 and over 5.00 km²) and ±11 % (0.10–0.20 km²) based on the size groupings presented in Table 4.

5.4 Lake surface temperature error estimation

Water surface temperature estimates were validated against all known and/or open-access in situ measurements of lake temperature in Greenland (Fig. 9). The only continuous/long-term in situ surface measurements (i.e. ≦2 m) are from six lake records in southwest Greenland – Kangerluarsunnguup Tasia (64°${\mathrm{07}}^{\prime }{\mathrm{50}}^{\prime \prime }$ N, 51°${\mathrm{21}}^{\prime }{\mathrm{36}}^{\prime \prime }$ W) and Qassi-Sø (64°${\mathrm{09}}^{\prime }{\mathrm{14}}^{\prime \prime }$ N, 51°${\mathrm{18}}^{\prime }{\mathrm{27}}^{\prime \prime }$ W) (Greenland Ecosystem Monitoring, 2024), Russell Lake (67°${\mathrm{13}}^{\prime }{\mathrm{77}}^{\prime \prime }$ N, 50°${\mathrm{07}}^{\prime }{\mathrm{63}}^{\prime \prime }$ W) (courtesy of Kristian K. Kjeldsen), and three lakes as part of the Asiaq Greenland Survey hydrological monitoring programme (Qassi-Sø 2024 measurements; Qamanersuaq, 63°${\mathrm{47}}^{\prime }{\mathrm{71}}^{\prime \prime }$ N, 50°${\mathrm{00}}^{\prime }{\mathrm{50}}^{\prime \prime }$ W; and an unnamed lake referred to as Asiaq station 924, 64°${\mathrm{12}}^{\prime }{\mathrm{99}}^{\prime \prime }$ N, 51°${\mathrm{36}}^{\prime }{\mathrm{39}}^{\prime \prime }$ W).

Figure 9 Comparison of in situ surface (≦2 m) water temperature measurements with remotely sensed temperature estimates (°C) from Kangerluarsunnguup Tasia (circle), Qassi-Sø (2011–2021) (triangle), Russell Lake (square), Qassi-Sø (2024) (cross), Qamanersuaq (star) and Asiaq station 924 (ST924) (diamond). The colour of each point corresponds to the month that the observation was collected. All observations are from afternoon acquisitions (between 13.00–15.00 UTC).

Comparison of the 133 coinciding in situ measurements with those estimated using the remote sensing approach adopted here exhibit a strong correlation (r²=0.87), with an RMSE of 1.68 °C, suggesting that the remotely sensed temperature estimates are reliable (Fig. 9). This trend appears to be consistent regardless of the time of year. An interesting cluster of data points is evident, originating from measurements taken at Qamanersuaq and Asiaq station 924 which could be related to specific lake characteristics, such as lake depth/morphology or suspended sediment concentration. An error estimation of ±1.2 °C is determined, based on the average difference from data points across all lake sites.

Ideally, a correction factor specifically for calibrating values to Greenland lakes would be adopted. However, in situ validation datasets in Greenland are sparse and the derived correction factor appears to agree well with the limited datasets available. In the future, more in situ observations would strengthen the assessment, with a possibility to derive a Greenland-specific correction scheme. In addition, the influence of ice presence is limited based on the selection of strictly summer scene acquisitions and removing temperatures below 0 °C. In future work, icebergs and ice cover could be explicitly removed before a temperature estimate is derived to reduce the influence of their presence. This is currently beyond the scope of the dataset given the challenges with classifying the presence of ice on lakes across large regions in the Arctic (e.g., Carrea et al., 2025; Dye et al., 2025).

6 Potential applications and future updates

6.1 Uses for the ice-marginal lake inventory series

The inventory series presented here is the first step to quantifying the terrestrial storage of meltwater, and how it changes over time, which would be highly valuable for refining estimations of the future sea level contribution of the Greenland Ice Sheet and surrounding PGICs. The dataset extends the value of a single, time-static inventory by providing a consistent, multi-year record of lake evolution. This temporal dimension enables analyses that capture lake variability and persistence through time. Tentative findings have been outlined, yet further analysis and evaluation against other datasets is needed to investigate causal links. For example, the inventory series could be used to address the drivers of change in lake area with comparison to potential influences such as meltwater flux, sedimentation rates, bedrock type, and GLOF magnitude and frequency (e.g. Veh et al., 2025). The inventory series could also be incorporated with mapping efforts of terrestrial lakes (i.e. no contact with the ice margin) to provide a detailed overview of dynamic and stable storage of water at the terrestrial margins of Greenland (e.g Danish Climate Data Agency, 2025). Additionally, the inventory series would be a valuable dataset for examining lacustrine terminus retreat dynamics, expanding investigations from a case study basis (e.g. Mallalieu et al., 2021; Langhamer et al., 2024) to a regional and/or national scale (e.g. Dye et al., 2022).

The ice-marginal lake inventory series is applicable to climate and cryosphere research, enabling inter-annual comparison of lake change (abundance, extent and surface temperature) over time, similar to inventories for other regions such as Svalbard (Wieczorek et al., 2023). Such inventories have been used to characterise ice dam types (e.g., Rick et al., 2022), monitor GLOFs (e.g., Lützow et al., 2023), and assess lake conditions in catchments of interest (e.g., Hansen et al., 2025). Lake conditions could also provide insights into glacier dynamics in lacustrine settings around Greenland; for example, to investigate submarine melting in lacustrine settings and its impact on glacier retreat (e.g., Mallalieu et al., 2021), and to explore regional variability and evolution of lake surface temperatures and their influence on downstream ecology (e.g., Fellman et al., 2014; Carrea et al., 2025). More widely, the lake changes documented in this inventory series would be valuable to studies of the redistribution of mass on the earth surface, affecting gravity, geodesy and lithospheric elastic response (e.g., Ran et al., 2024).

Beyond scientific research, the inventory series will also be a useful resource in Greenland's assessment of infrastructure, with hydropower being the main sector that could benefit. Given Greenland's commitment to the Paris Agreement strongly suggests the expansion of current hydropower infrastructure, the ice-marginal lake inventory series could be valuable in infrastructure assessments (Naalakkersuisut, 2023). For example, the inventory series can be used to distinguish glacier-fed lakes from catchment-fed lakes, identify draining lakes, and other characteristics that are useful to discern viable catchment regions.

6.2 The future of the ice-marginal lake inventory series

It is planned to update the ice-marginal lake inventory series annually with new inventory years, using the methodology and data sources outlined here. Further automated classification methods need to be explored and incorporated into the data production pipeline to address the under detection of ice-marginal lakes, as highlighted in the validation of automated classifications against manual detections (see Sect. 5.2). Methodologies such as Forel-Ule color indexing (FUI) have been successfully applied to the detection of lakes in the High Arctic (Urbański, 2022), and more widely in a global context (Wang et al., 2021), which could be suitable for applying to ice-marginal lakes in Greenland after testing their accuracy and feasibility. Another avenue to explore is the inclusion of manually delineated lake extents where lakes have not been identified with the automated classification approaches. However, this would add further labour to the manual curation of the inventory series. An alternative would be to look at implementing new automated classification methods with machine learning, using the existing lake classifications as the foundation of a training dataset. Lake classification aided by machine learning has been successfully used for supraglacial lake detection on the Greenland Ice Sheet, so the use of machine learning in ice-marginal lake detection is worth further investigation (Lutz et al., 2023; Melling et al., 2024).

One of the key limitations of this work to be addressed in the future is the reliance on static data products, in particular the static ArcticDEM 2 m mosaic for classification, and the MEaSUREs GIMP static ice margin for filtering. The use of static data products in the inventory series presented here highlights the importance of high-labour, time-consuming manual dataset curation. For the DEM classification, an alternative would be SAR-derived DEMs from the TanDEM-X mission or the ArcticDEM strip data product, which are both time variant, but data coverage is lacking currently and scenes covering all Greenland may not be possible from year to year (e.g., Lutz et al., 2024). Another option could be to use coarser spatial resolution DEM products, such as PRODEM (500 m) (Winstrup et al., 2024), however, smaller lakes would not be identifiable. For the ice margin filtering, machine learning ice margin products show promise in being used in future editions of the inventory series, such as AutoTerm (trained with the TeamPicks dataset) (Goliber et al., 2022; Zhang et al., 2023a). The use of dynamic ice margin datasets in the future could negate the need for generating a classification spatial buffer around the margin data and instead classify ice-marginal lakes directly from their intersection with the ice margin position.

Another opportunity for future iterations of the inventory series could be to include past years, prior to the Sentinel satellite era. However, this is limited by the open availability of SAR and multi-spectral satellite imagery at a high spatial resolution (i.e. 10 m). The addition of valuable metadata could also be explored, including information such as the characteristics and dynamics of each classified ice-marginal lake. The type of damming has been included in other inventories, proving to be useful for assessing present and future lake conditions under a changing climate (e.g., Rick et al., 2022). Incorporating known GLOFs and/or drainage periods for each lake would also provide insight into abrupt changes in terrestrial water storage and be highly valuable information for infrastructure assessments, such as hydropower utilities (e.g., Dømgaard et al., 2024).

7 Code and data availability

The dataset is openly available on the GEUS Dataverse at https://doi.org/10.22008/FK2/MBKW9N (How et al., 2025), distributed under a CC BY 4.0 license (https://creativecommons.org/licenses/by/4.0/, last access: 12 November 2025). If the dataset is presented or used to support results of any kind then we ask that a reference to the dataset be included in publications, along with any relevant publications from the data production team. If the dataset is crucial to the main findings, we encourage users to reach out to the authorship team as this will likely improve the quality of the work that uses this product. The production code for making the inventory series is openly available at https://doi.org/10.5281/zenodo.15933630 (GrIML v1.0.3) and https://github.com/GEUS-Glaciology-and-Climate/GrIML (How, 2025a, b). It is distributed as a deployable and version-controlled Python package, including Jupyter Notebook tutorials on how to run the pipeline and basic handling of the dataset. If the production code is used or adapted, then we ask for a reference to be included in publications.

8 Conclusions

Here, a series of annual inventories is presented that represent ice-marginal lake abundance, surface area extents, and surface temperature estimates across Greenland for the years 2016 to 2023. Ice-marginal lakes are mapped across the margin of the Greenland Ice Sheet and its surrounding PGICs. The dataset demonstrates lake change over the 8-year period, which can be assessed at various scales, from individual lake, to regional, to Greenland-wide change. The annual ice-marginal lake inventory series is openly available on the GEUS Dataverse with a cite-able DOI at https://doi.org/10.22008/FK2/MBKW9N (How et al., 2025) and an open, reproducible workflow (How, 2025b).

The dataset reveals small fluctuations in the abundance of lakes year on year, with the largest variations occurring at the NO (22 %) and NE (14 %) PGIC margins, and the SW Ice Sheet margin (8 %). The NE region holds the largest lakes, with a median lake area of 0.40 km² at the ice sheet margin and 0.24 km² at the PGIC margins (including Romer Sø, the largest lake in the inventory series). Between 2016 and 2023, 283 lakes grew in size, 240 shrank in size and 1373 remained the same size (±0.05 km²). A summer surface temperature estimate is provided for each lake across the inventory series, demonstrating an average temperature fluctuation between 4.3 °C (2018) and 5.9 °C (2023) and evident spatial and temporal trends influenced by lake setting and size.

The ice-marginal lake inventory series was validated against manually identified lakes to assess its accuracy in lake abundance, revealing that 64 % of manual identifications are adequately captured in the inventory series. This formed an error estimate of −809 lakes (36 %). SAR backscatter and multi-spectral indices classifications were compared to assess the uncertainty in detected lake extent, providing an error estimate of ±5 %. Lake surface temperature estimates were compared against existing in situ surface (≦2 m) lake measurements from Greenland, exhibiting a strong correlation (r²=0.87; RMSE =1.68 °C) and an error estimate of ±1.2 °C.

The annual ice-marginal lake inventory series is a valuable addition to addressing current limitations in terrestrial water storage and its influence on Greenland's future sea level contribution. This dataset is the first step towards quantifying meltwater storage at the margins of the Greenland Ice Sheet, and surrounding PGICs. It also provides insight into lake change over time, and the resulting impact on glacier dynamics, such as lacustrine frontal ablation (i.e. submarine melting and calving). Beyond the cryospheric science community, the dataset will be invaluable to related disciplines in biology and ecology, where changes in lake conditions shape Arctic ecosystems and biological activity. On a national level, the inventory series could be a useful resource in environmental management and infrastructure assessment, for instance in the expansion of hydropower utilities as suggested in Greenland's new commitments to the Paris Agreement.