Top-of-permafrost ground ice indicated by remotely sensed late-season subsidence

Abstract. Ground ice is foundational to the integrity of Arctic ecosystems and infrastructure. However, we lack fine-scale ground ice maps across almost the entire Arctic, chiefly because there is no established method for mapping ice-rich permafrost from space. Here, we assess whether remotely sensed late-season subsidence can be used to identify ice-rich permafrost. The idea is that, towards the end of an exceptionally warm summer, the thaw front can penetrate materials that were previously perennially frozen, triggering increased subsidence if they are ice rich.
Focusing on northwestern Alaska, we test the idea by comparing the Sentinel-1 Interferometric Synthetic Aperture Radar (InSAR) late-season subsidence observations to permafrost cores and an independently derived ground ice classification. We find that the late-season subsidence in an exceptionally warm summer was 4–8 cm (5th–95th percentiles) in the ice-rich areas, while it was low in ice-poor areas (−1 to 2 cm; 5th–95th percentiles). The distributions of the late-season subsidence overlapped by 2 %, demonstrating high sensitivity and specificity for identifying top-of-permafrost excess ground ice. The strengths of late-season subsidence include the ease of automation and its applicability to areas that lack conspicuous manifestations of ground ice, as often occurs on hillslopes. One limitation is that it is not sensitive to excess ground ice below the thaw front and thus the total ice content. Late-season subsidence can enhance the automated mapping of permafrost ground ice, complementing existing (predominantly non-automated) approaches based on largely indirect associations with vegetation and periglacial landforms. Thanks to its suitability for mapping ice-rich permafrost, satellite-observed late-season subsidence can make a vital contribution to anticipating terrain instability in the Arctic and sustainably stewarding its ecosystems.



