Seasonal Observations at 79°N Glacier (Greenland) From Remote Sensing and in situ Measurements

This study investigates seasonal ice dynamics of Nioghalvfjerdsfjorden or 79°N Glacier, one of the major outlet glaciers of the North East Greenland Ice Stream. Based on remote sensing data and in-situ GPS measurements we show that surface melt water is quickly routed to the ice-bed interface with a direct response on ice velocities measured at the surface. From the temporally highly resolved GPS time series we found summer peak velocities of up to 22% faster than their winter baseline. These average out to 9% above winter velocities when relying on temporally lower resolved velocity estimates from TerraSAR-X intensity offset tracking. From our GPS time series we also found short term ice acceleration after the melt season. By utilizing optical satellite imagery and interferometrically derived digital elevation models we were able to link the post melt season speed-up to a rapid lake drainage event (<24 h) with an estimated drainage volume of 28 × 106 m3. We further highlight that GPS measurements are needed to resolve short term velocity fluctuations with low amplitudes, whereas remote sensing estimates are rather useful for the calculation of general trends in velocity behavior.


INTRODUCTION
Each summer substantial amounts of melt water reach the base of the Greenland Ice Sheet through existing faults such as crevasses and moulins or through hydro-fracturing beneath supraglacial lakes. Once the water reaches the ice-bed interface it can act as a lubricant, seasonally enhancing ice flow in the ablation area (e.g., Zwally et al., 2002;Hoffman et al., 2011;Joughin et al., 2013;Tedesco et al., 2013). Typically, the increase in ice velocity is most pronounced at the beginning of the melt season due to an inefficient subglacial drainage system which is highly sensitive to an increase in water pressure (Bartholomew et al., 2010). During the melt season a more efficient channelized basal drainage system is developed resulting in a deceleration of ice velocities toward the end of the melt season. This general seasonal velocity pattern is particularly observed at glaciers in the south-eastern and western part of Greenland (Bartholomew et al., 2010;Moon et al., 2014;Doyle et al., 2015). In contrast, several glaciers, mostly located in the northern part of Greenland show a direct and strong correspondence between ice velocity and melt water availability suggesting that the development of an efficient basal drainage system is limited or shifted for these glaciers (Moon et al., 2014;Rathmann et al., 2017;Vijay et al., 2019).
So far, most studies investigating seasonal variations of surface velocities are concentrated to the western part of Greenland where supraglacial lakes form and drain frequently and surface velocities show the typical seasonal pattern described above (e.g., Zwally et al., 2002;Bartholomew et al., 2010;Joughin et al., 2013). However, it has been shown previously that the largest increase of supraglacial lakes will likely take place in north-eastern Greenland (Igneczi et al., 2016). Until now little is known about the influences of supraglacial lakes on ice dynamics in this area.
In this study we focus on Nioghalvfjerdsfjorden or 79 • N Glacier, located in north-eastern Greenland (Figure 1). Draining approximately 6% of the Greenland Ice Sheet by area (Krieger et al., 2020), 79 • N Glacier is one of the major outlet glaciers of the North East Greenland Ice Stream, a potentially large contributor to future sea level rise. In order to further investigate the seasonal behavior of 79 • N Glacier we combine remote sensing observations with in situ data acquired in 2017. Based on radar backscatter data from the Sentinel-1 satellite mission we link surface melt events to variations in glacier velocities. The latter are derived from in situ GPS measurements accompanied by estimates from intensity offset tracking on high resolution TerraSAR-X data. We further investigate the effect of a rapid supraglacial lake drainage event on the velocity dataset and give estimates on the timing and amount of water released during the lake drainage.

