Articles | Volume 16, issue 2
Research article
03 Mar 2022
Research article |  | 03 Mar 2022

Overestimation and adjustment of Antarctic ice flow velocity fields reconstructed from historical satellite imagery

Rongxing Li, Yuan Cheng, Haotian Cui, Menglian Xia, Xiaohan Yuan, Zhen Li, Shulei Luo, and Gang Qiao

Antarctic ice velocity maps describe the ice flow dynamics of the ice sheet and are one of the primary components used to estimate the Antarctic mass balance and contribution to global sea level changes. In comparison to velocity maps derived from recent satellite images of monthly to weekly time spans, historical maps, from before the 1990s, generally cover longer time spans, e.g., over 10 years, due to the scarce spatial and temporal coverage of earlier satellite image data. We found velocity overestimations (OEs) in such long-span maps that can be mainly attributed to velocity gradients and time span of the images used. In general, they are less significant in slow-flowing grounded regions with low spatial accelerations. Instead, they take effect in places of high ice dynamics, for example, near grounding lines and often in ice shelf fronts. Velocities in these areas are important for estimating ice sheet mass balance and analyzing ice shelf instability. We propose an innovative Lagrangian velocity-based method for OE correction without the use of field observations or additional image data. The method is validated by using a set of ground truth velocity maps for the Totten Glacier and Pine Island Glacier which are produced from high-quality Landsat 8 images from 2013 to 2020. Subsequently, the validated method is applied to a historical velocity map of the David Glacier region from images from 1972–1989 acquired during Landsat 1, 4, and 5 satellite missions. It is demonstrated that velocity overestimations of up to 39 m a−1 for David Glacier and 195 m a−1 for Pine Island Glacier can be effectively corrected. Furthermore, temporal acceleration information, e.g., on basal melting and calving activities, is preserved in the corrected velocity maps and can be used for long-term ice flow dynamics analysis. Our experiment results in the Pine Island Glacier (PIG) show that OEs of a 15-year span can reach up to 1300 m a−1 along the grounding line and cause an overestimated grounding line (GL) flux of 11.5 Gt a−1 if not corrected. The magnitudes of the OEs contained in both velocity and mass balance estimates are significant. When used alongside recent velocity maps of 1990s–2010s, they may lead to underestimated long-term changes for assessment and forecast modeling of the global climate change impact on the Antarctic ice sheet. Therefore, the OEs in the long-span historical maps must be seriously examined and corrected. We recommend that overestimations of more than the velocity mapping uncertainty (1σ) be corrected. This velocity overestimation correction method can be applied to the production of regional and ice-sheet-wide historical velocity maps from long-term satellite images.

1 Introduction

Ice flow velocity fields on the Antarctic ice sheet (AIS) have been mapped by using synthetic-aperture radar (SAR) and optical satellite images to study ice-sheet-wide ice flow dynamics and AIS responses to global climate changes (Rignot et al., 2011a; Gardner et al., 2018; Shen et al., 2018; Greene et al., 2020). An important direct use of such velocity fields is to estimate the ice discharge from the AIS to the Southern Ocean and perform mass balance analyses (Shepherd et al., 2012, 2018; Gardner et al., 2018; Shen et al., 2018; Rignot et al., 2019). Any errors in the reconstructed velocity fields can cause uncertainties in the estimated AIS mass balance and associated contribution to global sea level (GSL) change (Rignot et al., 2011b; Church et al., 2013; DeConto and Pollard, 2016; IPCC, 2019).

The state of the reconstructed velocity fields can be represented by a series of velocity maps that are derived from periodical satellite images. The variations of the available satellite images in both spatial and temporal coverage determine the extent and time span of the velocity maps. While ice-sheet-wide annual velocity maps, e.g., in the Inter-mission Time Series of Land Ice Velocity and Elevation (ITS_LIVE) project (Gardner et al., 2018, 2019), have been produced from recent satellite images (Storey et al., 2014; Lillesand, 2015; Li and Roy, 2017; Sabins and James, 2020), regional velocity maps at a seasonal or monthly scale have been generated from optical and SAR images (e.g., Landsat and Sentinel; Frezzotti et al., 1998, 2000; Nakamura et al., 2010; Zhou et al., 2014; Greene et al., 2017, 2018, 2020; Moon et al., 2021; Wulder et al., 2019). Furthermore, a weekly ice velocity mapping scheme based on multi-mission satellite images was proposed in Altena and Kääb (2017). However, owing to low image quality, large geolocation errors, and low temporal and spatial coverage, satellite images prior to 1990 are generally less available; appropriate images for ice-sheet-wide or large regional velocity mapping with a shorter time span (e.g., seasonal or annual), especially from 1960s and 1970s, are scarce. Therefore, earlier velocity maps have been produced with a limited extent and a longer time span. For example, combinations of first-generation film-based ARGON images of 1963 (Ruffner, 1995; Kim, 2004), early Landsat MSS images of 1970s, and TM images of 1980s (Chander et al., 2009) have been used to create regional velocity maps with a time span ranging from 1 to 23 years (Bindschadler and Scambos, 1991; Bindschadler et al., 1996; Wang et al., 2016; Cheng et al., 2019; Rignot et al., 2019), although the unique case of a 1963 velocity map of the Rayner Glacier from two ARGON stereo image pairs has been presented (Li et al., 2017). Such a long time span is problematic when we use the feature-matching technique for velocity mapping. For example, at time1 a feature, with an initial velocity v0 at the first location, is taken in the first image. The same feature is tracked in the second image taken at time2 after traveling at the velocity v0 and an acceleration a for a time span of Δt (time2–time1). Thus, the velocity v=v0+aΔt increases along with the time span Δt if acceleration a exists. Given a constant acceleration, the velocity can be overestimated if the time span is long or the velocity overestimation is proportional to the time span. With a 10-year span, an overestimation of  70 m a−1 (or  10 %) was estimated for an area near the grounding line of Totten Glacier (Fig. A1) (Chad A. Greene, personal communication, 2020). Berthier et al. (2003) compensated for the overestimations in a Mertz Glacier mapping study with an image span of 11 years by assigning the overestimated velocities to middle points of traveled segments. This technique should have corrected a large portion of the overestimation given the relatively weak spatial velocity gradient along the main trunk of the glacier. Although it has not been brought to further attention in publications, given its nature and magnitude, this velocity overestimation issue should be fully understood, and a comprehensive correction method should be developed so that corrected historical velocity maps can be analyzed alongside modern maps to create a long record of cohesive ice flow dynamics. Furthermore, this capability of building a long record of AIS ice flow dynamics is important for estimation of long-term AIS mass balance and prediction of the future GSL contribution (Rignot et al., 2019).

This study is a part of our efforts to develop an ice velocity map of East Antarctica for 1963 to 1989 (Li et al., 2017; Ye et al., 2017; Cheng et al., 2019). In this paper we prove the existence of ice velocity overestimation in long-term historical velocity maps and present an innovative correction method. The proposed correction method is based on the Lagrangian velocity that can be calculated from the overestimated velocity map itself and thus does not require any field observations or additional satellite images. We used a set of ground truth velocity maps with time spans of 1 to 7 years derived from recent Landsat 8 images of 2013 to 2020 from the Totten Glacier (TG), East Antarctica, to validate the correction method. We then applied the validated method to a historical velocity map in the David Glacier region, East Antarctica, which was produced from images from 1972 to 1989 acquired during satellite missions of Landsat 1, 4, and 5. We show that the 17-year velocity overestimation can be successfully corrected to within the uncertainty (1σ) of the 1-year map. Another experiment at the Pine Island Glacier (PIG), one of the most dynamic glaciers in Antarctica, was carried out to demonstrate that the overestimation correction method can effectively adjust the overestimated velocities as large as 195 m a−1, while preserving the velocity change signature caused by temporal accelerations due to basal melting and calving activities for long-term ice flow dynamics analysis.

2 Methods

2.1 The velocity overestimation issue

We describe an acceleration-induced overestimation using a typical scenario in AIS (Fig. 1a) where ice flow accelerates over a long slope from several glaciers originating from the inland interior, running through the main trunk, and discharging to the ocean (Bamber et al., 2000; Cuffey and Paterson, 2010; Rignot et al., 2011a). In order to quantify the velocity overestimation, for any point Po(xo,yo) on a glacier (Fig. 1a), we define its velocity in two different frameworks, Eulerian and Lagrangian (Chu and Fan, 2014; Chenillat et al., 2015; Altena and Kääb, 2017). First, the velocity in the Eulerian framework (E velocity) is defined as VE=D/dt, the straight-line distance (D) divided by the time span (dt). Note that for simplicity in equation derivation and discussion, we interchange a vector with its scalar; thus, velocity and speed are not strictly distinguished in this paper. Hence a velocity field described by a velocity map is also defined in the Eulerian framework. In reality, the start point of D is measured on the first image and end point on the second image; the two images are taken with a time span of dt apart (Li, 1998; McGlone, 2013). The reconstructed velocity map is stored as velocity components (VxE,VyE) in the x and y directions, from which the velocity field can be reconstructed. The velocity, represented by a velocity map, may indicate an instantaneous or average velocity depending on the time span dt (Cuffey and Paterson, 2010; Rignot et al., 2011a). Second, the velocity in the Lagrangian framework (L velocity) is defined as VL=S/dt, the curved distance traversed along the flow line (S) divided by the time span (dt; Altena and Kääb, 2017). Given a time span and an initial point Po(xo,yo) (Fig. 1a), its L velocity can be determined by tracking the point using a velocity map (Chu and Fan, 2014; Chenillat et al., 2015). Operations in the Lagrangian framework are often performed to estimate advected ice features (Altena and Kääb, 2017) or forecast future events in earth science applications, e.g., mud and debris flow of landslides (Debella-Gilo and Kääb, 2013; Feng et al., 2016), ocean currents at different depths (Glenn et al., 2016; van Sebille et al., 2018), and storm center motion of a typhoon or hurricane (Cram et al., 2007; Euler et al., 2019).