Introduction
Permafrost conditions are changing across the Arctic, as evidenced by widespread observations of ground warming and increasing terrain instability Jorgenson et al., 2015;Segal et al., 2016;Box et al., 2019). Susceptibility to terrain instability is primarily governed by the presence and abundance of excess ice in the permafrost, i.e. the volume of ice 20 in the ground which exceeds the total pore volume that the ground would have under natural unfrozen conditions (French and Shur, 2010;Kanevskiy et al., 2017). As excess ice melts and the meltwater drains, the ground will settle, slump or collapse (Morgenstern and Nixon, 1971;Kokelj and Jorgenson, 2013;Shiklomanov et al., 2013). Even though such thermokarst is an infrastructure hazard, we lack accurate fine-scale estimates of excess ground ice over most of the Arctic (Heginbottom, 2002; ice content profiles e(y)  Figure 1. Late-season subsidence is predicted to be closely related to the excess ice content at the top of permafrost in an exceptionally warm (and wet) summer. When the thaw front penetrates the permafrost in the late season, the melt of excess ice in the permafrost, where present, will give rise to increased subsidence. Early-season subsidence reflects excess ice at the top of the active layer, which may also be present in units with ice-poor permafrost (top row), such as young floodplains. The subsidence time series are cumulative and relative to the beginning of the thaw season. Melvin et al., 2017). This lack is a major limitation for sustainably planning in the Arctic and for accurately anticipating how 25 ecosystems and the hydrologic cycle will change (Prowse et al., 2009).
The paucity of fine-scale ground ice maps is largely due to the fact that permafrost ground ice is not directly observable from space (Heginbottom, 2002). Current approaches for mapping ground ice have significant shortcomings. Maps obtained from palaeogeographic modelling of ground ice aggradation and degradation are currently limited to coarse scales (Jorgenson et al., 2008;O'Neill et al., 2019). For localized maps, the standard approach is to upscale costly field observations and expert 30 interpretations based on imperfect indirect associations with vegetation cover and surficial geology (Pollard and French, 1980;Heginbottom, 2002;Reger and Solie, 2008;Paul et al., 2020). This works well where near-surface perennial ground ice can be reliably excluded (such as under active floodplains; Jorgenson et al. (1998); Reger and Solie (2008)), or where there are robust indicators of excess ground ice that can be recognized using -largely manual -image analysis. These include aggradational landforms such as palsas (Borge et al., 2017) and degradational features that include thermokarst lakes, thaw slumps and ice 35 wedge pits (Dredge et al., 1999;Farquharson et al., 2016;Kokelj et al., 2017;Zhang et al., 2018). Two problems with this approach are the paucity of reliable ground ice indicators in many areas (Mackay, 1990;Jorgenson et al., 2008;Reger and Solie, 2008), and that identifying ice-rich permafrost using degradational features works best when it is already too late. What is needed is a physically based observational strategy that can identify ice-rich permafrost across a wide range of conditions.
Here, we assess remotely-sensed late-season subsidence in an exceptionally warm summer as an indicator of excess ground 40 ice at the top of permafrost. The idea is that, towards the end of a hot summer, the thaw front can penetrate materials that were previously perennially frozen (Fig. 1). If these materials do not contain excess ground ice, we generally do not expect to observe elevated late-season subsidence. If they are rich in excess ground ice, the melt of this ice is predicted to induce a characteristic late-season acceleration of subsidence (Harris et al., 2011). Late-season subsidence is a sign that the thaw front has penetrated ice-rich materials at depth, which is difficult to observe in the field. As such, it may be a subtle ( 1 cm) precursor for long-45 term instability, or even the expression of rapid and pronounced ( 5 cm) thermokarst in an exceptionally warm summer. In addition to warm temperatures, wet conditions (Jorgenson et al., 2006;Shiklomanov et al., 2010;Douglas et al., 2020) and surface disturbance (Jorgenson et al., 2010a;Liu et al., 2015) can also promote deep thaw, triggering late-season subsidence in ice-rich permafrost.
Late-season subsidence is a physically based indicator of top-of-permafrost excess ground ice. To sketch the physical con-50 nection, we will make the simplifying assumption that even on submonthly time scales, thaw consolidation equals the melt of excess ground ice in a given period (Morgenstern and Nixon, 1971). We further neglect summer heave Mackay (1983). In a simple 1D scenario, the subsidence s(t 1 , t 2 ) between times t 1 and t 2 is then equal to the total excess ice that melts during that period:

55
where e(y) is the excess ice content [−] per unit depth dy. y f (t) is the depth of the thaw front relative to the surface at the beginning of the thaw season; it is assumed to be a monotonic function of t. By judicious choice of t 1 and t 2 , one can determine at which depth to probe the excess ice content. By focusing on the late season, we intend to isolate the excess ice at the base of the active layer and top of the permafrost (Harris et al., 2011;Bartsch et al., 2019). Sensitivity to the excess ice at the top of the permafrost -segregated and massive ice -is important because it is the permafrost ground ice that determines the susceptibility 60 to terrain instability (Kokelj and Jorgenson, 2013).
Remote sensing of late-season subsidence has become possible on a circumpolar scale with the advent of the Sentinel-1 InSAR satellites (Torres et al., 2012). The repeat interval of 12 days allows subseasonal variability in subsidence to be resolved (Chen et al., 2018;Bartsch et al., 2019;Rouyet et al., 2019). With other satellites, the longer repeat intervals necessitated a shallow upper limit y f (t 1 ), thus precluding subseasonal separation of deeper excess ice from that in the upper active layer (Liu et al., 2015;Chen et al., 2020). Consequently, past InSAR studies have largely focused on long-term subsidence (Liu et al., 2015;Iwahana et al., 2016). This observation strategy can be used to detect ice-rich permafrost indirectly through thermokarst but presupposes prolonged degradation (Liu et al., 2015;Iwahana et al., 2016). Conversely, Sentinel-1 may allow for the detection of ice-rich permafrost within a single season through late-season subsidence, but its suitability for identifying ice-rich permafrost is unknown.

70
To assess the suitability of remotely sensed late-season subsidence as a permafrost ground ice indicator, we use Sentinel-1 observations near Kivalina, northwestern Alaska, from 2017-2019. The summer of 2019 was the warmest since 1980 and also wet. To determine the specificity and sensitivity of the satellite observations, we compare them to ground ice cores and to an independent ground ice map, which we derived based on manual interpretation of high-resolution optical images and field observations. Based on these assessments, we appraise the suitability and discuss the limitations of late-season subsidence as a ground ice indicator. These findings will serve to enhance the automated mapping of ground ice and anticipating terrain instability on pan-Arctic scales.

Study area
Our study area is located in the northwestern Alaskan Arctic, near the town of Kivalina (Fig. 2a). The surficial geology and topography are varied (DOWL Engineers, 1994;Tryck Nyman Hayes, 2006). The spectrum includes marine deposits near the  Alaska. According to MERRA-2 reanalysis data (Gelaro et al., 2017;Global Modeling and Assimilation Office, 2020) shown in Fig. 2b), the thawing degree days in 2019 exceeded those of the years 2008-2017 by more than a third. Average daily temperatures in 2019 were particularly elevated in late June and early July (Fig. S1a), and they also remained consistently above 10 • C in late August and the first half of September. The summer of 2018 was also warm, but not as exceptional as that of 2019. Precipitation also varied, with 2019 being the wettest and 2018 the driest of the three (Fig. S1b).  (Torres et al., 2012;Copernicus Sentinel, 2020) were acquired at a frequency of 5 GHz with VV polarization in the Interferometric Wide mode (single-look resolution of ∼10 m). Acquisitions from an easterly look direction (37 degrees incidence angle; descending orbit; path 15, frame 367) were available at 12-day intervals.
There was one exception in 2017, during which a gap of 18 days occurred.