SEASONAL MELT PATTERNS
To obtain information as to when the melt season starts, for how long it lasts and to what altitude surface melt can be expected we extracted a dense time series of backscatter values from the Sentinel-1 archive. Sentinel-1A was launched on 3 April 2014, followed by its twin satellite Sentinel-1B on 25 April 2016. Both satellites are equipped with a C-band Synthetic Aperture Radar (SAR) sensor with a shifted repeat cycle of 12 days, resulting in a 6 day revisiting time when both satellites are employed. In this study we used numerous Level-1 Single Look Complex (SLC) products obtained in Interferometric Wide swath (IW) mode. The data were acquired between January 2015 and January 2019 in relative orbit 74 over the study area. Following the suggestions of Winsvold et al. (2018) we build our analysis on a time series of the terrain corrected backscatter coefficient gamma nought (γ 0 ). The processing includes calibrating, mosaicing and multi-looking of the TOPS SAR data followed by a terrain correction and geocoding step. Radiometric calibration was done within the GAMMA R software resulting in σ 0 backscatter values using the ellipsoid area as reference area for normalization. Subsequently all bursts of each acquisition were mosaiced into a continuous SLC image. In order to reduce noise, multi-looking was performed with 10 looks in range direction and 2 looks in azimuth direction. To obtain terrain corrected γ 0 values we applied an additional reference area normalization based on the ArcticDEM, gridded to 100 m spatial resolution (Small, 2011;Porter et al., 2018). Finally, all SAR images were geocoded and linear backscatter values were translated to decibels.
Backscatter values returned from glaciers and ice sheets heavily depend on instrument related parameters such as incidence angle and polarization as well as on characteristics of the upper surface layers such as roughness, water content, snow density, grain size and impurities (e.g., Rau and Braun, 2002;Winsvold et al., 2018). As wet snow absorbs most of the SAR signal, low surface scattering prevails in these areas, allowing to discriminate between wet and frozen zones on a glacier.
In order to estimate the onset and duration of surface melt we employed the binned backscatter values (γ 0 ) from the time series shown in Figure 2b. Following QSCAT scatterometer studies over Greenland (e.g., Wang et al., 2007) we compared all γ 0 values to the mean γ 0 of the respective altitude band for December, January and February of the previous winter (Wγ 0 ). If γ 0 drops below Wγ 0 − 3dB the respective bin is treated as surface melt and is marked by a black dot in  (Figure 2b and Video S1). In 2018 surface melt starts on 15 June and ends on 7 September (84 days). Similar to 2017 the maximum extent of surface melt is found on 2 August, while isolated patches are found at low elevations until September 2018 (Figure 2b and Video S1). Due to the 6-day repeat pass of Sentinel-1 these estimates might be underestimated by up to 10 days in the most unfavorable case.
As surface melt heavily depends on air temperature, we additionally employed temperature data from nearby Automatic Weather Stations (AWS) as well as from the recently released ERA5 reanalysis dataset (Hersbach and Dee, 2016). Here we used daily averages of near surface air temperature acquired by AWS KPC_L and KPC_U of the Programme for Monitoring the Greenland Ice Sheet (PROMICE) (Fausto and van As, 2019). KPC_L and KPC_U are located at an altitude of 370 and 870 m, respectively and their positions are marked in Figure 1. At the position of KPC_U we also extracted 2 m air temperature data from the ERA5 reanalysis dataset (Figure 2a). ERA5 reanalysis data are available at a spatial resolution of 0.25 • and it is worthwhile noting that PROMICE observations are not assimilated in the reanalysis dataset (Delhasse et al., 2020). Between ERA5 2 m air temperature data and KPC_U near surface air temperature data we found a squared correlation coefficient of 0.98 when investigating the time interval shown in Figure 2a ( Figure S1).