Figure 1Illustration of acceleration-induced overestimation in a velocity map derived from a satellite image pair with a long time span: (a) schematic scene of accelerated ice flow in a glacier–ice-shelf system in AIS, (b) calculation of E velocities at Po(xo,yo) with an equal increment in time span (1 year), and (c) increase in E velocity at the same point Po(xo,yo) as the time span increases – the cause of the acceleration-induced velocity overestimation.


Assume that we use a set of n+1 images (Imagei, i=0,1,n) that are taken with a time interval of 1 year to produce n velocity maps, V0−i, each derived from an image pair (Image0, Imagei). The time span of the maps increases from 1 to n years. These maps are defined in the Eulerian framework (Fig. 1b). For Po(xo,yo) its location after i years, Pi(xi,yi), is determined in Imagei; its Eulerian distance D0−i is measured using both Image0 and Imagei (i=1,2,n). Consequently, the E velocity of a i-year span (Δto−i) at Po(xo,yo) is

(1) V 0 - i E = D 0 - i Δ t 0 - i .

As the time span increases at a fixed rate of 1 year, the traversed straight-line distance D0−i (red lines in Fig. 1b), correspondingly E velocity V0-iE, increases rapidly because of the acceleration over the traverse (Fig. 1c). In principle, every V0-iE (i=1,2,n) value represents the velocity at the same point Po(xo,yo) (Fig. 1a) in these n velocity maps. In the cases where Imagei was not available and thus the maps V0−i (i=1,n-1) were not produced, we only had the map V0−n with the longest span of n years. It is obvious that at Po(xo,yo) its n-year velocity V0-nE is significantly larger than the 1-year velocity V0-1E (Fig. 1c). In general, we define the velocity overestimation of a i-year E velocity as

(2) OE 0 - i = V 0 - i E - V 0 - 1 E .

Here we use a velocity map of a 1-year span as a baseline (overestimation-free) throughout the paper for simplicity, which can be changed for glacier regions of different ice flow dynamics (spatial acceleration, mainly caused by bed topography and slopes). For example the baseline span is 1 year for TG in Experiment 1 and 3 months for PIG in Experiment 3. We require that the overestimation of the baseline map is negligible or smaller than σ (velocity mapping uncertainty).

2.2 Overestimation correction based on Lagrangian velocity

We propose a method for correction of overestimation in a long-term velocity map based on its Lagrangian velocities. The velocity field described by the baseline map is used to calculate trajectories of ice mass and L velocities. The objective is to calculate the overestimation OE0−i in V0-iE (i=2,n) so that these velocities of longer spans are corrected using Eq. (2) and describe the same velocity field as the 1-year V0-1E.

For point Po(xo,yo) in Fig. 1a, the E velocity V0-iE (i=1,2,n) is presented by individual bars in Fig. 2. Correspondingly its L velocity of i-year span can be calculated as (Halliday et al., 2013)

(3) V 0 - i L = S 0 - i L Δ t 0 - i ,

where S0−i is the Lagrangian trajectory (L trajectory) distance (i=1,2,n):

(4) S 0 - i L = t 0 t i V 0 - 1 E ( t ) dt .

Within a baseline time span (e.g., 1 year or shorter), the difference between the E distance and L distance may be considered within the mapping uncertainty (1 σ): V0-1EV0-1L. Beyond that span, the curved trajectory distance is longer than the straight distance (S0-iD0-i, see Fig. 1b); thus, we have V0-iLV0-iE (i=2,3,n). An approximation of the increased difference between V0-iL and V0-iE is presented by the trends of L velocity (blue line) and E velocity (red line), which start from the same velocity at V0-1L (or V0-1E) and reach the maximum difference between V0-nL and V0-nE at the end (Fig. 2).

Figure 2Derivation of equation for overestimation correction using L velocity. Eulerian velocities V0-iE (i=1,2,n) are represented as bars. The red line is the average Eulerian velocity V0-iE of V0-1E and V0-nE. The blue line is the average Lagrangian velocity V0-iL of V0-1L and V0-nL derived from V0-1E. The black line is the average Lagrangian velocity V0-iL(n) of V0-1L(n) and V0-nL(n) derived from V0-nE.​​​​​​​


Premise I. Within a baseline time span (e.g., 1 year or shorter) each segment (from Pi−1 to Pi in Fig. 1b) is relatively short, and the E and L velocity difference is smaller than σ. Furthermore, over the map span of n years (e.g., 5–10 years or longer) the accumulated E and L distances along the entire trajectory do not deviate significantly from each other, so that the maximum velocity difference in the Lagrangian and Eulerian frameworks (end points of blue and red lines in Fig. 2) is limited within a threshold (V0-nL-V0-nEkσ), where k is a constant and σ is the velocity mapping uncertainty.

In reality, the available historical images may only allow us to produce an E velocity map with the longest time span, i.e., V0−n, which leads to the maximum overestimation as defined in Eq. (2). Based on this map V0−n of n-year span, the i-year span L velocity (black line in Fig. 2) is defined as follow:


Consequently, V0-1L(n) (start point of black line in Fig. 2) is set to be equal to V0-nE.

Premise II. Within the time span of n years (e.g., over 5–10 years), the velocity field described by V0-1E and V0-nE does not change significantly, so that the line between V0-nL(n) and V0-1L(n) (black line in Fig. 2) and that between V0-nL and V0-1L (blue line in Fig. 2) are approximately parallel to each other. Accordingly, the difference between their simple averaged accelerations is within a threshold (aVL(n)-aVLkσΔtn​​​​​​​), where k is a constant and σ is the velocity mapping uncertainty, and

(7) a V L ( n ) = V 0 - n L ( n ) - V 0 - 1 L ( n ) Δ t n , a V L = V 0 - n L - V 0 - 1 L Δ t n .

Given the long-span velocity V0-nE we can calculate the n-year L velocity V0-nL(n). Premise II measures the degree of spatial acceleration invariance over the n-year span (parallelity between black and blue lines in Fig. 2). Thus, it is a necessary condition for using the difference of V0-nL(n)-V0-1L(n), or  (V0-nL(n)-V0-nE), from the n-year-span map to substitute for the difference of V0-nL-V0-1E from the 1-year-span map; furthermore, based on Premise I the overestimation can be computed as

(8) OE 0 - n = V 0 - n E - V 0 - 1 E = V 0 - n L - V 0 - 1 E = V 0 - n L ( n ) - V 0 - n E .

Consequently, using the map with the longest time span of n years, we can effectively go back to V0-1E using a correction term defined as correction =V0-nE-V0-nL(n).

In principle, Premise I presents a sufficient condition, V0-nL-V0-nE kσ. In case Premise I does not hold, for example, because of temporal accelerations induced by calving and basal melting activities (Experiment 1 and Experiment 3), we are still able to correct the OE portion induced by spatial acceleration and preserve the uncorrected OE portion induced by temporal acceleration in the residual ε for ice flow dynamics change analysis (Fig. A2).

Overestimation correction theorem. Assume that the necessary condition in Premise II is met; spatial-acceleration-induced overestimations in long-time-span velocities V0-nE can be corrected or reduced using the following correction term, regardless of temporal acceleration:


If Premise I holds, correction -OE0-n; otherwise, |correction|<|OE0−n|, preserving the velocity increases induced by temporal accelerations in the residual term ε (Fig. A2).

Additionally, given n discrete E velocities V0-iE (i=1,2,n; Fig. 2), the acceleration can be estimated using a linear regression model. Its least squares (LS) estimation is (Montgomery et al., 2021)

(11) a ( V E ) = n i V 0 - i E - i V 0 - i E n i 2 - i 2 .

On the other hand, this acceleration can also be approximated by an averaged acceleration using the initial and end velocities:

(12) a ( V E ) = V 0 - n E - V 0 - 1 E Δ t n .

The two acceleration estimates can be compared to determine if the averaged acceleration be used in the event that the intermediate velocities V0-iE (i=2,n-1) are not available. The red line is then used instead of V0-iE (i=1,2,n; Fig. 2).

2.3 Implementation aspects