105
To estimate a subseasonal displacement d time series from the Sentinel-1 observations Copernicus Sentinel (2020), we applied Short BAseline Subset (SBAS) processing (Berardino et al., 2002). The rationale of this Interferometric SAR (InSAR) approach is to derive displacement time series from redundant interferograms. We first formed interferograms with temporal separations of up to 24 days, using spectral diversity techniques for coregistration (Scheiber and Moreira, 2000) and removing the topographic phase contribution using the TanDEM-X DEM (Rizzoli et al., 2017;TanDEM-X DLR, 2020). After multi-110 looking to a resolution of 60 m, we unwrapped the interferograms using SNAPHU (Statistical-Cost, Network-Flow Algorithm for Phase Unwrapping; Chen and Zebker (2001)). Ionospheric phase corrections were deemed to be unnecessary. We then estimated the displacement time series, with a temporal sampling of 12 days, from the interferogram stack. We used a weighted least squares approach based on the singular value decomposition (Berardino et al., 2002), with the weights determined by the Cram'er-Rao phase variance estimate (Tough et al., 1995).

115
The time series are reported as displacements d along the line-of-sight direction, with positive values corresponding to increasing distance. Owing to the hilly terrain, we chose not to convert them to vertical displacements. However, we did not discern aspect-dependent trends that are associated with downslope movements. If downslope motion were comparable to surface-perpendicular motion, the west-facing slopes would show greater line-of-sight displacements than the east-facing slopes. For simplicity's sake, we will refer to displacements with increasing distance as subsidence. The time series record the 120 cumulative subsidence since the first radar observation.  Figure 3. Late-season and early-mid-season subsidence were reconstructed by spline fitting and their accuracy assessed at stable validation points. a) Three cubic spline basis functions, corresponding to decelerating, S-shaped and accelerating subseasonal subsidence; for display purposes, we normalized them so that the peak subsidence rate is 1 mm/day. b) Reconstructed displacements over the validation points during the early-to-mid and late season (de and d l , respectively; positive values corresponding to subsidence). Individual values for all points and years shown at the top, and a kernel density estimate of the distribution below.

Referencing and assessing its quality
To reduce long-wavelength atmospheric errors in the subsidence observations, we spatially referenced the raw time series at multiple locations with outcropping bedrock or a thin rubble veneer. These locations were assumed to be stable (Reger and Solie, 2008;Antonova et al., 2018).

125
To assess the quality of the referencing, we chose eight bedrock validation points distributed across the study region. Assuming these points to remain stable, the subsidence observations at these locations are a measure of the observational uncertainty.
This uncertainty estimate subsumes decorrelation errors and uncompensated atmospheric contributions under a single number.

Late-season subsidence from spline fitting
We estimated the late-season subsidence d l by fitting spline basis functions to the referenced subsidence time series. The 130 advantage of fitting a flexible and yet simple spline function is that measurement noise, such as residual atmospheric errors, can be reduced (Berardino et al., 2002).
To capture a range of subseasonal subsidence patterns, we used three cardinal quadratic B-splines for the subsidence rate (first derivative), corresponding to the cubic spline basis functions for the subsidence time series shown in Fig. 3a). For each pixel, we estimated the three coefficients by ordinary least squares.