SURFACE VELOCITIES
Horizontal surface velocities were obtained by GPS measurements and intensity offset tracking on high resolution TerraSAR-X stripmap data. In summer 2017 four autonomous GPS stations were installed near the hinge zone of 79 • N Glacier (Figure 1). All GPS stations included one state-of-the-art Novatel FlexPak6 L-Band GPS receiver, one Novatel choke ring GPS antenna, at least one data logger, two solar panels, one solar charge controller and two 105 Ah batteries. The antennas were mounted on aluminum tripods and fixed to the ice by additional weights and long nails. GPS data were logged every FIGURE 1 | Overview of the study area. Hatched polygon marks the area used for the backscatter analysis shown in Figure 2b. Catchment areas are from Krieger et al. (2020). Inlet shows the locations of GPS receivers shown in Figure 3 and Figure S2. KPC_U and KPC_L are the closest Automatic Weather Stations (AWS) from the PROMICE project (Fausto and van As, 2019). In the background is a Sentinel-2 acquisition from 28 July 2017. Hinge zone is determined from ERS-II SAR interferometry, data were acquired in April/May 2011. Upper limit of major crevasses is shown as dashed black line. 15 s and were processed using Precise Point Positioning (PPP) processing implemented in the Waypoint R software including final precise international GNSS Service ephemerides. Velocities were calculated based on horizontal displacements within a 24 h moving window over all measurements. Each velocity estimate was then allocated to the center time of the moving window. Finally, average velocities were calculated for each day to obtain an overview of the general velocity pattern shown in Figures 3b-e. In order to detect short term variations in ice velocities (<1 d) related to e.g., lake drainage we produced an additional dataset using a 4 h moving window with no daily averaging ( Figure S2). The standard deviation of the horizontal displacement is estimated with <5 mm for all stations. This results in a velocity error of up to 0.06 and 0.01 m d −1 for the 4 and 24 h moving window, respectively.
We further obtained velocity fields from TerraSAR-X intensity offset tracking for a similar time period. The processing includes co-registration of SLC data separated by 11 days followed by cross correlating the backscatter intensity of each SLC pair in a predefined moving search window (e.g., Strozzi et al., 2002). Following Rückamp et al. (2019) we chose a window size of 250 m in both range and azimuth direction with a step size of 50 m. Range and azimuth offsets were then projected into a polar stereographic coordinate system assuming surface parallel ice flow. Surface velocity fields were filtered following the three step filtering approach described by Lüttig et al. (2017). This approach is designed to efficiently detect outliers in surface velocity fields and is based on 1. a segmentation filter which is based on the assumption that velocity fields can be divided into segments of continuous ice flow (Rosenau et al., 2015), 2. anomalies to a running median and 3. variations of flow directions within a moving search window. Finally, all velocity fields were gridded to 250 m and small data gaps were interpolated with an inverse distance weighting scheme. The latter step was applied since surface velocities were calculated during the melt season where the presence of water strongly modifies the radar backscatter resulting in rather noisy velocity estimates and data gaps. However, as our GPS receivers are located in the central part of the ice stream we expect only small variations in the ice velocities of adjacent flowlines, suggesting a lower sensitivity toward grid resolution when compared to e.g., shear margins . As a measure of confidence we included Overall the TerraSAR-X derived velocities are in good agreement with the GPS derived velocities (Figures 3b-e). By combining both datasets we found that surface velocities at the four GPS sites start to increase from their winter baseline in early June and are peaking at the end of July / beginning of August 2017 (annotated as summer speed-up in Figure 3b). However, velocity estimates derived from GPS data are resolved to a higher temporal resolution than those from intensity offset tracking (1 vs. 11 d). Therefore, a greater amount of detail can be observed in the GPS derived dataset. While all GPS receivers show highest surface velocities on 1 August 2017 this peak is averaged out in the TerraSAR-X estimates leading to slower summer velocities. Depending on the dataset and the location, summer acceleration can be up to 22% for the GPS and 9% for the TerraSAR-X derived estimates, respectively. Since we lack winter velocities from the GPS stations due to power limitations, all relative velocity changes throughout this study are referenced to the corresponding pixel values of a TerraSAR-X winter velocity field with data acquired on 10 and 21 December 2017 (not shown in Figure 3b-e). In order to detect the upper spatial limit of seasonally enhanced flow velocities we calculated differences between all TerraSAR-X velocity fields and the winter 2017 reference field Figure 6). If possible we also computed GPS derived velocity estimates specific to each 11-day time period of the TerraSAR-X estimates. Here we found similar values as from the intensity offset tracking with a summer acceleration of up to 7% when compared to the December TerraSAR-X velocity field.
After returning to their approximate winter baseline in late August all GPS receivers show a simultaneous speed-up in late September which is significantly lower in the record of GPS 7 (annotated as lake drainage in Figure 3b, a close-up is shown in Figure S2). Directly after this second velocity peak, velocities remain under their winter baseline for several days, until they slowly start to increase again.