Trajectory and L velocity computation. The computation of the L trajectory distances S0-iL and S0-iL(n) involves the numerical implementation of the integral of the velocity field from Po(xo,yo) to Pi(xi,yi) along a flow line on the maps of V0−1 and V0−n using Eqs. (4) and (6), respectively. We select an area of the velocity map of 2013 of Pine Island Glacier, West Antarctica (Fig. 3a and b; Gardner et al., 2018), to show the details of integral implementation. The velocity map has a time span of 1 year and a resolution of 240 m. At each grid cell the velocity is stored according to its components (VxE,VyE). We preprocess each component separately by using a moving window smoothing filter to reduce noise. A 10-year trajectory goes through an ice flow profile from P0 at 2560 m a−1 to P10 at 3935 m a−1 with an average acceleration of 138 m a−2. As shown in the enlarged area for the annual trajectory from P2 to P3 (Fig. 3c), the integral is carried out by accumulation of sub-trajectories (between white dots) with a monthly increment. The sub-grid positions of the monthly segments are interpolated for integration of the overall L distance. Finally, we have the result of V0-1E= 2560 m a−1, S0-10L= 34 960 m, and V0-10L= 3496 m a−1, compared to V0-10E= 3491 m a−1. The 10-year overestimation (V0-10E-V0-1E) at P0 is approximately 931 m a−1 (36 %). Examples of overestimation correction in PIG are given in Experiment 2 and the Discussion section.

Denoising and map quality. The quality of the reconstructed velocity field is important for L trajectory calculation and subsequent overestimation correction. Compared to the above Landsat 8 velocity map, velocity maps derived from historical images may have higher uncertainties due to the low image quality that is inherent in the satellite images of early missions (Kim, 2004; Ye et al., 2017). In some cases, there may be even gaps in the velocity map because image features may not be tracked by using the optical flow or feature-matching techniques when ice features disappear through glaciological processes or large advection motion in ice shelf fronts (Scambos et al., 1992; Schenk, 1999; McGlone, 2013). Therefore, preprocessing of the velocity maps should be performed to eliminate outliers; then interpolation may be applied to fill small gaps so that the integral of L trajectories and the calculation of L velocities can be realized. Despite the sub-pixel accuracy of the orthorectification of historical images and the sub-pixel precision of the cross-correlation-based feature-matching method, the reconstructed velocity maps may still be subject to a systematic bias called pixel locking, which results in matched features moving close to integer pixel positions (Shimizu and Okutomi, 2005). This pixel level bias may introduce additional uncertainties to areas of large velocity gradients. The pixel-locking effect can be reduced by using a coarse-to-fine hierarchical feature-matching method as demonstrated in the Antarctic and planetary environments (Debella-Gilo and Kääb, 2011; Li et al., 2011; Heid and Kääb, 2012; Li et al., 2017).

Figure 3Implementation of numerical integral and L velocity calculation in Pine Island Glacier, West Antarctica: (a) black rectangle indicates the area of 10-year trajectory in panel (b), numbered yellow rectangles with crosses are five selected areas in Experiment 3 (Fig. A3), and background is a 1-year velocity map (Gardner et al., 2018); (b) enlarged trajectory area with annual points from P0 to P10; and (c) numerical integral of the third-year sub-trajectory from P2 to P3 with a monthly increment.

3 Results

3.1 Experiment 1: validation of the correction method at Totten Glacier, East Antarctica

3.1.1 Preparation of validation velocity maps

This experiment is designed to validate the proposed velocity overestimation correction method using the velocity map of Totten Glacier (TG) (Fig. 4a; Gardner et al., 2019), which is one of the most active glaciers in East Antarctica and has been experiencing significant mass loss since 1989 (Li et al., 2017; Gardner et al., 2018; Shen et al., 2018). Figure 4b shows an ice flow velocity map of TG from 2013 from ITS_LIVE (Gardner et al., 2019). This 240 m resolution velocity map was derived from the 15 m resolution Landsat 8 images with a 1-year span (January to December 2013). Under the primary forces of gravity, margin shear, basal drag, and others (Bamber et al., 2000; Cuffey and Paterson, 2010; Rignot et al., 2011a), the ice flow accelerates over a path of 130 km along the centerline at a velocity of  600 m a−1 upstream of the grounding line to  2300 m a−1 at the ice shelf terminus.

Figure 4Totten Glacier as an example for the validation of the velocity overestimation correction method: (a) the Totten Glacier region and five areas (white rectangles) selected in different parts of the ice shelf; background is Landsat 8 image of 3 December 2013; black line is grounding line from Rignot et al. (2011c); and (b) velocity map of Totten Glacier of 2013 (Gardner et al., 2018); red and yellow lines are the shelf front of 6 April and 9 September 2015, respectively; light blue line is the PSI boundary (Fürst et al., 2016).

In order to validate the velocity overestimation correction method, we need to determine if premises I and II are fulfilled. This further requires us to have E velocity from the map series of V0−i (i=1,2,n), which are not readily available except the 1-year map of 2013 in Fig. 4b. To avoid lower-quality historical velocity maps that may influence the effectiveness of the validation, we use the earliest available high-quality Landsat 8 images from 2013 to 2020 to produce velocity maps V2013–i (i= 2014, … 2020​​​​​​​). The validation should also consider velocities of areas with different ice flow dynamics. Thus, we select five areas (white rectangles in Fig. 4) to produce smaller maps ( 10 km × 20 km), which contain 7-year trajectories and represent different ice flow dynamic characteristics of the ice shelf, i.e., Area 1 close to the grounding line, areas 2 and 4 along the main trunk, Area 5 on the tributary flow, and Area 3 near a suture zone between the main and tributary glaciers.

Using the matching method adapted for ice surface features (Li et al., 2011, 2017), we produced 35 smaller-sized velocity maps V2013–i (i= 2014, … 2020) for seven time spans and five areas. Information about the Landsat 8 images used and the acquisition dates are given in Table A1. The OE-free velocity maps of 1-year V2013–2014 of all five areas are illustrated in panels 1a–5a of Fig. 5. The corresponding overestimated maps of 7-year span V2013–2020 are shown in panels 1b–5b.

Figure 5Velocities in five areas (rectangles in Fig. 4) of the Totten Glacier are used to validate Premises I and II. Panels (1a)(5a)​​​​​​​ show the reconstructed 1-year velocity maps V2013–2014 in the areas; matched points (red triangles) are used to map E velocity V2013–iE (i= 2014, … , 2020) (red lines in panels 1c5c); points along the flow line (blue dots) are tracked from the 1-year maps and used to calculate L velocity V2013–iL (blue lines in panels 1c5c). Similarly, panels (1b)(5b) illustrate the reconstructed 7-year velocity maps V2013–2020 in the areas with the matched points (red triangles) for E velocity V2013–iE7; points along the flow line (black crosses) are tracked from the 7-year-span velocity map and used to calculate L velocity V2013–iL7 (black lines in panels 1c5c). In each area (panels 1c5c) the difference between red and blue lines increases with time but is limited within kσ at the end (Premise I), except a large k in Area 4 because of the effect of a calving event; the black line is above the blue line due to the spatial-acceleration-induced OE over 7-year span, but they are relatively parallel (Premise II), and thus a correction can be estimated.


The uncertainty of the maps is estimated from two error sources: the geolocation error of the satellite images σRef and the error of feature matching σMatch:

(13) σ Velocity = 1 Δ t σ Ref 2 + σ Match 2 ,

where σRef for Landsat 8 images is 18 m (Storey et al., 2014) and σMatch is set to 8 m,  0.5 pixel (Li, 1998; Li et al., 2017). Thus, the uncertainty of the measured E distance is 20 m. Consequently, the velocity uncertainty σVelocity becomes 20 m a−1 and 3 m a−1 for a 1-year map and a 7-year map, respectively.

In each area, we select an ice feature in the upper stream part in 2013. Its locations in the following 7 years are determined from the satellite images of 2014–2020 using the feature-matching method and shown as red triangles in the maps (Fig. 5). Correspondingly, we track its annual trajectory locations from the 1-year and 7-year maps using Eqs. (4) and (6) and plotted them as blue dots and black crosses, respectively, on the maps. Using these trajectories in the Eulerian and Lagrangian frameworks we calculate V2013–iE, V2013–iL, and V2013–iL7 and plot them as red, blue, and black curves in panels 1c–5c, respectively.

3.1.2 Validation of premises I and II

In general, the L velocities from the 1-year map (blue curves in Fig. 5) are mostly faster than the corresponding E velocities (red curves). To validate Premise I we further examine the L velocity V2013–2020L and E velocity V2013–2020E of the 7-year span, i.e., differences between the end points of blue and red curves in all five areas (panels 1c–5c of Fig. 5). They are less than or equal to 8 m a−1 ((VLVE)2013–2020 in Table 1), except −31 m a−1 in Area 4, which was influenced by a calving event in 2015 as mentioned in the Discussion section. The averaged difference of all five areas is −4± 16 m a−1, that is, within 20 m a−1, the uncertainty (1σ) of the 1-year velocity map (2013–2014) that is considered OE-free. Consequently, the conditions in Premise I are met.

Table 1Velocity and acceleration in Eulerian and Lagrangian frameworks used for validation of the overestimation correction method. “Actual E velocity and OE” lists actually mapped 1-year and 7-year E velocities and their differences as overestimations in all five areas. Premise I contains 7-year L velocities computed from the 1-year velocity map and corresponding L and E velocity differences, which are used for validating Premise I. Premise II illustrates averaged L accelerations computed from the 7-year and 1-year velocity maps, respectively, as well as their differences, which are used for validating Premise II. “Overestimation correction” presents 7-year L velocities computed from the 7-year map, overestimation corrections, and E velocities and residuals (or errors) after correction.

Download Print Version | Download XLSX