135
The late-season subsidence d l was defined to be the cumulative subsidence between 10 August (t 1 ) and 10 September (t 2 ) for all years. We chose t 1 such that the thaw front in an exceptionally warm summer could plausibly have reached the average long-term permafrost table. To test the sensitivity of the d l estimates to this choice, we also computed d l starting 10 days earlier or later. The end point t 2 was chosen to minimize the confounding impact of diurnal frost heave (Chen et al., 2020). We expect the thaw front to be within <5% of its seasonal maximum at t 2 . For the specified values of t 1 and t 2 , we derived d l from the 140 spline fit.
We also estimated total subsidence during the early and mid thaw season (10 June to 10 August), d e . While d e is not expected to contain direct information about permafrost ground ice, being sensitive to soil texture and interannual variations in active layer moisture (Lewkowicz, 1992;Harris et al., 2011), it provides a simple reference for assessing the late-season speed-up characteristic of melting top-of-permafrost excess ground ice.

145
3.2 Assessment with independent ground ice data

Manual mapping
We assessed the suitability of late-season subsidence as an indicator of permafrost ground ice content by comparing the InSAR observations to an independently derived ground ice map, published as Zwieback (2020a). The map comprised two primary categories: ice rich and ice poor. To assess the sensitivity and specificity, we computed the observed late-season subsidence 150 distributions conditional on the location being ice rich and ice poor, respectively.
The independently derived ground ice map was obtained using manual interpretation and expert knowledge, drawing on field observations and high-resolution (∼1 m) satellite imagery. The mapped focus area, 8 km by 8 km in size, was chosen because of the wide range of ecotypes and the availability of field observations. A drawback of the map are its gaps, because we classified the ground ice content of 19% of the area as indeterminate. Furthermore, six percent of the area were discarded 155 in the comparison to avoid unrepresentative values over lakes and infrastructure.
The key considerations in the manual mapping were: 1. Ice-rich permafrost was assigned to ecotypes where high-resolution imagery and field observations revealed direct indicators of excess ice at the top of permafrost. Ice wedge polygons -some only barely visible or inconspicuous except along lakeshores and beaded streams, others visibly in an advanced state of degradation -are widespread (Fig. 4). They 160 are abundant in old alluvial and lacustrine deposits (DOWL Engineers, 1994;Shannon & Wilson, Inc., 2006), as well as in colluvial sediments that are rich in retransported silt (DOWL Engineers, 1994).

2.
A major caveat is that the presence of ice wedges provides no direct information on the ice content in the polygon interiors. However, previously taken cores (by Shannon & Wilson, Inc. (2006); to be discussed later) from centres in various terrain units were ice rich, and the presence and ongoing expansion of thermokarst ponds and lakes also provides 165 support for this assumption.
3. Pingos in drained lake basins are also direct indicators of excess ground ice (Mackay, 1973). 4. Exposed bedrock and rubble-covered surfaces (Fig. 4b) were classified as ice poor. This classification is generally supported by geotechnical investigations and field evidence from the study area (Pewe et al., 1958;DOWL Engineers, 1994).
Deviations from this general pattern cannot be excluded (Robinson and Pollard, 1998;French and Shur, 2010). 5. Active and recent inactive floodplains were categorized as ice poor ( Fig. 4a; DOWL Engineers (1994)) because nearsurface permafrost, where present, is too young for abundant ground ice to have aggraded (Jorgenson et al., 1998). We did not observe any indications of abundant excess ground ice in these ecotypes.
6. Such indirect indicators as vegetation cover and surficial geology were used to classify areas as ice rich where similar adjacent areas were clearly ice rich. For instance, if polygons were visible over 80% of an alluvial deposit, the entire 175 deposit was classified as ice rich (Fig. 4a).
7. The ice content of the remaining locations was deemed to be indeterminate. The majority of these locations are in uplands and on hillslopes, where ice content is inherently variable (Morse et al., 2009). In Fig. 4b, the areas above the clearly ice-rich toe and midslope deposits (Fig. 4b)  abandoned floodplains) were also classified as indeterminate (Jorgenson et al., 1998).
The mapping process is inherently subjective. The discretization of top-of-permafrost excess ground ice into a limited number of categories presents challenges (Tryck Nyman Hayes, 2006;Paul et al., 2020). The most difficult decisions were about where to draw the line between the ice-rich/ice-poor and the indeterminate category.