LAKE DRAINAGE
To explain the unusual speed-up observed in September 2017, we carefully checked radar and optical satellite acquisitions in this time period. Between 18 and 19 September 2017 we found a sudden drainage of supraglacial lake 11 (name originates from arbitrary indexing of supraglacial lakes in the area, lake position is at −22.978 • E, 79.338 • N, Figure 1) on cloud free Sentinel-2 imagery (Figures 4b,c). In order to obtain an estimate of the drainage volume we further generated two Digital Elevation Models (DEMs) from bistatic TanDEM-X data acquired on 13 and 24 September 2017 following the methods described by Neckel et al. (2013). Both DEMs were vertically adjusted to the same reference DEM by removing a linear trend estimated over stable bedrock. Pixel-by-pixel DEM differences were calculated on the same spatial grid with a resolution of 10 m. The elevation differences are shown in Figure 4d and show a subsidence of >50 m in the center of the lake. Finally, elevation changes were Frontiers in Earth Science | www.frontiersin.org FIGURE 3 | Relation between near surface air temperature measured at PROMICE AWS KPC_U and KPC_L (a) (Fausto and van As, 2019), horizontal surface velocities as measured at the four GPS sites indicated in Figure 1 (b-e) and surface melt estimates from the Sentinel-1 time series shown in Figure 2b (f). TerraSAR-X derived surface velocities are centered at the mean date between data acquisitions with light blue bars in x-direction indicating the 11-day period between the dates of data acquisition. Light blue bars in y-direction represent the weighted standard deviation of the pixel values employed in the interpolation of the shown velocity estimates. GPS estimates are averaged over one day. Note the different data ranges in (b-e). Scatter plots (g-j) show the relation between GPS measured surface velocities relative to a TerraSAR-X derived winter velocity field and the maximum elevation of surface melt from the Sentinel-1 time series shown in Figure 2b and (f).
Frontiers in Earth Science | www.frontiersin.org translated into drainage volume by where h lake are the elevation changes within manually digitized shorelines, and p is the pixel size which was set to 100 m to reduce noise. In order to solely employ h lake values of elevation loss, we introduced a threshold of -1 m in the calculation of V drain . Following Equation 1 we estimated a volume of 28×10 6 m 3 for this drainage event. Almost 2 month before on 21 July 2017 we installed a data logger at this specific lake with the intention to autonomously record water pressure and temperature. Unfortunately, the data logger could not be recovered in 2018 but was probably lost during the rapid drainage event in September 2017. However, when installing the logger we measured a lake depth of 10.80 m at the location marked by the yellow dot in Figure 4a ( Figure S2). This value fits well to the TanDEM-X derived elevation difference of -11.05 m at the exact location.