For each area we calculated a least squares acceleration estimate a(V2013–iL7) of the 7-year L velocity (black curves in panels 1c–5c of Fig. 5) over all years (i= 2013, 2014, … 2020) using Eqs. (11); then we calculated an average acceleration aV2013–2020L7 using Eq. (12) based only on the velocities of the beginning and end years (2013 and 2020). In all five areas their differences (Table A2) are less than 1.3 m a−2, which is less than σa (3 m a−2), the acceleration uncertainty estimated from the velocity uncertainty of 20 m a−1 (1σ). That means that within the uncertainty of the 1-year velocity map we can substitute the more rigid LS acceleration a represented by a black curve in Fig. 5 with the average acceleration a estimated from the straight line between the beginning and end velocities. Similarly, the differences between accelerations a and a of the 1-year L velocity V2013–iL in all five areas (Table A2) are also less than 1.5 m a−2. Within the uncertainty of the 1-year velocity map we can then substitute a (blue curve in Fig. 5) with a (straight line) for the 1-year L velocity V2013–iL.

Validation of Premise II requires us to examine the accelerations represented by the 7-year L velocity V2013–iL7 (black curves in Fig. 5) and the 1-year L velocity V2013–iL (blue curves). The relatively well maintained parallelity of the black and blue curves in all five areas shows that the velocity trend in the Lagrangian framework has not changed significantly over the 7 years. This is further quantified by the differences between the averaged 7-year and 1-year accelerations aVL(7)-aVL in Table 1. The differences are within ±2.1 m a−2 in all five areas, with an overall difference of -0.8±1.3 m a−2, which are all within σa (3 m a−2). Thus, the conditions in Premise II are also fulfilled.

3.1.3 Overestimation correction

As shown in the first part of Table 1, we compare the overestimated E velocity V2013–2020E in the 7-year map V2013–2020 with V2013–2014E in the 1-year map V2013–2014. Their differences make the actual overestimation OEActual in five areas. The mean overestimation of 34 m a−1 with a range from 3 to 69 m a−1 is significant in comparison to σ (20 m a−1). Thus, the overestimation should be corrected. It should be noted that a higher overestimation is expected for historical velocity maps of longer time spans, e.g., 10–15 years.

Using our correction method in Eqs. (9)–(10), we estimated the correction term (V2013–2020E-V2013–2020L(7)) and applied it to the overestimated V2013–2020E(Table 1). After the overestimation correction, the corrected velocity VCorr.E from the 7-year span V2013–2020E is agreeable to the 1-year span V2013–2014E within an average difference of 14 m a−1 (ε, residual), which is less than σ (20 m a−1). Overall, the overestimations are effectively corrected.

3.2 Experiment 2: velocity overestimation correction in the David Glacier region, East Antarctica

Experiment 2 in the David Glacier region on the Scott Coast, East Antarctica, demonstrates the applicability of the introduced method to correct the overestimation in historical multi-year span velocity maps derived from long-term images from 1972 to 1989. Velocities in this region from 1988 to 1992 were mapped by using GPS and image-feature-tracking techniques (Frezzotti et al., 1998, 2000). A new GPS campaign was carried out in the region during 2005–2006 (Danesi et al., 2008). Velocity changes from 2016 to 2020 were detected using Sentinel-1A SAR images (Moon et al., 2021). In this experiment we produced a velocity map of the region from 64 Landsat images collected from 1972 to 1989 (Table A3). The image pairs used for velocity mapping have mainly four time spans, namely 1, 15, 16, and 17 years, with their footprints illustrated in Fig. 6a. The Landsat images were first orthorectified using ground control points (GCPs), which were selected at outcrops, blue ice, and ice rises (Ye et al., 2017). Then the velocity field was reconstructed in three steps: (a) measurement of a set of manually selected seed points for representing the initial structural velocity information in stationary, low velocity, and dynamic flow areas; (b) velocity point generation by the ice-surface-feature-based matching method controlled by a triangular irregular network (TIN) model and initiated by the seed points; and (c) dense grid matching for generation of the gridded ice velocity field with a mixed time span from 1972 to 1989. This feature and grid image matching method has been developed for surface mapping from optical satellite images in challenging planetary and polar environments and successfully applied to the reconstruction of the ice velocity fields in the Rayner and Fimbul ice shelves in East Antarctica from historical optical satellite images of 1963–1987 (Li et al., 2011, 2017). The produced velocity map has a grid spacing of 500 m after a kriging interpolation from all matched feature and grid points (Fig. 6b). The uncertainty (1σ) of the velocity points is 4, 31, and 62 m a−1 for time spans of 15–17 years, 1 year (TM), and 1 year (MSS), respectively.

Figure 6Application of the overestimation correction method to a historical velocity map of 1972–1989 of the David Glacier region, East Antarctica. (a) Footprints and time spans of the Landsat image pairs used for velocity mapping with a background of Landsat Image Mosaic of Antarctica (Bindschadler et al., 2008), (b) produced velocity map of 1972–1989 and selected points for examination of velocity overestimation and correction, (c) overestimation values of all velocity points generated using the image-feature-matching technique, and (d) velocity points with OE  1σ (red points).

At all velocity points, including 8564 matched features from the feature-matching process and 59 783 grid points from the dense-grid-matching process, we calculated the velocity overestimation values (Fig. 6c). The OEs are generally low on the grounded ice, mostly within ±4 m a−1, which are less than the mapping uncertainties (1σ) of different time spans (Fig. 6d and Table 2). The higher OEs ranging from  4 to  50 m a−1 are on the ice tongue, ice shelves, and glaciers. We selected a total of 16 points from four glaciers and the inland region to examine detailed OEs in different areas of ice flow dynamics (Fig. 6b and c). The OEs at three points in the inland region (IL1–IL3) are negligible (Table 2). Points DG1–DG6 along the main trunk of the David Glacier present an increasing trend of OEs from 1 to 4 m a−1 on the ice tongue up to 36 m a−1 at the ice shelf outlet, as the average Lagrangian acceleration aVL(n) increases from 0.1 to 2.2 m a−2. However, DG7 is located in a relatively steady velocity area and thus has negligible overestimation. Furthermore, the OE values at the selected points on the smaller glaciers, Reeves Glacier (RG1 and RG2), Priestley Glacier (PG1 and PG2), and Campbell Glacier (CG1 and CG2), range from 0 to 12 m a−1 and are generally greater than those on the grounded ice. At the velocity points where the OEs are greater than or equal to the velocity mapping uncertainties (1σ), the corrected velocities are calculated by applying the overestimation corrections.

Table 2Ice flow dynamics parameters, computed OEs, and corrected velocities at 16 selected points.

Download Print Version | Download XLSX

Overall, for a historical velocity map of 1972–1989 with time spans from 1 to 17 years in the David Glacier region, the overestimations and corrections mainly occur on the David Glacier and smaller glaciers where accelerations in the spatial domain (not the temporal domain) exist, while those on the grounded ice are non-significant. We recommend that the overestimations that are greater than or equal to the mapping uncertainties (1σ) be corrected.

3.3 Experiment 3: velocity overestimation correction at Pine Island Glacier, West Antarctica

We further apply the OE correction method at the Pine Island Glacier, one of the most dynamic glaciers in Antarctica. Four areas (1, 2, 3, and 4) along the main trunk and one area (5) in the margin in PIG are selected (Fig. 3a). The longest time span of 7 years (2013–2020) is determined, restricted by the very high velocity (up to  4500 m a−1) on the ice shelf, untraceable surface features over long spans, and availability of the high-quality Landsat 8 images. Hence, most of the ice shelf cannot be covered by a 7-year span velocity map, and the five selected areas are near and upstream from the grounding line. The five areas were all mapped with a baseline span of 3 months (panels 1a–5a in Fig. A3) but different longest spans (n years, n=6, 7, 7, 7, and 3 for the five areas) (panels 1b–5b in Fig. A3) because of difficulties in image feature tracking. Accordingly, in each area one surface feature was selected and tracked for n years to estimate E velocity V2013–iE; its L velocities of V2013–iL and V2013–iL(n) were computed from the 3-month-span and n-year-span maps, respectively. Consequently, the resulting E velocity (red lines), L velocity from the 3-month-span maps (blue lines), and L velocity from the n-span maps (black lines) are illustrated in panels 1c–5c in Fig. A3.

The mapped ground truth data of 3-month- and 7-year-span velocities in five selected areas showed actual OE values in PIG from 18 m a−1 in Area 5 of the inland interior to 626 m a−1 in Area 2 near the grounding line at the glacier outlet (Table 3). Much of these OEs are attributed to the temporal acceleration patterns in panels 1c–5c in Fig. A3 where the tracked E velocities V2013–iE of red lines are, in the later years, consistently above the L velocities V2013–iL of blue lines from the 3-month-span maps of 2013. These temporal accelerations were mainly caused by the drastic calving activities in and after 2017 in PIG and continuous basal melting in the Amundsen Sea sector during the study period (Joughin et al., 2020, 2021). Furthermore, the black and blue lines are approximately parallel (panels 1c–5c in Fig. A3). Hence, the proposed method corrected on average 97 m a−1 of the spatial-acceleration-induced OE, which is  40 % of the overall OE (245 m a−1) in five areas (Table 3). The average residual of 148 m a−1, accounting for another 60 % of the overall OE, represents the OE portion induced by the temporal acceleration in Eq. (9).