Cores
We also compared the late-season subsidence observations to geotechnical assessments of permafrost ground ice content from >3 m deep cores. The cores were taken in 2005 by Shannon & Wilson, Inc. to investigate the suitability of three proposed relocation sites for Kivalina (Shannon & Wilson, Inc., 2006). Where ice wedge polygons were discernable at the surface, cores were taken from the polygon centres (in one case: centre and wedge).
Ice-rich permafrost was encountered underneath all sites, but the Tatchim Isua site (Fig. 5) also encompassed an ice-poor estimates of the visible ice content. All cores but one were ice rich with visible ice contents >30% in the forms of segregated and wedge ice. The only ice-poor core, from the bench, was a gravelly soil grading into weathered bedrock (5% visible ice content in joints). Two cores taken from within the same ice wedge polygon were combined for the comparison exercise. Table   S1 summarizes the locations and ground ice properties of all 13 coring locations.

205
4.1 Spatial and temporal variability of late-season subsidence 4.1.1 Regional variability The observed late-season subsidence in our study area showed two distinct modes in 2019 (Fig. 5a). One corresponded to no or very small subsidence (-1 -1 cm). The other corresponded to regions with elevated subsidence (4-8 cm) in 2019. In 2017, the subsidence in these regions was lower (1-4 cm).

210
There were no notable instances of pronounced negative late-season subsidence estimates in any of the years. Negative estimates correspond to heave during the late season (10 August to 10 September). Small negative values (∼-1 cm) occurred in all years. These need to be interpreted in relation to the observational accuracy. Observations at the validation points (bedrock locations shown as circles in Fig. 5) correspond to a root-mean-square accuracy of 0.6 cm. The distribution shown in Fig. 3b encompasses estimates of either sign with a magnitude that is small compared to the observations in the subsiding areas.  Late-season subsidence varied with topographic position and surficial geology. Elevated late-season subsidence exceeding 3 cm was found in low-lying areas including old alluvial deposits and drained lake basins (Fig. 5a-b). Many hillslopes also featured elevated late-season subsidence in 2019. Conversely, late-season subsidence was consistently low in elevated uplands (e.g. bedrock or talus without vegetation cover) and along rivers (recent floodplains without or with dense vegetation cover.

Temporal variability 220
The late-season subsidence exhibited more year-to-year variability than the early-mid-season subsidence (Fig. 5, S2). Three examples in Fig. 6 illustrate the inter-annual variability of the subsidence time series.
A late-season acceleration of subsidence in 2019 was observed in two locations that were independently determined to be ice rich (Fig. 6 a and b). According to the radar observations, the rate of subsidence more than doubled in the late season. The rate and the total subsidence of ∼ 5cm during the late season is approximately a factor of three smaller in the years 2017 and 225 2018.
For the ice-poor floodplain in Fig. 6c, the interannual variability in late-season subsidence was less than 1 cm. The magnitude of the late-season subsidence was always small, comparable to the observational uncertainty in all years.

Suitability in an exceptionally warm summer 230
Late-season subsidence in the hot and wet summer of 2019 was markedly different for ice-rich and ice-poor areas (Fig. 7a). It exceeded 4 cm at all cored locations rich in top-of-permafrost ground ice. For the ice-rich areas, as determined by independent manual mapping , d l varied between 4 and 8 cm (5th and 95th percentile; Fig. 7c-d). For ice-poor areas the 5th-95th percentile range was -1 to 2 cm. Only in the extreme tails (2%) do the distributions overlap, ensuring a robust separability of the two classes. Based on the distributions, we applied a threshold of 2.5 cm to separate ice-rich from ice-poor permafrost (Fig. 7e).  Figure 8. Contours plot of a kernel density estimate of the early-mid-season and late-season subsidence for ice rich (purple) and ice poor (grey), as determined independently by manual mapping (Fig. 7d). The markers correspond to the values observed at the location of the ice cores (triangles in Fig. 5). Negative values, corresponding to heave, are of comparable magnitude to the referencing accuracy of ∼ 1 cm.
[cm] The separability based on the 2019 late-season subsidence d l is better than that based on the early-mid-season subsidence d e . The d e distributions of ice-rich and ice-poor areas overlap (Fig. 8c), whereas the d l distributions are concentrated around two separate peaks (see also Fig. 7a). When the start of the late-season period is pushed backward, the separability based on d e improves whereas that based on d l decreases (Fig. S3).
The candidate relocation site Tatchim Isua was characterized by a narrow zone with low late-season subsidence (Fig. 9).

240
This ∼100 m wide zone largely coincides with the gravel-covered bench that a single core from 2005 (Tab. S1) indicates to be ice poor (Shannon & Wilson, Inc., 2006). Late-season subsidence was elevated (∼7 cm) at the location of the seven cores  Fig. 7c and 5d). First row a-d): points that were independently determined to be ice rich; second row e-h) indeterminate according to manual mapping, but d l > 2.5cm indicates they are ice rich; third row i-l): independently determined to be ice poor; fourth row m-p): indeterminate according to manual mapping, but d l indicates they are ice poor.
further downslope or upslope. All cores contained ice-rich materials at the top of permafrost, but in the field the ice-rich nature was not readily apparent at the surface near the proximal coring locations (Shannon & Wilson, Inc., 2006).
Further examples from a range of geologic settings serve to illustrate the suitability of d l for identifying ice-rich permafrost.

