Seasonal Patterns of Greenland Ice Velocity From Sentinel-1 SAR Data Linked to Runoff

Accurate projections of the mass loss from the Greenland Ice Sheet (GrIS) require a complete understanding of the ice-dynamic response to climate forcings on seasonal and interannual timescales and would greatly benefit from more observational evidence. Here, we analyze a 5-year, high-resolution data set of ice velocities of the GrIS using K -means, an unsupervised clustering algorithm, to identify ice-sheet wide characteristic seasonal flow patterns


Introduction
The Greenland Ice Sheet (GrIS) is losing mass and has been a major contributor to sea-level rise since year 2000 (Fox-Kemper et al., 2021), and the contributions from surface melting and discharge at the marine outlet glaciers contributed equally to the mass loss in the period 1992-2018 (Shepherd et al., 2020).The dynamic response of the marine outlet glaciers to atmospheric and oceanic forcings remains a key uncertainty in determining the future ice-sheet mass-balance and projecting future sea level changes (Bevis et al., 2019;Fox-Kemper et al., 2021;Goelzer et al., 2020;Slater et al., 2020).Fast flowing ice feeds into the marine outlets and carries ice accumulated in the interior to the ocean.The discharge from these marine outlet glaciers have generally increased since the 2000s (Mankoff et al., 2020), and the glaciers further respond to climate on seasonal timescales (Bartholomew et al., 2010;Davison et al., 2020;Joughin et al., 2008;Moon et al., 2014).The fast flow in the marine terminating glaciers is maintained by the high basal water pressure.Surface meltwater penetration to the bed and reduction of basal friction is one mechanism that facilitates fast seasonal flow in these outlets (Larsen et al., 2016;Rathmann et al., 2017;Zwally et al., 2002), and changes in resistive stress at the calving front can also induce faster seasonal flow in the summer (Moon et al., 2014;Rathmann et al., 2017).However, the link between surface melting, seasonal variability in ice flow and inter-annual trends is poorly understood.This is critical in a future warming climate as the area with surface melt will increase, risking an increased effect on ice velocities and the dynamic mass loss.
Abstract Accurate projections of the mass loss from the Greenland Ice Sheet (GrIS) require a complete understanding of the ice-dynamic response to climate forcings on seasonal and interannual timescales and would greatly benefit from more observational evidence.Here, we analyze a 5-year, high-resolution data set of ice velocities of the GrIS using K-means, an unsupervised clustering algorithm, to identify ice-sheet wide characteristic seasonal flow patterns.We include all areas flowing >0.3 m/d and obtain an ice-sheet wide overview of the seasonality and the interannual variability.It shows both a spatial and interannual variability in seasonal flow patterns, both along individual glaciers and between glaciers.We compare with runoff from a regional climate model and infer that the ice-sheet wide patterns are linked to the availability of water penetrating to the base of the ice.