Table 3Application of the OE correction method in PIG. “Actual E velocity and OE” includes E velocities of 3-month and n-year spans and their differences as OEs. “Overestimation correction” presents n-year L velocities from n-year-span map, OE corrections, corrected E velocities, and residuals (or errors) after correction.

Download Print Version | Download XLSX

In summary, OEs exist in historical long-span velocity maps of Antarctica, from smaller glaciers, such as David Glacier, to the fast-flowing glacier of PIG, where spatial accelerations are produced by, e.g., bed topography and slopes. In general, they are less significant in slow-flowing grounded regions with low spatial accelerations. Instead, they take effect in places of high ice dynamics, for example, near grounding lines and often in ice shelf fronts. Velocities in these areas are important for estimating ice sheet mass balance and contribution to global sea level changes, as well as for analyzing ice shelf instability. For instance, in the David Glacier the large OE corrections (up to 36 m a−1) occur on the ice shelf (Fig. 6). The OEs of a 7-year span, up to 69 m a−1 (Table 1), would be about  50 % of the velocity increase detected during 1989–2015 in the Totten Glacier (Li et al., 2016). The PIG experiment showed an extreme case in Antarctica where the OEs of a 7-year span can go as large as 626 m a−1 ( 20 %) near the grounding zone (Table 3); furthermore, the OEs of a 15-year span can reach up to 1300 m a−1 along the grounding line and cause an overestimated GL flux of 11.5 Gt a−1 if not corrected (Fig. A4). Therefore, the magnitude of the OEs contained in the long-span historical velocity maps is significant. When overestimated historical maps of 1960s–1980s are used alongside recent maps of 1990s–2010s for assessment of the long-term global climate change impact on the Antarctic ice sheet and for forecast modeling, the overestimated historical states may lead to underestimated long-term changes. Furthermore, compromised forecasting may result. Thus, the OEs in the long-span historical maps must be seriously examined and corrected.

4 Discussion

In this section we discuss the applicability of the proposed method in terms of overestimation-free time span, influence of complex glacier geometry, overestimation in fast flowing glaciers, and comparison with the “midpoint” method.

4.1 Threshold of the overestimation-free time span for trajectory segments

The choice of a short time span ΔtRef (e.g., a few months to 1 year) for the OE-free segments along a trajectory in Premise I makes sure that the difference between the E and L velocities within the span is negligible, or less than σ (velocity mapping uncertainty). It is also the baseline time span of the initial OE-free velocity map. Determination of this threshold has an implication on validation of Premises I and II, as well as the integration period of the trajectory segments from Pi to Pi+1 (Fig. 3b). Estimation of ΔtRef can be performed in each glacier by a linear regression between the E or L velocity V and time span Δt, V=KΔt+b. Given b, K, and σ (20 m a−1), ΔtRef can be calculated as σΔtRef=σ-bK. Our experiment results show ΔtRef as 3.2 years and 3 months for TG and PIG, respectively. Thus, the estimated OE-free time spans appear to be related to ice flow dynamics of the glaciers. In Experiment 1 and Experiment 3, we used 1 year for TG and 3 months for PIG, respectively. We suggest that a linear regression for ΔtRef estimation be performed before extensive historical velocity mapping is carried out.

4.2 Glaciers with complex geometry

Within a long time span (e.g., over 5–10 years) in Premise I the difference between L and E velocities accumulated over all segments along the entire trajectory, i.e., the end-point deviation between the red and blue lines in Fig. 2, is measured with a more tolerable threshold: (V0-nLV0-nE)kσ. Although the OEs of the trajectory segments are controlled by ΔtRef, ice mass moving along a curved flow line over this long span may result in an accumulative discrepancy. In Experiment 1 we showed that the average difference (V0-nLV0-nE) in TG is within 1σ. Our further experiments in five smaller Antarctic glaciers with complex geometry, including the George VI, Abbott, Dotson, Crosson, and Getz glaciers, resulted in (V0-nLV0-nE) values that are negligible (smaller than 1σ) in all five glaciers. We further investigated Area 2 of PIG (Fig. 3a), a section with a high velocity of 3278 m a−1 and the highest curvature along the main trunk (7 %, namely, a 1305 m cross flow deviation from the 19 720 m long straight line). Based on the feature-tracking result using eight annual Landsat 8 images, this curvature-induced E and L velocity difference is calculated as 206 m a−1 over a 7-year span, among which 195 m a−1 (95 %) was corrected (correction for Area 2 in Table 2). The results in these glaciers suggest that the velocity overestimation in glaciers with complex geometry can be mostly corrected and should not affect the applicability of the proposed method.

4.3 OE correction in fast-flowing glaciers

The Totten Glacier and Pine Island Glacier are among the most dynamic glaciers in Antarctica. The accelerated ice mass loss in these two fast-flowing glaciers has been realized by ice shelf basal melting and front calving (Li et al., 2015; Joughin et al., 2020, 2021). In Experiment 1 and Experiment 3 we show that such temporal accelerations do not alter the necessary condition in Premise II significantly, so that aVL(n)-aVLkσΔtn holds with k=1 and 3 on average for all experiment areas in TG (Table 1) and PIG (Table A4), respectively. The spatial-acceleration-induced OEs can be corrected. On the other hand, the impact of temporal accelerations is represented in the sufficient condition in Premise I, (V0-nL-V0-nE) kσ. For example, Area 4 of TG is located in a dynamic area of a velocity of  1400 m a−1, about 30 km from the ice shelf front (Fig. 4). Its E velocity (red curve in panel 4c in Fig. 5) exceeded the L velocity (blue curve) in 2016 and showed a significant acceleration thereafter. This resulted in the largest actual overestimation OEActual of  69 m a−1 (Table 1), which may be attributed to large shelf front calving occurring during the period between 6 April and 20 September 2015, with an area loss of 136 km2 (2 %; Fig. 4b). The L velocities V2013–iL (blue curve in panel 4c of Fig. 5) are derived from the 1-year-span map V2013–2014 before 2015 and thus do not reflect this temporal acceleration; furthermore, in L velocities V2013–iL(7) (black curve) derived from the 7-year-span map V2013–2020, the effect of the calving-induced acceleration is averaged out over the 7-year span. Thus, the back and blue lines are approximately parallel with an average acceleration difference of only −0.5 m a−2 (Table 1). Consequently, it appears that the calving event did not significantly reduce the buttressing potential of the ice shelf because the ice shelf front retreat occurred outside of the passive shelf ice (PSI) boundary (light blue line in Fig. 4b; Fürst et al., 2016); the other four areas are located farther away, inland, and away from the PSI boundary and thus were not influenced as much as in Area 4. Therefore, the calving-induced temporal acceleration only occurred in Area 4 (k= 1.5 in Premise I) but not in the four other areas (k≤0.5). Hence, the calculated correction of −19 m a−1 in Area 4 is at a similar level to that in the four other areas and has only adjusted the OE caused by the spatial acceleration along the trajectory. The velocity increase due to the temporal acceleration by the calving event remains in the adjusted map as a large residual of 50 m a−1, which can serve as a signature to study the relationship between calving activities and ice flow dynamics.

With k= 1, 2, 2, 10, and 5 in Premise I for the five selected areas from the inland interior to the grounding zone in PIG (Table A4), the temporal accelerations caused by climate warming impact are indicated by the high values of k (10 and 5) near the grounding line at the glacier outlet. Overall, the proposed method was able to correct  40 % of the total OEs, leaving another  60 % in residuals as the signature of the velocity changes induced by the temporal accelerations.

In addition, based on an annual E velocity map of 2013 in PIG from ITS_LIVE (Fig. A4a), L velocities along the grounding line (GL; Gardner et al., 2018) with time spans of 1 to 15 years (Fig. A4b) and the associated ice discharge were computed. The actual flux gates were set with nodes separated every 240 m, which were located up to 13 km upstream of the GL to reduce the uncertainty of ice thickness data (Gardner et al., 2019) from BedMachine (Morlighem et al., 2020). The results show that the velocity OEs of the 15-year span can reach up to  1300 m a−1 in the GL region. Such OEs in velocity can further cause an overestimation in ice discharge across flux gates upstream the GL, which increases rapidly by  6.3 Gt a−1 within the first 4-year span (Fig. A4c); thereafter, the flux overestimation slows down until a maximum of 11.5 Gt a−1 is reached at the 15-year span. This suggests that in fast-flowing glaciers like PIG, OEs in velocity maps with a time span of greater than 0.5–2 years should be corrected. The fact that the OE effect in PIG appears to level off after a couple of years is mainly because the velocity in the majority of the ice shelf is approximately leveled at  4000 km a−1 and the acceleration is thus very small. This makes the average L velocities along the GL increase in a much-reduced pace where an integral of E velocities is performed in the leveled velocity area a few years after crossing the GL (Fig. A4d). Consequently, the flux OE indicated in Fig. A4c showed a similar trend. Extension of this exercise to the entire continent involves a complex situation, including deceleration caused by ice rises, mapping errors near ice shelf shear margins, GL data points that would calve into the ocean within n-year span, and other potential factors. It should be noted that given a velocity range in flux computation, the ice flow angle between a flux gate and ice flow controls the flux magnitude. In principle, at the same gate this angle changes with velocities of different time spans and flow line patters, resulting in an overestimation or underestimation in flux for individual glaciers. A systematic and in-depth investigation should be carried out to handle the overestimation issue in GL flux of all of Antarctica.