245
Fig. 10a-d shows instances of ice-rich permafrost with ice wedge polygons. They all exhibited elevated d l of 4-8 cm, corresponding to an increased subsidence rate during the late season. Conversely, the observed subsidence was below 1 cm for the points shown in Fig. 10i-l, which were independently determined to be poor in ground ice.
The most interesting cases are those 19% of the area where the manual ground ice mapping was indeterminate because there was no strong evidence for either category. They exhibited a bimodal distribution of d l (Fig. 7a). The larger mode d l ∼ 5 cm, 250 comparable to those of ice-rich areas. Examples include colluvial hillslope deposits without conspicuous ground ice indicators ( Fig. 10e-h). The characteristic late-season acceleration in subsidence leading to elevated d l ∼ 5 cm uniquely indicates that the top of the permafrost is ice rich. The smaller mode d l ∼ 1 cm roughly matches the observations over ice-poor terrain. The negligible subsidence observed at the locations shown in Fig. 10m-p suggests that the materials that thawed late in summer contained little excess ice.

Suitability in cooler years
In the average year of 2017, the suitability was reduced because the late-season subsidence distributions of ice-rich and icepoor regions overlapped substantially (Fig. 8a). The late-season subsidence d l of 80% of the terrain that was mapped to be ice rich and all the ice-rich coring locations was less than 2 cm. On average, it was 4 cm, or a factor of 5 smaller, than during 2019 (standard deviation: 1 cm).

260
In the warm summer of 2018, the separability based on the d l distributions was intermediate. The distributions overlapped at the 10% level (Fig. 8b), compared to 2% in 2019. This is largely due to the smaller late-season subsidence of ice-rich terrain compared to the exceptionally warm summer of 2019.

265
Comparing the late-season subsidence to the independently determined ground ice map and the ice cores, we note that ice-poor and ice-rich permafrost are well distinguishable in the exceptionally warm summer of 2019 (Fig. 7). Ice-poor areas were stable, whereas ice-rich areas exhibited pronounced (∼ 5 cm) late-season subsidence (Fig. 10).
The late-season speed-up in subsidence is interpreted as the thaw front penetrating ice-rich materials at the top of permafrost (Shur et al., 2005;French and Shur, 2010). A characteristic feature is that the subsidence rate increased up to fivefold in late 270 summer (Fig. 10). Harris et al. (2011) interpreted the late-season acceleration observed in-situ in Svalbard by the sharp contrast in excess ice between the ice-rich layer and the materials above. In practical terms, the acceleration increases the robustness of the separability with respect to the chosen starting day of the late season period (Fig. S3), thus facilitating the identification of ice-rich materials at the top of the permafrost. If the ice-rich layer is thick, the late-season acceleration in subsidence can be a precursor for long-term terrain instability (Kanevskiy et al., 2017).

275
Inter-annual variability in late-season subsidence of ice-rich areas poses challenges for ground ice mapping. Potential drivers of inter-annual variability and trends include surface changes (e.g., soil moisture, disturbance, snow) and variable meteorological conditions such as precipitation (Shiklomanov et al., 2010;Bartsch et al., 2019;Douglas et al., 2020). In contrast to the exceptionally warm and wet summer of 2019, the limited late-season subsidence of ∼ 2 cm in the warm summer of 2018 in the ice-rich area shown in Fig. 6a indicates the thaw front did not sufficiently penetrate into the ice-rich materials. That the 280 identification strategy presupposes an initial degradation of excess ground ice constitutes its biggest limitation.

Limitations
The excellent separability in our study area in the exceptionally warm summer does not imply that interferometric observations of late-season subsidence are a universally applicable basis for mapping ice-rich permafrost. Limitations arise from observational uncertainties and from the imperfect sensitivity and specificity of subseasonal subsidence as an indicator of 285 top-of-permafrost excess ground ice.