Plain Language Summary
The Greenland Ice Sheet (GrIS) is currently losing mass but the response from the marine outlet glaciers to atmospheric and oceanic forcings is uncertain and limiting the accuracy of the estimated future mass loss.Here, we analyze a 5-year satellite based ice velocity data set with an unprecedented spatial and temporal resolution in order to investigate the variability on all fast flowing areas of the GrIS.We use a machine learning algorithm to sort the annual time series into characteristic seasonal patterns, and we compare the results to modeled runoff from a climate model.We find that individual glaciers are not classified to be one type, but the response depends on the availability of drained meltwater, and that the seasonal pattern of ice velocity varies spatially and temporally, both along individual glaciers and between neighboring glaciers.We conclude that the seasonal pattern of response to runoff provide insights to the evolution of the subglacial hydrological system during the runoff season.Understanding the response of ice flow to meltwater and how it is linked to the subglacial hydrological system is crucial for understanding the dynamic response of the ice sheet to future climate warming.
Previous studies have mostly focused on the response to surface runoff of either land-terminating (e.g., Andrews et al., 2014;Bartholomew et al., 2010;Stevens et al., 2016;van de Wal et al., 2008van de Wal et al., , 2015;;Zwally et al., 2002) or marine-terminating glaciers (e.g., Davison et al., 2020;Sakakibara & Sugiyama, 2020).Seasonal responses of the ice velocity are observed in both systems to be driven by surface melt reaching the bed of the ice, where it modulates the water pressure in the subglacial hydrological system leading to changes in basal friction and ultimately the flow of ice.The subglacial hydrological system evolves through the melt season depending on the availability of runoff (Hewitt, 2013;Schoof, 2010;Zwally et al., 2002).The coupling between the drainage of surface runoff, subglacial hydrology and ice flow dynamics is complex and related to the variability of the water pressure at the base and the efficiency of the drainage system (Bartholomew et al., 2010;Davison et al., 2019;van de Wal et al., 2015).Basal water drains through an inefficient drainage system of subglacial cavities under high pressure or through subglacial channels under low pressure that develop when the melt water exceeds a critical level (Schoof, 2010).Basal friction is linked to the basal water pressure, and a switch from inefficient to efficient basal drainage is generally associated with ice flow deceleration and depends on the amount and variability of drained water (van de Wal et al., 2015).The drainage systems may evolve throughout the year, due to increased melting and drainage during summer or due to ice flow that alters the meltwater channels.
At marine-terminating glaciers, seasonal flow patterns are not only controlled by runoff through the subglacial hydrological system, but also by terminus position and the rigidity of the proglacial ice mélange through the effects on back pressure and changes in glacier geometry (e.g., G. Cheng et al., 2022;Joughin et al., 2008;Moon et al., 2014Moon et al., , 2015)).In recent years, studies of seasonal patterns have been carried out covering larger areas using satellite observations (e.g., Moon et al. (2014) and Vijay et al. (2019Vijay et al. ( , 2021))).Moon et al. (2014) studied seasonal flow patterns near the terminus of 55 marine terminating glaciers in Greenland based on 3-6 ice velocity estimates by year, and identified 3 distinct types of seasonal flow patterns by manual inspection of the time series.Type 1 is closely linked to the position of the calving front with speed-up related to terminus retreat, while the other two are connected to runoff and related changes in basal friction.Type 2 is characterized by a speed-up in early summer and slow down to pre-melt season velocity starting in mid-summer while type 3 has mid-summer slow down to a clear late summer minimum with velocities rebounding over winter and spring.Vijay et al. (2019Vijay et al. ( , 2021) carried the analysis forward in time and included more glaciers using the newly available satellite data, but still focusing on the glacier fronts.
In this study, we extend the above analyses to all fast flowing areas of the GrIS (Figure 1) using 5 years of ice velocity mosaics at 12 days temporal resolution based on data from ESA's Sentinel-1 satellites.Due to the enormous and unprecedented number of time series, we employ an unsupervised clustering algorithm to determine ice-sheet wide characteristic seasonal flow patterns.We compare the resulting patterns to runoff from a regional climate model (RCM).This enables us to identify distinct modes of seasonal responses and to explore the large-scale response of ice flow to variations in melt water, one of the main drivers of seasonal variations in ice flow.

Data and Methods
Surface ice velocity is obtained from the PROMICE ice velocity product (Solgaard & Kusk, 2021;Solgaard et al., 2021), which is a time series of ice velocity mosaics covering the GrIS.It is derived using offset tracking on synthetic aperture radar (SAR) data from the European Space Agency's Sentinel-1 satellites.The time series is posted at a spatial resolution of 500 m and each mosaic spans two Sentinel-1 cycles (24 days).A new mosaic is available every 12 days and thus a given mosaic overlaps by 12 days with the previous and subsequent maps.The time series starts September 2016 and includes ∼30 mosaics per year.Here, we use data spanning the period 1 January 2017 to 31 December 2021.
The 5 years long time series in each grid point is split into annual series.In a given grid point, there are thus 5 time series, that is, one for each year in the range [2017,2021].In order to reduce the impact of noise, we exclude annual time series where any of the time steps has? a velocity <0.3 m/d or the estimated uncertainty exceeds >20% of the observed velocity.Only complete annual time series are included in the analysis.Due to these filtering steps, the number of time series included varies from year to year.We normalize each annual time series by its maximum value that year.This is done to focus on the annual pattern rather than the velocity magnitudes.No further filtering is performed and our results show the dominating patterns resulting from both long and short term trends for the years [2017,2021].Following these preprocessing steps, we have created a data set consisting of normalized time series from all 5 years, in total 632,835 time series.
We apply K-Means clustering (Pedregosa et al., 2011) on this ensemble of annual time series of ice velocity to identify characteristic Greenland wide seasonal patterns of flow.The K-means algorithm groups time series with similar signal shapes into a specified number of clusters by minimizing the within-cluster variance.K-Means is an unsupervised machine learning algorithm, that is, there is no information provided to guide the algorithm.The algorithm follows these main steps (Pedregosa et al., 2011): 1. Initial centroids of the clusters are chosen with the K-means++ function, which distributes the initial centroids at equal distances instead of the random initialization used in basic K-means.2. All the time-series are sorted to the closest centroid according to the Euclidean distance to each centroid.3. New centroids are derived as the mean of the time-series assigned to each of the current centroids.
By iterating through steps 2 and 3 until the centroids remain constant within a threshold, the algorithm minimizes the within-cluster variance (or inertia) given by where μ j is the mean of the observations in cluster j, C is the set of clusters and x i is time series i.The process is repeated several times varying the initial centroids, and the final centroids are chosen as the run with the lowest inertia (Equation 1).
The K-means algorithm requires the number of clusters, k, to be specified a priori.A reasonable value of k can be determined using a so-called "elbow plot" showing the inertia (Equation 1) of the final centroids as a function of k (Figure S1 in Supporting Information S1).A good choice of k is when the inertia decreases linearly with increasing k.
The K-means algorithm was run for k ∈ [2, 14] and the resulting elbow plot shows the sum of the squared distance from each time series to the nearest cluster center (Figure S1 We compare the annual time series of ice velocity to annual time series of estimated runoff from the RCM RACMO2.3p2, in order to investigate the link between seasonal flow and runoff.RACMO2.3p2 is run at 5.5 km spatial resolution for the period 1958-2021, over a domain including the GrIS, its peripheral ice caps, and surrounding Arctic glaciers (Noël et al., 2019).The daily runoff product used here was further statistically down scaled from RACMO2.3p2 at 5.5 km onto a 1 km grid (Noël et al., 2019).Detailed model description and evaluation using in-situ and remote sensing measurements are available in Noël et al. (2018) and Noël et al. (2019).Due to the difference in spatial resolution and grid, the ice velocity is compared to the runoff time series in the nearest grid point.This direct comparison assumes that runoff immediately enters the subglacial hydrological system.Access to the bed varies across the ice sheet and thus adds to the uncertainty in the timing of the peak runoff accessing the bed.Each annual runoff time series is normalized by the maximum value that year.

Results
The K-means algorithm groups the time series into overarching patterns.For k = 4 (Figure 2), cluster 2 is characterized by an early summer peak in velocity that slows down to an annual minimum in mid-summer/mid-runoff season, followed by a speed-up over winter and spring.Cluster 3 time series has no early summer peak but slows down over summer.An annual minimum is reached at the end of summer followed by an increase in velocities over the late summer, fall and early winter.Cluster 1 and 4 both have a mid summer peak in velocity followed by a slow down.For cluster 4 the amplitude is small on a background of accelerating flow while the velocities have a significantly larger amplitude and no trend in the case of cluster 1.The shape of the speed-up in clusters 1 and 4 largely matches the shape of the runoff curve, while the late summer minimum for cluster 2 and 3 coincides with In all 6 basins, 2019 was the year with the highest amount of runoff (Figure 3).In the northern basins (NO, NW, and NE) there is a tendency for the number of cluster 2 time series to correlate with runoff.Thus in years with higher runoff larger areas have a slow-down over summer.In contrast, the distribution of clusters in the southern basins (SE and SW) is largely insensitive to the variability in the amount of runoff over the study period.
Three things are worth noting for increasing k.First, as expected from the elbow plot (Figure S1 in Supporting Information S1) the standard deviation bands around the mean clusters become narrower.Second, as k increases, the clusters branch out as more specialized versions with different timings or amplitudes or annual trends.For example, the two clusters for k = 2 with cluster 1 having a slow down over the runoff season reaching an annual low late in the season and an increase in velocity during winter, and with cluster 2 having a velocity variation which is synchronous with melt, are present for all k = [2, 7] with narrower standard deviation bands as k increases (see Figures S2-S6 in Supporting Information S1).Finally, k = 6 includes a cluster with minimal seasonal change (Figure S5 in Supporting Information S1); this cluster includes ∼17% of the time series.For k < 6, these time series are not identified separately, but contribute to the standard deviation of the other clusters.

Discussion
Using an unsupervised clustering algorithm, we are able to identify the general seasonal patterns of flow governing large areas of the GrIS.Independent of the number of clusters, k, that the glaciers are distributed between, the clusters share the characteristics of the two types identified previously by Moon et al. (2014) related to runoff and the development of the subglacial hydrological system.If we only sort our data points into two clusters (k = 2) (Figure S2 in Supporting Information S1), cluster 1 can be recognized as type 3 flow and cluster 2 as type 2 flow, but the standard deviation bands are wide as all time series are forced into the two clusters.In the case of four clusters, both cluster 2 and 3 show type 3 behavior, while cluster 1 and 4 have the characteristics of type 2. Our analysis thus statistically substantiate, that the general patterns identified by Moon et al. (2014) and Vijay et al. (2019Vijay et al. ( , 2021) ) are the dominant patterns of flow throughout the faster flowing areas of the GrIS.
In agreement with Moon et al. (2014) and Vijay et al. (2019), we find that the slow down during summer in cluster 2 and 3 (for k = 4) largely coincides with peak runoff and that the increase in velocity from the seasonal low starts at the end of the runoff season.Likewise, for cluster 1 and 4 the acceleration of flow occurs largely at the same time as the onset of runoff, while deceleration starts as the runoff season declines.Moon et al. (2014) and Vijay et al. (2019) link these differences in seasonal behavior to the state of the subglacial system.The slow down during high runoff to a seasonal minimum is indicative of a switch from an inefficient subglacial drainage system to an efficient system, while the system remains inefficient or limited in the case of clusters 3 and 4. The development of the subglacial drainage system over the runoff season depends on the timing and duration of the melt, the volume of runoff reaching the bed, basal melt rates, geology, topography, ocean interaction and internal instabilities (e.g., Bartholomew et al., 2010;Chu et al., 2016;Doyle et al., 2018;Schoof, 2010;Zwally et al., 2002).Variability in the distribution of seasonal flow characteristics is therefore expected -also between neighboring glaciers -and has been observed previously in studies of near terminus flow (Davison et al., 2020;Moon et al., 2014;Vijay et al., 2019Vijay et al., , 2021)), of specific regions (Larsen et al., 2016;Lemos et al., 2018;Rathmann et al., 2017;Sakakibara & Sugiyama, 2020) or predicted by modeling (Hewitt, 2013).Here, we provide the first GrIS wide temporal and spatial overview.
In our analysis all 632,835 time series are assigned a cluster.A closer look at the individual time series (e.g., Figure S12 and S13 in Supporting Information S1) within each cluster reveals that they are all variations of, but rarely identical to the cluster centroid, and are assigned that particular cluster as the closest match.This adds significant noise to each cluster especially for low values of k.Moon et al. (2014) notes that the patterns identified in their study should be viewed as end members of the response to runoff.Thus not all types of flow are distinguishable by the four clusters and all time series within each cluster may not result from the same combination of processes.Previous studies have shown that Jakobshavn Isbrae and Helheim Glacier mainly are controlled by processes at the marine terminus (Ashmore et al., 2022;G. Cheng et al., 2022;Joughin et al., 2008).However, in their analyses of Jakobshavn Isbrae, Joughin et al. (2012) and Ashmore et al. (2022) found that other processes contributes as well but to a smaller degree.For Helheim Glacier in southeast Greenland, Ultee et al. (2022) concluded that both terminus position and runoff are important controls on the seasonal variations in flow.This is most likely the case for most glaciers and must be taken into account when interpreting the clusters.For example, the inland area of Petermann Glacier and the downstream area of Jakobshavn Glacier are both assigned a mixture of cluster 3 and 4. However it is expected that the processes driving the variations are different at the two locations due to the different amounts of runoff, proximity to the front and the magnitude of the flow speed.Figures S15-S20 in Supporting Information S1 in the Supplementary give information on absolute flow speeds and the amplitude of the seasonal variations.
Some regions or glaciers are clearly dominated by one type of seasonal pattern while others vary between neighboring glaciers or along individual glaciers.For the large marine outlets in north Greenland, the prevailing seasonality is cluster 1 (type 2), not only near the termini as was also found by (Vijay et al., 2021) but over large parts of the downstream areas.In southern Greenland, many of the large outlets are also mainly dominated by one type, in this case cluster 2 (type 3).
Our analysis shows that especially for the northern basins (NO, NW, and NE; Figure 1) the volume of runoff generally correlates with the number of cluster 2 (type 3) time series (Figure 3).We note an increase in the number of cluster 2 (type 3) glaciers in 2019, which was a high melt year.This increase was also noted by Vijay et al. (2021).During our study period, especially the NW basin is sensitive to runoff.In the NO basin, Ryder Glacier (Figure 1) is a large contributor to the increase in cluster 2 grid cells switching from cluster 1 (type 2) in 2017 and 2018 to cluster 2 (type-3) in 2019.Also areas of Humboldt Glacier and Academy Glacier switch to cluster 2 in 2019 indicative of a more efficient drainage system.In 2020, Ryder Glacier is a mixture of cluster 2 and 1, but in 2021 it is back to cluster 2 like in the high melt year 2019.Conversely, the number of time series classified as cluster 2 in the SW basin is insensitive to available runoff over the study period although 2019 also stands out is a high melt year in this region.
In land terminating areas, it has been both observed and predicted by modeling (e.g., Davison et al., 2019;Hewitt, 2013;van de Wal et al., 2015) that near the margin, melt discharge to the bed is sufficient for developing and sustaining of an efficient drainage system.This will lead to a late season minimum (typically cluster 2).In areas further away from the margin, where the melt season is shorter, runoff less abundant, and an efficient system does not develop, water pressures are thought to remain high in this weakly connected system throughout the melt season (Andrews et al., 2014).In these areas a late seasonal minimum is not observed (typically cluster 1).The interannual response of the seasonal cycle to variations in melt water is, however, less established.Williams et al. (2020) studied the response of flow at the western Greenland margins to multi-annual trends in available runoff and observed a slow down in reaction to increased melt.In years with higher melt, the efficient system is sustained for a longer time and expands further inland draining water from the weakly connected system, thereby lowering the water pressure in these areas and thus decreasing velocities.In years with lower melt, the weakly connected system is not drained as much by the low pressure in the efficient system and can recharge from basal melting thereby increasing the water pressure inland.Hewitt (2013) phrases this as the distributed system having a long-term 'memory' compared to the efficient, channelized system in the sense that its return to its initial state depends on the level and extent of the drainage.An example of this is Ryder Glacier where the subglacial hydrological system does not necessarily return to the same state immediately after a high runoff year (Figure 1 and Figure S10 in Supporting Information S1).Thus, the impact of the sub-glacial 'memory' must be explored further.Conversely, large areas in southern Greenland already develop an efficient system at the average level of runoff (cluster 2/type 3) and thus the seasonality in these areas is less sensitive to further runoff.
A caveat in our analysis is the varying number of time series each year due the constraint that only complete annual time series are included.An example of this is the large fraction of the total number of time series constituted by cluster 3 and 4 (k = 4) in 2017, 2018, and 2021 in the NE basin.This is due a to a large increase in the number of time series from the inland part included those years and not due to changes in flow.The downstream areas of the glaciers in this basin are dominated by cluster 1 as in the NO basin.Where our analysis includes times series near the glacier fronts we see good agreement between our results and those of Vijay et al. (2021) for the years (2017)(2018)(2019) where the studies overlap.There is also agreement between our results and those of Davison et al. (2020) who studied the evolution of the sub-glacial drainage system near the front of three glaciers in southwest Greenland for the years 2016-2018 (see Figure 1 for location).In that study, Kangiata Nunaata Sermia (KNS) was classified as cluster 1 or 4 (type 2) while its neighbors, Narsap Sermia (NS) and Akullersuup Sermia (AS), exhibited cluster 2 (type 3) behavior.The glacier front of NS is not included in our analysis, but our results align with the findings for KNS and AS.However, we see that the cluster 1 or 4 (type 2) behavior of KNS is only near the front, the upstream part of the glacier continuing far inland is cluster 2 (type 3).The exception is 2017, where a large part of the glacier is accelerating and is classified as cluster 4 (type 2).

Conclusions
We have applied K-means, an unsupervised clustering algorithm, to group annual time-series of ice velocity by seasonal flow patterns.This provides us with a temporal and spatial overview of the location of similar types of flow over all the fast flowing areas of the GrIS -not only near the outlet glacier termini.The patterns identified by the algorithm share the same characteristics as those identified by Moon et al. (2014) near the terminus of marine terminating glaciers and can be related to the state of the subglacial system.We show here that the characteristics documented in Moon et al. (2014) are the prevailing patterns across the entire fast-flowing part of the GrIS.
The results for four different clusters provide a spatial overview of the most common seasonalities in flow across the ice sheet.Furthermore, our results show that for an increasing number of clusters, k, the emerging clusters are special cases or combinations of those for lower values of k.In future studies, this could be investigated further and exploited to give a more complete picture of seasonal flow patterns and drainage system evolution.
In the northern basins, we see a positive correlation between the number of time series each year with an efficient drainage system classified as cluster 2 (type 3) and variations in runoff, while the distribution of clusters in the southern basins is largely insensitive to these variations.We suggest, that this is due to the differences in average runoff.However, a closer look reveals that there is no direct link everywhere in the north for example, at Ryder Glacier.Thus, the "memory" of the sub-glacial drainage system needs further investigation.
We show here that it is possible to obtain insights to the interannual variability of seasonal flow patterns in satellite-derived velocities by using the unsupervised clustering algorithm, K-means.This can be exploited in future studies exploring the coupling between ice flow and external forcings on multi-seasonal scales.

Figure 1 .
Figure 1.(a) Overview of the Greenland Ice Sheet and a view of the spatial distribution of the clusters in 2019.Sectors of the ice sheet is indicated in full black lines, glacier catchments in thin gray.Specific glacier catchments discussed in the text: HuG: Humboldt Glacier, AG: Academy Glacier, RG: Ryder Glacier, DJG: Daugaard-Jensen Glacier, HG: Helheim Glacier, KNS: Kangiata Nunaata Sermia, AS: Akullersuup Sermia, NS: Narsap Sermia, JI: Jakobshavn Isbrae and KG: Kangilernata Sermia.(b-d) are zoom-ins on the three areas marked in (a) and show the cluster distribution in 2018 and 2019.The legend applies to all subfigures.
) in Supporting Information S1.The inertia decreases steeply between k = 2 and k = 4, and decreases linearly in the range k = [4, 5, 6], indicating that more of the variability in the data is explained by the clusters for k = 4 than for k = 2, but less is gained for k > 4. We thus base our analysis on k = 4, but include results for k = [2, 7] in the discussion.Figures S2-S6 in Supporting Information S1 show how the centroids and standard deviations change with increasing k.
the end of the runoff season.In general, clusters 2 and 4 contain more fast flowing glaciers compared to clusters 1 and 3 (FigureS14in Supporting Information S1).The distribution of the clusters reveals a large spatial and interannual variability (Figure1and Figures S7-S11 in Supporting Information S1 for k = 4).We use the basins and glacier catchments defined by Mouginot and Rignot (2019) (Figure1).Most basins and glaciers are characterized by a spatial and temporal variability of belonging to several clusters.The large outlets in north Greenland, however, generally exhibit one type of flow on their fast flowing parts -usually cluster 1. Examples of how the seasonal variability can vary along a glacier are shown in FiguresS12 and S13in Supporting Information S1 of ice velocity in points along a flow line for Kangilernata Sermia in west Greenland and Daugaard Jensen Glacier in east Greenland (for locations see Figure1).Calving front positions from TermPicks and CALFIN (D.Cheng et al., 2021;Goliber & Black, 2021) are also shown in the figures.Kangilernata Sermia is dominated mainly by cluster 2 downstream and cluster 1 upstream, while for Daugaard Jensen Glacier the prevailing patterns are cluster 3 and 4.

Figure 2 .
Figure 2. k = 4: The mean and standard deviation of the clusters shown together with the mean of the normalized runoff in the nearest grid point of each time series included in the specific cluster.

Figure 3 .
Figure 3. Barplots of the annual distribution of clusters for each of the sectors displayed in Figure 1.The total runoff for each year in each sector is indicated in black.