4.4 Comparison with the “midpoint” method

The “midpoint” method presented in Berthier et al. (2003) compensates for overestimations by assigning the overestimated velocities to middle points of the trajectories. We use the velocity measurements in Experiment 1 to compare the performances of these two OE correction methods. Since Area 4 was affected by a calving event during the time span, here we use the other four areas (areas 1, 2, 3, and 5, Fig. 4a). We estimated the E velocities of 7 years (2013–2020) V2013–2020E and assigned them to the midpoints of the trajectories in the four areas (Table A5). They were then compared to the 1-year-span velocity V2013–2014E at the midpoints to calculate the bias εM. Similarly, the same overestimated 7-year-span E velocities V2013–2020E were corrected using the OE correction method of this paper and assigned to the start points of the trajectories as Vcorrected (Table A5). They were then compared to the 1-year-span V2013–2014E also, but at the start points to calculate another set of bias εS. As shown in Table A5, the proposed OE correction method achieved a higher overall accuracy of 4 ± 10 m a−1 compared to 12 ± 14 m a−1 of the midpoint method.

In summary, our validation and application experiment results demonstrated that the proposed OE correction method can be applied to different types of glaciers in Antarctica, regardless of their ice flow dynamics. While the OEs are effectively corrected, temporal velocity increases caused by global-warming-induced activities can be preserved for long-term change studies. The implication is that, when using newer velocity maps of short spans along with historical maps of long spans produced in previous studies over the past few decades, the overestimation of historical velocities could have caused an underestimation in the long-term acceleration magnitude. On the other hand, new efforts in historical velocity mapping at an ice-sheet-wide or a large regional scale should be made with a full consideration of OE corrections. The applicability of the method to glaciers in Greenland should be further investigated in future research.

5 Conclusions

Velocity overestimation exists in Antarctic ice flow velocity maps produced from optical satellite images of long time spans. Such overestimations are inevitable for historical velocity maps due to the poor availability of earlier satellite images in Antarctica, especially before 1990. The results in this study show that the overestimations can reach up to  69 m a−1 in the Totten Glacier, East Antarctica, over a 7-year span and  931 m a−1 in the Pine Island Glacier, West Antarctica, over a 10-year span. The overestimated historical velocity maps should be adjusted before they are combined with recent velocity maps to build a long-term record for monitoring and forecasting the Antarctic ice flow dynamics and impact of global climate changes on the ice sheet. We used a set of ground truth velocity maps in the Totten Glacier produced from Landsat 8 images of 2013 to 2020 to validate the proposed innovative method for velocity overestimation correction. Based on the validated method, we successfully corrected the overestimations of a velocity map of 1972–1989 with time spans from 1 to 17 years in the David Glacier region, East Antarctica. Another experiment at the Pine Island Glacier (PIG) was carried out to demonstrate that the method can effectively adjust overestimated velocities as large as 195 m a−1, while preserving the velocity change signature for long-term ice flow dynamics analysis. In summary we draw the following conclusions.

  1. The proposed Lagrangian velocity-based method is effective and easy to implement because the overestimation corrections are calculated by using the long-term velocity map itself only, without field observations or additional image data.

  2. The premises of the correction method are validated by using a set of ground truth velocity maps developed from high-quality Landsat 8 images from 2013 to 2020 to show the rigidity of the method. The velocity overestimations of up to a 7-year span are proven to be corrected effectively to within the uncertainty (1σ) of the 1-year map.

  3. The validated correction method is then successfully applied to correct overestimations in a historical velocity map of 1972–1989 with time spans from 1 to 17 years.

  4. It is proven that ice velocity change information of temporal acceleration events, e.g., caused by shelf front calving and basal melting in the Totten Glacier and Pine Island Glacier, is preserved after the correction and can be used for long-term ice flow dynamics analysis.

The magnitude of the OEs contained in the long-span historical velocity maps is significant. When overestimated historical maps are used alongside recent maps for assessment of the long-term global climate change impact on the Antarctic ice sheet and for forecast modeling, the underestimated long-term changes and compromised forecasting may result. Thus, the OEs in the long-span historical maps must be seriously examined and corrected. It is recommended that OE values be computed for long-term historical velocity maps and corrections be made when OEs are more than the velocity uncertainty (1σ). This velocity overestimation correction method can be applied to adjust velocity fields for the production of regional and ice-sheet-wide historical velocity maps from long-term satellite images before the 1990s.

Appendix A

Figure A1Velocity map of the Totten Glacier from ITS_LIVE (left) used to explain the concept of velocity overestimation caused by acceleration (right): “Over 16 years, that parcel of ice travels about 13 km downstream (red path). It begins at a velocity of about 720 m/yr, and in the first 8 months it travels at an average rate very close to 720 m/yr. But then the ice picks up speed as it moves downstream, so in the first 10 years it does not travel just 7200 m – it actually travels about 7900 m, or an average speed of 790 m/yr” (Chad A. Greene, personal communication, 2020).

Figure A2Overestimation correction in case of temporal accelerations.


Figure A3Velocities in five areas (numbered yellow rectangles in Fig. 3) of PIG are used to demonstrate application of the OE correction method.​​​​​​​ Panels (1a)(5a) show the reconstructed 3-month velocity maps V2013–2014 in the areas; matched points (red triangles) are used to map E velocity V2013–iE (i= 2014 (3 months), 2015, … , n) (red lines in panels 1c5c); points along the flow line (blue dots) are tracked from the 3-month maps and used to calculate L velocity V2013–iL (blue lines in panels 1c5c). Similarly, panels (1b)(5b) illustrate the reconstructed n-year velocity maps V2013–n in the areas with the matched points (red triangles) for E velocity V2013–iEn; points along the flow line (black crosses) are tracked from the n-year-span velocity map and used to calculate L velocity V2013–iLn (black lines in panels 1c5c). In each area (panels 1c5c) the difference between red and blue lines increases with time; the black line is above the blue line due to the spatial-acceleration-induced OE over n-year span, but they are relatively parallel, and thus a correction can be estimated. Background in panels (1a)(5a) and (1b)(5b) is Landsat 8 image of 12 February 2014.

Figure A4(a) One-year-span ice velocity map of PIG (2013) from ITS_LIVE (Gardner et al., 2019), (b) L velocities of 1- to 15-year span V0-iL (i= 1, 2, … 15) along GL calculated from the 2003 velocity map, (c) ice discharge across flux gates upstream the grounding line estimated with time spans from 1 to 15 years, and (d) L velocities V0-iL (averaged over the GL portion across main trunk) vs. time span.

Table A1Information for the Landsat 8 images used for velocity mapping in Experiment 1.

Download Print Version | Download XLSX

Table A2Acceleration estimates using least squares regression a and simple average a applied to L velocities of 7 years VL(n) and 1 year VL.

Download Print Version | Download XLSX

Table A3Information for the Landsat images used for velocity mapping in Experiment 2.​​​​​​​

Download XLSX

Table A4Application of the OE correction method in five selected areas in PIG. Premise I contains n-year tracked E velocities and L velocities computed from the 3-month velocity map, as well as their differences. Premise II illustrates averaged L accelerations computed from the n-year and 3-month velocity maps, respectively, as well as their differences.​​​​​​​

* With an average span of 6 years for all areas and σ=20 m a−1, the mean of aVL(n)-aVL in Premise II is 8 m a−2 (3σΔtn).

Download Print Version | Download XLSX

Table A5Comparison of the proposed OE correction method with the midpoint method in Berthier et al. (2003).

Download Print Version | Download XLSX

Code availability

Code is available from the corresponding author upon request.

Data availability

ITS_LIVE ice flow velocity maps used in this study are available at (last access: 26 February 2022) (Gardner et al., 2019).

Author contributions

RL led the study and developed the overestimation correction model. YC, MX, XY, ZL, and SL carried out velocity mapping. HC did the programming. RL, YC, and GQ were involved in data analysis and presentation.

Competing interests

The contact author has declared that neither they nor their co-authors have any competing interests.


Publisher’s note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.


We thank the editor and two reviewers for their constructive comments and suggestions. We also thank the United States Geological Survey (USGS) for the Landsat images.

Financial support

This research has been supported by the National Natural Science Foundation of China (grant no. 41730102), the National Key Research and Development Program of China (grant no. 2017YFA0603100), the Chinese Arctic and Antarctic Administration (grant no. CXPT2020017), and the Fundamental Research Funds for the Central Universities.

Review statement

This paper was edited by Etienne Berthier and reviewed by Chad Greene and one anonymous referee.


Altena, B. and Kääb, A.: Weekly Glacier Flow Estimation from Dense Satellite Time Series Using Adapted Optical Flow Technology, Front. Earth. Sci., 5, 53,, 2017. 

Bamber, J. L., Vaughan, D. G., and Joughin, I.: Widespread complex flow in the interior of the Antarctic ice sheet, Science, 287, 1248–1250,, 2000. 

Berthier, E., Raup, B., and Scambos, T.: New velocity map and mass-balance estimate of Mertz Glacier, East Antarctica, derived from Landsat sequential imagery, J. Glaciol., 49, 503–511,, 2003. 