Observational uncertainties
Observational uncertainties chiefly arose from errors in the referencing, dominated by uncompensated atmospheric contributions, and from location-specific systematic and random errors.
The errors due to imperfect referencing were determined to be substantially smaller (1 cm; 3b) than the typical difference in 290 late-season subsidence between ice-rich and ice-poor permafrost terrain (5 cm; Fig. 7a).
This uncertainty metric does not account for systematic biases associated with changes in soil or vegetation moisture, as it was obtained from rocky, non-vegetated surfaces. Soil moisture commonly increases toward the end of the thaw season, which would correspond to a spurious subsidence signal (De Zan et al., 2015). However, the worst-case estimates of the bias (∼ 1 cm at C-band) are a factor of five smaller than the late-season subsidence observations (Zwieback et al., 2017). Vegetation moisture 295 in shrubs will decrease with senescence, corresponding to a spurious heave signal of ∼1 cm (Zwieback and Hajnsek, 2016).
The persistently small magnitude of the displacement estimates over shrub-covered inactive floodplains (Fig. 5) indicates that these systematic errors were not a major confounding factor in the present study. However, dedicated in-situ observations are needed to accurately characterize the observational uncertainties. We expect these errors to be greater in densely vegetated, often discontinuous, permafrost (Wang et al., 2020). The ample suitability of late-season subsidence for mapping ice-rich permafrost that we identified in the study area was arguably promoted by the exceptionally warm and wet summer of 2019. In less propitious circumstances, its sensitivity and specificity may be reduced.
Its sensitivity is impaired when excess ground ice does not manifest as elevated late-season subsidence. Such false negatives 305 occurred for instance in the warm summer 2018 but not in 2019 (Fig. 8b), indicating that in 2018 the thaw front did not penetrate substantially into the ice rich materials in the upper permafrost. Limited late-season subsidence in warm summers is promoted by the low to moderate ice contents above the ice-rich materials, i.e., the lower active layer and the transient layer (Shur et al., 2005). Consequently, identifying ice-rich permafrost using late-season subsidence calls for exceptionally warm summers such as 2019.

310
False negatives -even in an extremely warm summer -are expected to occur most commonly in the discontinuous and sporadic permafrost zone. There, disturbances such as forest fires are more likely to have obliterated the ecosystem-protected or ecosystem-driven perennial ground ice near the surface (Jorgenson et al., 2010b;Kanevskiy et al., 2012Kanevskiy et al., , 2014Paul et al., 2020), but see Burn (1997). Ice-rich permafrost may occur at depth, perhaps under a thick talik. The thermal regime of such permafrost is more complex (Jorgenson et al., 2010b;Connon et al., 2018), and melt-induced subsidence does not necessarily 315 occur primarily at the end of the thaw period. The sensitivity may further be diminished by the retardation of thaw consolidation due to inefficient drainage of the excess melt water (Morgenstern and Nixon, 1971); or by sinkhole formation and piping underneath cohesive materials such as peat (Osterkamp et al., 2000). In such challenging environments, subsidence-based ground ice mapping will benefit most from synergies with complementary mapping approaches.
The specificity is impaired when late-season displacements are not due to melt of top-of-permafrost ground ice. Elevated 320 late-season subsidence may reflect ice content at the bottom of the active layer. Cold permafrost and ample moisture supply promote ice segregation at the base of the active layer, with a seamless transition in ground ice content to the ice-rich upper permafrost (Mackay, 1981). Problems could further be caused by late-season deformation related to processes such as organic soil compaction and mass loss, or downslope transport (Stephens et al., 1984;Roulet, 1991;Matsuoka, 2001;Gruber, 2020;Zhang et al., 2020). It will be beneficial to incorporate expert knowledge and to exploit the characteristic late-season acceleration 325 shown in Fig. 6a-b in order to robustly identify ice-rich permafrost.
A final limitation of this study and the preceding discussion is that the complexity of ground ice was simplified to just two categories: ice rich and ice poor. In reality, however, excess ground ice content is a continuous parameter (Morse et al., 2009;Kanevskiy et al., 2012;Paul et al., 2020) whose magnitude could be constrained using late-season subsidence observations.
A big challenge, the sub-resolution spatial variability of ground ice, is exemplified by ice-wedge polygons. To what extent do 330 the satellite observations reflect the wedges and the polygon interior, and how does it vary with factors such as ice content of the centres and the thickness distribution of the protective layer above the wedges? Quantitative answers will require densely sampled ground ice cores (Morse et al., 2009;Jorgenson et al., 2010b).