DISCUSSION
By analyzing a dense time series of Sentinel-1 γ 0 backscatter coefficients in relation to surface elevation we obtained a good overview of when the melt season starts and to what altitude surface melt can be expected (Figure 2b). Also a clear correlation between surface melt and near surface air temperature is found when including additional meteorological datasets in the analysis (Figure 2a). Focusing on the time with surface velocity measurements, we find that in June, July and August 2017 surface melt is detected generally to elevations between 600 and 1,000 m with short term excursions reaching elevations >1,250 m. During the melt onset on 2 June 2017 surface melt reaches rather low elevations of 580 m, which is also reflected in the AWS data shown in Figure 3a. While the higher elevated KPC_U station measured the first positive degree day on 10 June 2017, i.e., 8 days after the detected melt onset, KPC_L recorded the first positive degree day contemporaneous on 1 June 2017. Having in mind that our GPS receivers are located at elevations between 70 and 390 m we conclude that melt water was available during the entire 3 month period in this region. A general increase in surface velocities is observed at all four GPS stations prior to August 2017 (Figures 3b-e). Also from our backscatter analysis we detected the largest extent of surface melt at this time reaching elevations of almost 1750 m on 1 August 2017 (Figures 2b, 3f). The temporal coexistence of peaks in surface velocities of all GPS receivers to this date is quite striking and suggests a strong and direct reaction of surface velocities to the availability of melt water. At this time also the highest daily mean near-surface air temperatures of +3.1 and +6.1 • C were measured at KPC_U and KPC_L, respectively, underpining that surface melt is largely driven by air temperature (Figure 3a). This is further supported by ERA5 reanalysis data where a significant anomaly is found when comparing the average 2 m air temperature of 1 August 2017 to the monthly average of July 2017 (Figure 5a, Hersbach and Dee, 2016). Also the 0 • C isotherm on 1 August 2017 matches quite well with the upper limit of surface melt as detected from the backscatter time series (Figure 5a). When investigating the scatter plots shown in Figures 3g,h it is likely that the larger the area affected by surface melt is the stronger the velocity response downstream the glacier. However, the question remains how the melt water is speeding up the glacier. In line with recent work by Rathmann et al. (2017) and Vijay et al. (2019) who also found a direct velocity response of 79 • N Glacier to the occurrence of surface melt, we suggest that surface melt water reaches the ice-bed interface quite quickly through predefined pathways such as cracks and crevasses. This is supported when investigating the spatial distribution of seasonal glacier speed-up as shown in Figure 6. Here the dashed black line marks the upper limit of major crevasses which was delineated with the help of several TerraSAR-X acquisitions.
With an average altitude of 680 m this upper crack limit lies at the lower limit of the estimated melt extent shown in Figure 2b. Even in the time period of the largest melt extent (Figure 2b) and the highest ice acceleration (Figures 6D-F) significant changes in surface speed stayed well below this line. This suggests that the melt water from higher elevated regions is routed on the ice surface until it reaches the crevassed zone at lower elevations.
Here the water is quickly routed to the ice-bed interface through numerous cracks and crevasses and is acting as a lubricant for the overlaying ice column causing a seasonal acceleration of the ice stream. For the Greenland Ice Sheet such mechanisms were first described by Zwally et al. (2002) but were previously known also for mountain glaciers (e.g., Iken et al., 1983). After the velocity peak in the beginning of August 2017 velocities started to decrease again accompanied by a decrease in melt extent and air temperature (Figures 3a-f). Then a second speed-up was recorded by our GPS receivers in September 2017 where no surface melt can be detected from the backscatter time series and near-surface air temperatures measured at KPC_U and KPC_L stayed well below 0 • C (Figure 3a). This second speed-up can clearly be attributed to a rapid drainage event of lake 11 between 18 and 19 September 2017. While the lake was completely covered by lake ice on 18 September 2017 no ice layer was present on 19 September 2017 and the lake area significantly decreased in less than 24 h (Figure 4b,c). We therefore conclude, that the lake drainage started at some point close to 18 September 2017 and the lake ice broke up due to the drainage event. Considering the short time period of the lake drainage, water-driven fracture propagation seems to be the most likely mechanism to establish a connection to the bed (e.g. Das et al., 2008). While GPS 1, 2, and 8 received a similar signal on 19 September only a minor speed-up was recorded by GPS 7 which is located upstream of lake 11 (Figures 1, 3b-e). This implies that the water, with an estimated volume of 28×10 6 m 3 was not routed through established channels but overwhelmed the subglacial hydraulic system and spread along the ice-bed interface downstream of lake 11. Here ice velocities increased by up to 24% when compared to pre-drainage velocities ( Figure S2). This is in contrast to other studies which report post-drainage ice velocities of 4-10 times larger than before a rapid lake drainage event (e.g., Hoffman et al., 2011;Tedesco et al., 2013). These  (a), lake is mostly ice free, yellow dot marks the location of lake depth measurement on 21 July 2017 ( Figure S2). Sentinel-2 acquisition on 18 September 2017 at 15:49 is shown in (b), lake is completely covered by lake ice. Additional pre-drainage shorelines are shown. Note the difference in illumination between August and September. Sentinel-2 acquisition on 19 September 2017 at 15:18 is shown in (c), the lake has drained and is not covered by ice anymore. TanDEM-X derived surface elevation changes between 13 and 24 September 2017 are shown in (d). differences are most likely attributed to the comparatively large distances of our GPS receivers to the lake. While Tedesco et al. (2013) report distances within 2 km of the lake and receiver locations along flowlines through the lake, our GPS receivers are located between 6.5 and 12.5 km from the lake and do not match any flowline through the lake (Figure 1). Further, it has been shown that changes in water supply are key for driving ice accelerations rather than the absolute amount of water input (e.g., Schoof, 2010;Bartholomew et al., 2012). This might explain why the relatively large drainage volume has only a limited effect on the ice dynamics in the area.
When focusing our backscatter time series on the melt season 2018, we find indirect evidence for a cold summer with high levels of precipitation as described elsewhere (Polar Portal Season Report, 2018). Compared to 2017 the melt season started almost 2 weeks later on 15 June 2018. The largest extent of surface melt is found on 2 August 2018 where it reaches a similar altitude as on 1 August 2017 (Figures 2b, 5c). While compared to 1 August 2017 lower values of daily mean near-surface air temperature of +1.4 and +3.2 • C were measured at KPC_U and KPC_L, respectively, precipitation might have enhanced surface melt and/or accumulated in the upper snow layers. At this day no helicopter supported field work was possible due to rainy conditions with low clouds and also ERA5 reanalysis data suggests precipitation in the area (Figure 5d). This is in contrast to 1 August 2017 were dry conditions can be assumed (Figure 5b). So far no GPS measurements are available for 2018 but it has been shown previously, that precipitation can have a direct influence on glacier velocities (Doyle et al., 2015).
The results presented in this study are in general agreement with the findings of Rathmann et al. (2017) who also suggest that the seasonal speed-up of 79 • N Glacier is largely driven by surface melt water lubricating the bed. By conducting modeling experiments, Rathmann et al. (2017) rule out other potential drivers such as seasonally enhanced sliding along the side walls of the floating tongue. For recent ice flow conditions Mayer et al. (2018) asses a strong buttressing effect of the floating ice tongue, especially within the narrower part of Nioghalvfjerd fjord (i.e., for the first ∼30 km downstream the floating part of the ice tongue). This implies a low sensitivity of seasonal ice velocities toward varying sea ice conditions, an assumption which can not be made for the neighboring Zachariae Isstrom. For the latter  Rathmann et al. (2017) find that the break up of the ice mélange coincides with the onset of surface melt upstream the glacier, making a distinction of processes governing the seasonal ice acceleration rather difficult. When investigating 45 Greenlandic glaciers, Vijay et al. (2019) could link the seasonal varying surface velocities of roughly half of them to surface melt induced changes of the subglacial hydrology. For less than a quarter of their study glaciers, Vijay et al. (2019) found a correlation between seasonal ice velocities and terminus changes. Even though surface melt water might also play a role in the seasonal velocity evolution of the latter group, Moon et al. (2014) suggest that fluctuations in terminus position are the primary controlling factor for these glaciers. According to the classification scheme introduced by Moon et al. (2014) and Vijay et al. (2019), 79 • N Glacier fits well within their type 2 behavior by showing a strong correlation between glacier velocity and runoff. This also implies that a seasonal development of an efficient channelized subglacial drainage system is rather limited for 79 • N Glacier.
While the above studies are solely based on remotely sensed surface velocities, we show that by the additional use of GPS measurements short term velocity variations are preserved which would be averaged out if only remote sensing estimates were available. Based on intensity offset tracking on Sentinel-1 data, Lemos et al. (2018) and Rathmann et al. (2017) found a seasonal acceleration of 79 • N Glacier by 10 and 11%, respectively. These estimates match comparably well with our estimate of 9% based on TerraSAR-X data, even though sampling location, time stamp and reference data were different in all three studies. Such velocity estimates might be sufficient for monthly averaged mass balance estimates or general trends in seasonal velocity behavior, but short term velocity fluctuations with rather low amplitudes are difficult to quantify from the current satellite missions with temporal baselines of several days.