Bindschadler, R., Vornberger, P., Blankenship, D., Scambos, T., and Jacobel, R.: Surface velocity and mass balance of Ice Streams D and E, West Antarctica, J. Glaciol., 42, 461–475,, 1996. 

Bindschadler, R., Vornberger, P., Fleming, A., Fox, A., Mullins, J., Binnie, D., Paulsen, S. J., Granneman, B., and Gorodetzky, D.: The Landsat Image Mosaic of Antarctica, Remote Sens. Environ., 112, 4214–4226,, 2008. 

Bindschadler, R. A. and Scambos, T. A.: Satellite-Image-Derived Velocity Field of an Antarctic Ice Stream, Science, 252, 242–246,, 1991. 

Chander, G., Markham, B. L., and Helder, D. L.: Summary of current radiometric calibration coefficients for Landsat MSS, TM, ETM+, and EO-1 ALI sensors, Remote Sens. Environ., 113, 893–903,, 2009. 

Cheng, Y., Li, X., Qiao, G., Ye, W., Huang, Y., Li, Y., Wang, K., Tian, Y., Tong, X., and Li, R.: Ice flow velocity mapping of East Antarctica from 1963 to 1989, International Society for Photogrammetry and Remote Sensing (ISPRS) Geospatial Week 2019, Enschede, the Netherlands, 10–14 June 2019, XLII-2/W13,, 2019. 

Chenillat, F., Blanke, B., Grima, N., Franks, P. J. S., Capet, X., and Riviere, P.: Quantifying tracer dynamics in moving fluids: a combined Eulerian-Lagrangian approach, Front. Earth Sci., 3, 43,, 2015. 

Chu, P. C. and Fan, C.: Accuracy Progressive Calculation of Lagrangian Trajectories from a Gridded Velocity Field, J. Atmos. Ocean. Tech., 31, 1615–1627,, 2014. 

Church, J. A., Monselesan, D., Gregory, J. M., and Marzeion, B.: Evaluating the ability of process based models to project sea-level change, Environ. Res. Lett., 8, 014051,, 2013. 

Cram, T. A., Persing, J., Montgomery, M. T., and Braun, S. A.: A lagrangian trajectory view on transport and mixing processes between the eye, eyewall, and environment using a high-resolution simulation of Hurricane Bonnie (1998), J. Atmos. Sci., 64, 1835–1856,, 2007. 

Cuffey, K. M. and Paterson, W. S. B.: Physics of Glaciers, 4th edn., Academic Press, Burlington, USA, ISBN-13 978-0123694614, ISBN-10 0123694612, 2010. 

Danesi, S., Dubbini, M., Morelli, A., Vittuari, L., and Bannister, S.: Joint geophysical observations of ice stream dynamics, in: Geodetic and Geophysical Observations in Antarctica, Springer, Berlin, Heidelberg,, 2008. 

Debella-Gilo, M. and Kääb, A.: Sub-pixel precision image matching for measuring surface displacements on mass movements using normalized cross-correlation, Remote Sens. Environ., 115, 130–142,, 2011. 

Debella-Gilo, M. and Kääb, A.: Monitoring slow-moving landslides using spatially adaptive least squares image matching, in: Landslide Science and Practice, Springer, Berlin, Heidelberg, (last access: 26 February 2022), 2013. 

DeConto, R. M. and Pollard, D.: Contribution of Antarctica to past and future sea-level rise, Nature, 531, 591–597,, 2016. 

Euler, C., Riemer, M., Kremer, T., and Schömer, E.: Lagrangian Description of Air Masses Associated with Latent Heat Release in Tropical Storm Karl (2016) during Extratropical Transition, Mon. Weather Rev., 147, 2657–2676,, 2019. 

Feng, T., Mi, H., Scaioni, M., Qiao, G., Lu, P., Wang, W., Tong, X., and Li, R.: Measurement of Surface Changes in a Scaled-Down Landslide Model Using High-Speed Stereo Image Sequences, Photogramm. Eng. Rem. S., 82, 547–557,, 2016. 

Frezzotti, M., Capra, A., and Vittuari, L.: Comparison between glacier ice velocities inferred from GPS and sequential satellite images, Ann. Glaciol., 27, 54–60,, 1998. 

Frezzotti, M., Tabacco, I. E., and Zirizzotti, A.: Ice discharge of eastern Dome C drainage area, Antarctica, determined from airborne radar survey and satellite image analysis, J. Glaciol., 46, 253–264,, 2000. 

Fürst, J. J., Durand, G., Gillet-Chaulet, F., Tavard, L., Rankl, M., Braun, M., and Gagliardini, O.: The safety band of Antarctic ice shelves, Nat. Clim. Change, 6, 479–482,, 2016. 

Gardner, A. S., Moholdt, G., Scambos, T., Fahnstock, M., Ligtenberg, S., van den Broeke, M., and Nilsson, J.: Increased West Antarctic and unchanged East Antarctic ice discharge over the last 7 years, The Cryosphere, 12, 521–547,, 2018. 

Gardner, A. S., Fahnestock, M. A., and Scambos, T. A.: ITS_LIVE Regional Glacier and Ice Sheet Surface Velocities, National Snow and Ice Data Center [data set], (last access: 26 February 2022​​​​​​​), 2019. 

Glenn, S. M., Miles, T. N., Serokal, G. N., Xu, Y., Forney, R. K., Yu, F., Roarty, H., Schofield, O., and Kohut, J.: Stratified coastal ocean interactions with tropical cyclones, Nat. Commun., 7, 10887(2016)​​​​​​​,, 2016. 

Greene, C. A., Blankenship, D. D., Gwyther, D. E., Silvano, A., and van Wijk, E.: Wind causes Totten Ice Shelf melt and acceleration, Science Advances, 3, e1701681,, 2017. 

Greene, C. A., Young, D. A., Gwyther, D. E., Galton-Fenzi, B. K., and Blankenship, D. D.: Seasonal dynamics of Totten Ice Shelf controlled by sea ice buttressing, The Cryosphere, 12, 2869–2882,, 2018. 

Greene, C. A., Gardner, A. S., and Andrews, L. C.: Detecting seasonal ice dynamics in satellite images, The Cryosphere, 14, 4365–4378,, 2020. 

Halliday, D., Resnick, R., and Walker, J.: Fundamentals of physics, 10th edn., John Wiley & Sons, ISBN 111823071X, 9781118230718, 2013. 

Heid, T. and Kääb, A.: Evaluation of existing image matching methods for deriving glacier surface displacements globally from optical satellite imagery, Remote Sens. Environ., 118, 339–355,, 2012. 

IPCC: The Ocean and Cryosphere in a Changing Climate: A Special Report of the Intergovernmental Panel on Climate Change, Working Group II Technical Support Unit, 756 pp., (last access: 26 February 2022​​​​​​​), 2019. 

Joughin, I., Shean, D. E., Smith, B. E., and Floricioiu, D.: A decade of variability on Jakobshavn Isbræ: ocean temperatures pace speed through influence on mélange rigidity, The Cryosphere, 14, 211–227,, 2020. 

Joughin, I., Shapero, D., Smith, B., Dutrieux, P., and Barham, M.: Ice-shelf retreat drives recent Pine Island Glacier speedup, Science Advances, 7, eabg3080,, 2021. 

Kim, K.: Satellite mapping and automated feature extraction: geographic information system-based change detection of the Antarctic coast, Doctoral dissertation, The Ohio State University, Columbus, Ohio, 171 pp., (last access: 26 February 2022​​​​​​​), 2004. 

Li, J. and Roy, D. P.: A Global Analysis of Sentinel-2A, Sentinel-2B and Landsat-8 Data Revisit Intervals and Implications for Terrestrial Monitoring, Remote Sens., 9, 902,, 2017. 

Li, R.: Potential of high-resolution satellite imagery for national mapping products, Photogramm. Eng. Rem. S., 64, 1165–1169, 1998. 

Li, R., Hwangbo, J., Chen, Y., and Di, K.: Rigorous Photogrammetric Processing of HiRISE Stereo Imagery for Mars Topographic Mapping, IEEE T. Geosci. Remote, 49, 2558–2572,, 2011. 

Li, R., Ye, W., Qiao, G., Tong, X., Liu, S., Kong, F., and Ma, X.: A New Analytical Method for Estimating Antarctic Ice Flow in the 1960s From Historical Optical Satellite Imagery, IEEE T. Geosci. Remote, 55, 2771–2785,, 2017. 

Li, X., Rignot, E., Morlighem, M., Mouginot, J., and Scheuchl, B.: Grounding line retreat of Totten Glacier, East Antarctica, 1996 to 2013, Geophys. Res. Lett., 42, 8049–8056,, 2015. 

Li, X., Rignot, E., Mouginot, J., and Scheuchl, B.: Ice flow dynamics and mass loss of Totten Glacier, East Antarctica, from 1989 to 2015, Geophys. Res. Lett., 43, 6366–6373,, 2016. 

Lillesand, T. R. W. K.: Remote sensing and image interpretation, 7th edn., John Wiley & Sons, ISBN-13 978-1118343289, ISBN-10 111834328X, 2015. 

McGlone, J. C.: Manual of Photogrammetry, sixth edn., ASPRS Publications, Maryland, MD, USA, ISBN-10 1570830991, ISBN-13 978-1570830990, 2013. 