Enhancing automated ground ice mapping
Late-season subsidence may enhance the automated mapping of permafrost ground ice. Remotely sensed late-season subsi-335 dence can be mapped on pan-Arctic scales, thanks to the global availability of Sentinel-1 data (Torres et al., 2012). A further practical advantage is that it lends itself to automation, as no manual interpretation and no calibration using in-situ cores are required. To automate the specification of the late-season period, globally available reanalysis data could be considered. Despite the potential for automation, we believe that the greatest potential lies in its synergistic use with expert knowledge, field observations and complementary mapping approaches.

340
Incorporating external constraints and expertise will be essential to counteract weaknesses of remotely sensed late-season subsidence. Site knowledge is indispensable for interpreting the stratigraphic complexity and sub-resolution variability of ground ice conditions. Most importantly, independent observations will be needed to estimate excess ice at depth and thus constrain total ice contents. External data can further assist in densely vegetated areas with larger errors (Wang et al., 2020) and mitigate false negatives when the thaw front fails to penetrate deep into the ice-rich materials (Shur et al., 2005).

345
Late-season subsidence is complementary to state-of-the-art mapping approaches based on visible manifestations of ground ice and indirect associations with vegetation cover or topographic variables. In addition to validating existing maps such as the coarse-scale product by O'Neill et al. (2019), late-season subsidence can be used synergystically in future ground ice mapping efforts. The greatest contribution of late-season subsidence will likely be where field observations are sparse or where ground ice is indistinct and poorly correlated with surface characteristics. Examples include uplands and hillslopes (Fig. 9), as well as 350 areas underlain by ice wedges from the early Holocene or the Pleistocene (Pewe et al., 1958;Dredge et al., 1999;Reger and Solie, 2008;Morse et al., 2009;Jorgenson et al., 2015;Farquharson et al., 2016).

Conclusions
We set out to assess the suitability of remotely sensed late-season subsidence for identifying ice-rich materials at the top of permafrost. We predicted that ice-rich permafrost would become detectable by enhanced subsidence toward the end of an 355 exceptionally warm summer. Our assessment was based on the comparison of Sentinel-1 satellite observations from 2017-2019 with independent ground ice data in northwestern Alaska. Our principal findings and conclusions are: 1. In the exceptionally warm and wet summer of 2019, the InSAR-derived late-season subsidence observations were large (∼ 4-8 cm) in areas that were independently determined to be rich in top-of-permafrost ground ice. Conversely, the observed late-season subsidence was consistently small (-1-2 cm) in ice-poor areas.

360
2. Distinguishing ice-rich from ice-poor terrain worked best in the exceptionally warm and wet summer, as the respective late-season subsidence distributions overlapped by less than 2%. In the preceding summers the overlap was 10% or larger.
3. The suitability of late-season subsidence as an indicator of top-of-permafrost excess ice will not always be as high as in this study. The greatest drawback is the lack of sensitivity when the thaw front does not penetrate deep into the 365 ice-rich layers, but in suitable terrains this limitation can be mitigated through the use of imagery acquired during an exceptionally warm and wet summer.
4. Remotely sensed late-season subsidence complements established techniques for estimating ground ice contents, such as manual mapping approaches that exploit visible manifestations of ground ice and indirect associations with surface characteristics. Because late-season subsidence is insensitive to excess ice at depth, it will be essential to incorporate 370 geological reasoning and indirect associations established using field observations for estimating total ice contents.
5. In summary, late-season subsidence can enhance the mapping of permafrost excess ground ice and the susceptibility to terrain instability. Its practical advantages include the pan-Arctic availability of data, the ease of automation, the independence from costly in-situ observations, and its applicability in absence of visible manifestations of ground ice.
Late-season subsidence in a warm and wet summer may be relatively small, but it can serve as an observable precursor 375 for much larger terrain changes.
Pan-Arctic expansion of ground ice mapping using late-season subsidence is timely and societally relevant. It is timely because of the widespread warming and accelerating degradation of permafrost. It is relevant because we lack accurate, fine-scale ground ice maps over essentially the entire Arctic. Late-season subsidence observations can make a vital contribution to anticipating terrain instability in the Arctic and to sustainably stewarding its vulnerable ecosystems.