CONCLUSIONS
By combining remote sensing observations with in situ GPS measurements we were able to link the occurrence of extensive surface melt to a direct velocity response of 79 • N Glacier. While we found a maximum seasonal velocity increase of 22% from GPS measurements, these average out to 9% when velocities are calculated over longer time intervals such as the 11-day repeat pass of TerraSAR-X. We therefore conclude, that surface velocities derived by remote sensing give a good overview over entire seasons or major speed-up events, but are not able to fully resolve short term velocity fluctuations triggered by e.g., rapid lake drainage events. For the latter we found an example even after the melt season in mid September when temperatures were well below 0 • C. By utilizing digital elevation models from the TanDEM-X satellite mission we were able to estimate a drainage volume of 28×10 6 m 3 which resulted in a minor speed-up of 24% when compared to pre-drainage ice velocities estimated from 3 GPS stations downstream the lake. Compared to studies in western Greenland the observed speed up is a magnitude lower. On the one hand this can be attributed to the large distances of our GPS receivers to the lake, on the other hand this could possibly support findings that not the amount of water is decisive for ice acceleration but rather the timing and the changes in water supply in combination with the effectivity of the hydraulic system.

DATA AVAILABILITY STATEMENT
The datasets generated for this study are available on request to the corresponding author.

AUTHOR CONTRIBUTIONS
NN processed and analyzed the remote sensing data and wrote the initial draft of the manuscript. OZ and VH processed the GPS data. NN, OZ, and DS conducted field work. AH designed the project, lead the field program and supervised the study. All authors continuously discussed the results.

FUNDING
This project has received funding from the European Union's Horizon 2020 research and innovation programme under grant agreement No 689443 via project iCUPE (Integrative and Comprehensive Understanding on Polar Environments). OZ received founding through BMBF Project GROCE (03F0778A).