Montgomery, D. C., Peck, E. A., and Vining, G. G.: Introduction to Linear Regression Analysis, 6th edn., John Wiley & Sons, ISBN-10 1119578728, ISBN-13 978-1119578727, 2021. 

Moon, J., Cho, Y., and Lee, H.: Flow Velocity Change of David Glacier, East Antarctica, from 2016 to 2020 Observed by Sentinel-1A SAR Offset Tracking Method, Korean Journal of Remote Sensing, 37, 1–11​​​​​​​,, 2021. 

Morlighem, M., Rignot, E., Binder, T., Blankenship, D., Drews, R., Eagles, G., Eisen, O., Ferraccioli, F., Forsberg, R., Fretwell, P., Goel, V., Greenbaum, J. S., Gudmundsson, H., Guo, J., Helm, V., Hofstede, C., Howat, I., Humbert, A., Jokat, W., Karlsson, N. B., Lee, W. S., Matsuoka, K., Millan, R., Mouginot, J., Paden, J., Pattyn, F., Roberts, J., Rosier, S., Ruppel, A., Seroussi, H., Smith, E. C., Steinhage, D., Sun, B., Broeke, M. R. V. D., Ommen, T. D. V., Wessem, M. V., and Young, D. A.: Deep glacial troughs and stabilizing ridges unveiled beneath the margins of the Antarctic ice sheet, Nat. Geosci., 13, 132–137,, 2020. 

Nakamura, K., Doi, K., and Shibuya, K.: Fluctuations in the flow velocity of the Antarctic Shirase Glacier over an 11-year period, Polar Sci., 4, 443–455,, 2010. 

Rignot, E., Mouginot, J., and Scheuchl, B.: Ice Flow of the Antarctic Ice Sheet, Science, 333, 1427–1430, 2011a. 

Rignot, E., Velicogna, I., van den Broeke, M. R., Monaghan, A., and Lenaerts, J.: Acceleration of the contribution of the Greenland and Antarctic ice sheets to sea level rise, Geophys. Res. Lett., 38, L05503,, 2011b. 

Rignot, E., Mouginot, J., and Scheuchl, B.: Antarctic grounding line mapping from differential satellite radar interferometry, Geophys. Res. Lett., 38, L10504,, 2011c. 

Rignot, E., Mouginot, J., Scheuchl, B., van den Broeke, M., van Wessem, M. J., and Morlighem, M.: Four decades of Antarctic Ice Sheet mass balance from 1979–2017, P. Natl. Acad. Sci. USA, 116, 1095–1103,, 2019. 

Ruffner, K. C.: CORONA: America’s First Satellite Program, (last access: 26 February 2022​​​​​​​), 1995. 

Sabins Jr., F. and James, F. A.: Remote Sensing: Principles, Interpretation, and Applications, 4th edn., Waveland Press, ISBN-10 1478637102, ISBN-13 978-1478637103, 2020. 

Scambos, T. A., Dutkiewicz, M. J., Wilson, J. C., and Bindschadler, R. A.: Application of image cross-correlation to the measurement of glacier velocity using satellite image data, Remote Sens. Environ., 42, 177–186,, 1992. 

Schenk, T.: Digital Photogrammetry, TerraScience, Laurelville, Ohio, ISBN-10 0967765315, 1999. 

Shen, Q., Wang, H., Shum, C. K., Jiang, L., Hsu, H. T., and Dong, J.: Recent high-resolution Antarctic ice velocity maps reveal increased mass loss in Wilkes Land, East Antarctica, Sci. Rep.-UK, 8, 4477,, 2018. 

Shepherd, A., Ivins, E. R., Geruo, A., Barletta, V. R., Bentley, M. J., Bettadpur, S., Briggs, K. H., Bromwich, D. H., Forsberg, R., Galin, N., Horwath, M., Jacobs, S., Joughin, I., King, M. A., Lenaerts, J. T. M., Li, J., Ligtenberg, S. R. M., Luckman, A., Luthcke, S. B., McMillan, M., Meister, R., Milne, G., Mouginot, J., Muir, A., Nicolas, J. P., Paden, J., Payne, A. J., Pritchard, H., Rignot, E., Rott, H., Sorensen, L. S., Scambos, T. A., Scheuchl, B., Schrama, E. J. O., Smith, B., Sundal, A. V., van Angelen, J. H., van de Berg, W. J., van den Broeke, M. R., Vaughan, D. G., Velicogna, I., Wahr, J., Whitehouse, P. L., Wingham, D. J., Yi, D., Young, D., and Zwally, H. J.: A Reconciled Estimate of Ice-Sheet Mass Balance, Science, 338, 1183–1189,, 2012. 

Shepherd, A., Ivins, E., Rignot, E., Smith, B., van den Broeke, M., Velicogna, I., Whitehouse, P., Briggs, K., Joughin, I., Krinner, G., Nowicki, S., Payne, T., Scambos, T., Schlegel, N., Geruo, A., Agosta, C., Ahlstrom, A., Babonis, G., Barletta, V., Blazquez, A., Bonin, J., Csatho, B., Cullather, R., Felikson, D., Fettweis, X., Forsberg, R., Gallee, H., Gardner, A., Gilbert, L., Groh, A., Gunter, B., Hanna, E., Harig, C., Helm, V., Horvath, A., Horwath, M., Khan, S., Kjeldsen, K. K., Konrad, H., Langen, P., Lecavalier, B., Loomis, B., Luthcke, S., McMillan, M., Melini, D., Mernild, S., Mohajerani, Y., Moore, P., Mouginot, J., Moyano, G., Muir, A., Nagler, T., Nield, G., Nilsson, J., Noel, B., Otosaka, I., Pattle, M. E., Peltier, W. R., Pie, N., Rietbroek, R., Rott, H., Sandberg-Sorensen, L., Sasgen, I., Save, H., Scheuchl, B., Schrama, E., Schroeder, L., Seo, K., Simonsen, S., Slater, T., Spada, G., Sutterley, T., Talpe, M., Tarasov, L., van de Berg, W. J., van der Wal, W., van Wessem, M., Vishwakarma, B. D., Wiese, D., and Wouters, B.: Mass balance of the Antarctic Ice Sheet from 1992 to 2017, Nature, 558, 219–222,, 2018. 

Shimizu, M. and Okutomi, M.: Sub-pixel estimation error cancellation on area-based matching, Int. J. Comput. Vision, 63, 207–224,, 2005. 

Storey, J., Choate, M., and Lee, K.: Landsat 8 Operational Land Imager On-Orbit Geometric Calibration and Performance, Remote Sens.-Basel, 6, 11127–11152,, 2014. 

van Sebille, E., Griffies, S. M., Abernathey, R., Adams, T. P., Berloff, P., Biastoch, A., Blanke, B., Chassignet, E. P., Cheng, Y., Cotter, C. J., Deleersnijder, E., Doos, K., Drake, H. F., Drijfhout, S., Gary, S. F., Heemink, A. W., Kjellsson, J., Koszalka, I. M., Lange, M., Lique, C., MacGilchrist, G. A., Marsh, R., Adame, C. G. M., McAdam, R., Nencioli, F., Paris, C. B., Piggott, M. D., Polton, J. A., Ruehs, S., Shah, S. H. A. M., Thomas, M. D., Wang, J., Wolfram, P. J., Zanna, L., and Zika, J. D.: Lagrangian ocean analysis: Fundamentals and practices, Ocean Model., 121, 49–75,, 2018. 

Wang, S., Liu, H., Yu, B., Zhou, G., and Cheng, X.: Revealing the early ice flow patterns with historical Declassified Intelligence Satellite Photographs back to 1960s, Geophys. Res. Lett., 43, 5758–5767,, 2016. 

Wulder, M. A., Loveland, T. R., Roy, D. P., Crawford, C. J., Masek, J. G., Woodcock, C. E., Allen, R. G., Anderson, M. C., Belward, A. S., Cohen, W. B., Dwyer, J., Erb, A., Gao, F., Griffiths, P., Helder, D., Hermosillo, T., Hipple, J. D., Hostert, P., Hughes, M. J., Huntington, J., Johnson, D. M., Kennedy, R., Kilic, A., Li, Z., Lymburner, L., McCorkel, J., Pahlevan, N., Scambos, T. A., Schaaf, C., Schott, J. R., Sheng, Y., Storey, J., Vermote, E., Vogelmann, J., White, J. C., Wynne, R. H., and Zhu, Z.: Current status of Landsat program, science, and applications, Remote Sens. Environ., 225, 127–147,, 2019.  

Ye, W., Qiao, G., Kong, F., Ma, X., Tong, X., and Li, R.: Improved geometric modeling of 1960s KH-5 ARGON satellite images for regional antarctica applications, Photogramm. Eng. Rem. S., 83, 477–491, 2017. 

Zhou, C., Zhou, Y., Deng, F., Ai, S., Wang, Z., and E, D.: Seasonal and interannual ice velocity changes of Polar Record Glacier, East Antarctica, Ann. Glaciol., 55, 45–51,, 2014. 

Short summary
Historical velocity maps of the Antarctic ice sheet are valuable for long-term ice flow dynamics analysis. We developed an innovative method for correcting overestimations existing in historical velocity maps. The method is validated rigorously using high-quality Landsat 8 images and then successfully applied to historical velocity maps. The historical change signatures are preserved and can be used for assessing the impact of long-term global climate changes on the ice sheet.