Detailed ice loss pattern in the northern Antarctic Peninsula: widespread decline driven by ice front retreats

. The northern Antarctic Peninsula (nAP, < 66 ◦ S) is one of the most rapidly changing glaciated regions on earth, yet the spatial patterns of its ice mass loss at the glacier basin scale have to date been poorly documented. We use satellite laser altimetry and satellite stereo-image topography spanning 2001–2010, but primarily 2003–2008, to map ice elevation change and infer mass changes for 33 glacier basins covering the mainland and most large islands in the nAP. Rates of ice volume and ice mass change are 27 . 7 ± 8 . 6 km 3 a − 1 and 24 . 9 ± 7 . 8 Gt a − 1 , equal to − 0 . 73 m a


Introduction
The northern Antarctic Peninsula (nAP) is one of two areas of the Antarctic ice sheet showing major mass loss, the other being the Amundsen Sea coast of West Antarctica's ice sheet.
Previous studies have shown large negative mass imbalances and significant elevation losses for the nAP (Ivins et al., 2011;Shepherd et al., 2012;Luthcke et al., 2013;Sasgen et al., 2013;McMillan et al., 2014). However, in general these studies have not resolved the spatial distribution of mass imbalance in detail. A mapping of these patterns, over the entire region but at the glacier basin scale, can provide more insight into processes responsible for the ice loss. Studies based on gravitational change detection using the Gravity Recovery and Climate Experiment satellite system (GRACE) have an inherent spatial resolution of roughly 250 km scale (Ivins et al., 2011;Luthcke et al., 2013;Sasgen et al., 2013), far larger than the scale of the nAP individual glacier basins and islands. Past altimetry-based studies (Pritchard et al., 2009;Flament and Rémy, 2012;Shepherd et al., 2012;McMillan et al., 2014) suffer from either sparse coverage or slope correction issues, or both, due to the steep terrain in the nAP. In the published assessments based on laser altimetry (Shepherd et al., 2012), broad assumptions and large extrapolations are required to interpolate the data across the dissected and rugged peninsula region. Mass budget methods (Rignot et al., 2004(Rignot et al., , 2008Rott et al., 2011;Shepherd et al., 2012), which aim to difference outflowing ice flux and surface mass balance (SMB) for each glacier basin, have to date shown results that are difficult to reconcile with other studies of the same glaciers (Shuman et al., 2011;Berthier et al., 2012). This is primarily due to spatially coarse SMB estimates from Published by Copernicus Publications on behalf of the European Geosciences Union. models or field measurements, difficulties in estimating the cross-sectional area of the glaciers, and differences in the span of time used to estimate ice flux changes (Berthier et al., 2012).
The goal of this study is to determine the spatial pattern of ice elevation changes in the nAP, improve estimates of mass balance for the region, and study the relationship of mass balance with ice shelf collapse and ice front retreats in the area. In light of known climate-related changes in the region, such as increasing surface air temperatures and surface melting, regional sea ice decline, and increasing accumulation (e.g., Mulvaney et al., 2012;Zagorodnov et al., 2012;Stammerjohn et al., 2008;Lenaerts et al., 2012;Barrand et al., 2013), our study reveals a pattern of ice mass loss in space and (we infer) in time that may be similar to the characteristics of mass loss in other areas of Antarctica in the coming century.

Methods
Our study combines satellite stereo-image digital elevation model differencing (dDEM) with repeat-track laser altimetry from the Ice, Cloud and land Elevation Satellite (ICE-Sat; Schutz et al., 2005), with the objective of providing an assessment of surface elevation change resolved at the scale of the major glacier catchments. We used stereo-image data from the Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER; Fujisada et al., 2005) and Satellite Pour l'Observation de la Terre 5 (SPOT5; Korona et al., 2009). Eight satellite stereo-image data sets from the ASTER sensor, and six from the SPOT-5 Haute Résolution Stéréoscopique (HRS) sensor (Table S1 and Fig. S1 in the Supplement) were processed using previously published methods (Shuman et al., 2011;Berthier et al., 2012;Gardelle et al., 2013).
For the ICESat repeat-track data (Release 633), we used 26 ground tracks from the 91-day-repeat orbit crossing the nAP and major ice-covered islands for the high-energy laser campaigns (ICESat Laser 2A through Laser 3J, September 2003-March 2008Shuman et al., 2006;Zwally et al., 2012). Cross-track elevation adjustment and along-track filtering are used to improve measurement quality, based on surface slopes (not elevations) derived from a recent Antarctic Peninsula DEM (Cook et al., 2012). We first eliminated ICESat profile tracks more than 300 m from the reference track position, and sections where the absolute slope from the gridded DEM was > ± 10 % slope (or ±5.7 • ) for the reference track or measurement track location, or areas where the absolute difference between along-track slopes of the measurement track and reference track exceeded 5 % (or ±2.9 • ). We further required the ICESat elevation data to be within 50 m (vertically) of the corresponding interpolated DEM elevation. All elevations are referenced to the Earth Geopotential Model 1996 (EGM96) geoid datum. To migrate the measurement track data to the reference track and compare ele-vations, we identified reference track stations every 43.75 m along the reference track (one-fourth the distance between ICESat altimetry shot locations along track). We then applied an elevation correction based on the difference between the interpolated gridded DEM elevation at the nearest reference track station and the ICESat track data point. ICESat campaign data were compared by differencing their migrated elevations, divided by the time in years between dates of track acquisition. To reduce effects of possible seasonal variations in elevation, we compared only near-integer-year separated repeat profiles, e.g., data from campaigns 2A to 3A (∼ October, ∼ 1 year apart) or 3B to 3H (∼ March, ∼ 2 years apart).
To evaluate different processes in elevation and ice mass change, we treat regions above and below 1000 m above sea level (a.s.l.) separately for each of 33 drainage basins. This is the approximate elevation of the crest of an extensive escarpment in the nAP separating plateau areas from individual glacier cirques. Above 1000 m a.s.l., and for islands without sufficient dDEM coverage (Robertson Islands, Snow Hill Island, and Joinville, Dundee, and D'Urville islands; Fig. 1), the rate of elevation change (dH /dt) is determined from satellite laser altimetry alone. In smooth high-elevation areas, correlation of satellite stereo-images often fails due to a lack of high-contrast surface features of sufficient horizontal scale (tens of meters). Below 1000 m a.s.l., a hypsometric interpolation method was applied to individual glacier basins to extend dDEMs and ICESat dH /dt measurements to areas not directly measured. ICESat dH /dt was weighted 10-fold relative to dDEM dH /dt to prevent small dDEM data areas from dominating the weighted dH /dt mapping, and to better utilize the higher accuracy of individual ICESat-based measurements. We used the relationship where N is the number of measurements within an elevation band (i.e., the number of 50 m grid cells for the dDEMs, or reference track site locations at 43.75 m spacing for ICE-Sat) and e is an inverse weighting of the measurement methods. For dDEMs we used a weight of 1, and for ICESat we used 0.1. This allowed the fewer but more accurate ICESatbased measurements to contribute to the final result in basins with extensive dDEM coverage. In several areas, ICESat data were available in regions not well covered by dDEM results (see Fig. S2). We also estimate the above-flotation mass loss of grounded-ice areas that retreated at least 2 km 2 during the study interval (2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010), as identified by image mapping (Cook et al., 2005;Cook and Vaughan, 2010). To estimate the volume and mass loss represented by these areas, we mapped the area of retreat during our study period (2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010) and half the mean elevation loss rate observed just above the area of grounded-ice retreat. This represents an assumption that the vertical elevation change rate of the retreated ice was identical to the region just upstream of the loss area, and that  Tables S2 and  S3. Major drainage basins evaluated by the study are outlined in white, sub-basins are indicated in blue. Base image is the MODIS Mosaic of Antarctica . Inset, location of the study area shown in Fig. 2. the time of ice front retreat (e.g., when the ice calved and drifted away) was midway through the study period. Errors for our assessment of dH /dt (Tables 1 and S2), are based on past analysis of the dDEM method (Shuman et al., 2011;Berthier et al., 2012), on inter-comparisons of the two methods at sites having both dDEM and ICESat measurements, and on crossover analysis of ICESat crosstrack-corrected data (Table S4). Past analysis for this region suggests that dDEM methods using mixed ASTER and SPOT5 imagery can have a ±5 m uncertainty for individual glacier basins, i.e., ∼ 1 m a −1 given a 5 yr time separation between DEMs. However, examining our ICESat and dDEM dH /dt at sites with both measurements (6158 sites) shows that the methods differ by just ∼ 0.3 m a −1 overall, a difference that ranges between 0.07 and 0.75 m a −1 over various sub-sets of our measurements (Table S4). This is in agreement with a previous study that showed reduced errors when dDEM results are averaged over basin-scale areas (Berthier et al., 2012). Seven crossover sites with slope-corrected ICE-Sat dH /dt measurements show good agreement with the dDEM measurements at the same locations (mean offset of +0.05 m a −1 ).
Errors in the ICESat cross-track correction for dH /dt are more dependent on slope errors in the Cook et al. (2012) DEM and not its absolute elevation accuracy. Assuming our selection criteria eliminated regions of significant error in the DEM, we estimate that across-track or along-track slopes in the Cook et al. (2012) DEM are accurate to within ±0.5 • over a length scale of 300 m, or ±8.7 m km −1 . A test of this was conducted by comparing the Cook et al. (2012) DEM slopes with a DEM acquired in 2009 by the NASA Land, Vegetation, and Ice Sensor (LVIS) airborne laser altimeter, covering about 20 % of the study regions (http://lvis.gsfc. nasa.gov/Crane.html). This showed that the mean difference in along-track slope in the overlap region was 0.06 ± 1.2 • when our criteria are applied to both data sets. For laser altimetry measurements alone, our inferred mean slope error of ±0.5 • implies a mean laser measurement pair cross-trackcorrection error of ±1.31 m (assuming a mean cross-track distance of 150 m). We assume this error is randomly distributed when averaging over a glacier basin. Thus, for the average of 20 measurement sites, the mean error is < 25 cm. Since laser measurement pairs may have 1 to 4 years separation in time, and many have multiple measurements at a single site (see Fig. S2), our overall mean error in elevation change rate is significantly less than this. Additionally, the majority of the basins we consider have many more than 20 measurement sites (Table S2).
Considering all sources of error, and variations in the time span of measurements for dDEM and ICESat measurements, data density variations for the basins, and the strong agreement between these independent altimetric methods, we adopt a mean error of ±0.15 m a −1 for regions of laser altimetry measurement alone (above 1000 m a.s.l.), and ±0.3 m a −1 for our dDEM plus altimetry measurements (below 1000 m a.s.l.) and the glacier basins, islands, and subbasins without laser altimetry. Errors for volume and mass change estimates scale with basin area.

Results
An overview of our results is shown in Fig. 2 and Table 1, and detailed basin-by-basin values are provided in Table S2. The results show that basins impacted by recent ice shelf loss and ice front retreat have very high rates of change, but also indicate that few areas -high or low, east or west -have positive dH /dt. Recent ice-shelf loss (ISL) basins (losses since 1986, and particularly in 1995 and 2002), all on the eastern side of the nAP, and four smaller glaciers with recent grounded-ice front loss (IFL; losses since 2000) on the western and northeastern side of the nAP, show a characteristic pattern of very high elevation loss rates just upstream of the ice front but far lower elevation loss rates at high elevation. Abbreviations for place names: nAP, northern Antarctic Peninsula; ISL, ice shelf loss; IFL, ice front loss. a Assuming mean density of 900 kg m −3 for all dV /dt measurements. Errors for these values are 0.9 times the sum of errors for dV /dt for each row; b Area determined from additional ASTER, SPOT, and Landsat images, spanning 2000-2002-2010 Rate of elevation loss measured for the first 50 m elevation band above area of grounded-ice retreat; d Volume loss assumes flotation was reached midway between 2001 and 2010 (period of observations); e Percent area covered by differential DEM satellite stereo-image data; f Number of repeat-track point measurements used. If < 20 ICESat dH /dt measurements are available, the regional mean ICESat dH /dt for areas > 1000 m (−0.31 m a −1 ) or, for sub-basins, the main basin mean is used; g Hypsometric weighting for areas below 1000 m elevation; h Errors on dV /dt can be determined by ±0. change for areas below 1000 m a.s.l. at 12 eastern-side ISL glacier basins (or sub-basins) is −2.6 m a −1 (range, +0.4 to −5.8 m a −1 ) and −2.2 m a −1 (−2.0 to −2.7 m a −1 ) for the four western-side and northeastern IFL sub-basins experiencing recent ice front retreat (> 2 km 2 since 2000). At elevations > 1000 m, elevation loss in the eastern ISL basins is small (mean, −0.10 m a −1 ; range +0.35 to −0.54 m a −1 ). Glacier systems on the western nAP coast and the western islands below 1000 m a.s.l., excluding the recent IFL regions, are changing at various rates (typically ∼ −0.15 m a −1 , range +0.7 to −1.6 m a −1 ). However, western-side basins are losing elevation at significant rates above 1000 m a.s.l. (mean, −0.59 m a −1 ; range −0.25 to −1.30 m a −1 ). In terms of mean water equivalent layer, the overall study area has changed by an average of −0.73 m a −1 w.e. during the study period; the western glaciers (basins 1-11) averaged 0.33 m a −1 w.e., and eastern glacier systems impacted by ice shelf loss (basins 19, 21-25, 26b, 27-30, and 31a) averaged 1.48 m a −1 w.e. elevation loss. We examine the rates of surface elevation change and cumulative ice volume change as they vary with altitude for three sub-regions of the study area in Fig. 3. The patterns of elevation change with altitude illustrate the differences between the western-side glacier and island regions and the eastern-side ISL areas, and also highlight the bi-modal hypsometry pattern characteristic of the nAP. Eastern-side ISL areas show dramatically decreasing elevation, and large volume changes at low elevations, but relatively small elevation losses in the upper-most catchment areas (Fig. 3c-d). Western-side glaciers show negative rates of elevation change at all elevations, and a steady cumulative volume decrease rate with altitude. The major glaciers of Scar Inlet ice shelf, the lone remaining large (> 50 km 2 ) ice shelf in the study area with significant tributary glaciers, show a unique pattern of ice loss at low elevation and some areas of thickening at altitude. We believe this is likely the pattern of elevation change present for the eastern nAP ISL glacier systems in the years immediately prior to shelf disintegration.

Discussion
The widespread elevation losses suggested here for both sides of the nAP at high elevations, and especially for the western side of the divide, have significant implications for the region's recent mass change history. Moreover, the detailed mapping on a basin-by-basin scale supports model and GPS studies of local bedrock uplift. A recent study (Neild et al., 2014), using an earlier (near-final) version of our presented data combined with continuous GPS uplift measurements at sites in the peninsula, modeled both the elastic response and long-term isostatic rebound in the region, showing that the nAP is underlain by very low viscosity mantle.
Previous observational studies have shown that the elevation decline pattern for ISL or IFL glaciers here and in other similar glaciated regions migrates upstream and diffuses on a scale of years to decades (Howat et al., 2007;Joughin et al., 2008;Shuman et al., 2011;Berthier et al., 2012), consistent with kinematic wave models of glacier response to ice front stress changes for tidewater glaciers (Pfeffer, 2007;Nick et al., 2009;Favier et al., 2014). In past work in this area (Shuman et al., 2011;Berthier et al., 2012), and with comparison to these results, we observe that eastern ISL glaciers are currently propagating kinematic waves upstream from their lower trunk areas, but this process has not yet had a  (Cook and Vaughan, 2010) are indicated in gray-blue, with years of major collapse events and the limit of extensive grounded-ice loss shown. Ice edge is from a 2009 MODIS mosaic ). significant impact on higher elevations for the eastern glacier basins. Western-coast nAP glacier front retreats, elevation losses, and accelerations have been documented (Cook and Vaughan, 2005;Pritchard and Vaughan, 2007;Kunz et al., 2012;Christ et al., 2014), with a major pulse of retreat beginning in the 1970s. Moreover, our work here shows that ongoing ice front losses within the study period behave much like smaller versions of the eastern-side glaciers impacted by ice shelf and glacier front retreat (Table 1; Fig. 2). The earlier losses inferred for the western-side fjord glaciers (e.g., Christ et al., 2014) appear to have now propagated throughout the entirety of the western basins, leading to significant and widespread surface lowering in the western upper catchment areas (> 1000 m a.s.l.) at greater rates than for the eastern side on average (Tables 1 and S2) 1960-1969to 2000, respectively (Potocki et al., 2011Goodwin, 2013). Models of precipitation input for the region (Saha et al., 2014;Dee et al., 2011;Lenaerts et al., 2012) also show a strong overall increase for the most recent decades, but some indicate a slight decline in the last decade, covering our dH /dt measurement period (Saha et al., 2010;Lenaerts et al., 2012;Shepherd et al., 2012). The large increase and later reduction in accumulation are associated with multidecadal warming (Barrand et al., 2013) and with reductions in sea ice extent northwest of the nAP (Stammerjohn et al., 2012), recently moderated by a slight cooling trend (Fogt and Scambos, 2013;Zagorodnov et al., 2012).
Multi-decadal accumulation, temperature, and snowmelt trends cause changes in the compaction rate of snow and firn, and can potentially impact measurements of surface elevation change (Ligtenberg et al., 2011). Using a model climate time series (based on reanalysis of weather data) spanning the period of our measurements (Regional Atmospheric Climate Model, RACMO-2.1/ANT; Lenaerts et al., 2012), a dH /dt for the firn column at 27 km spatial scale is obtained similar to that used in previous related analyses Gardner et al., 2013). The modeled inter-annual variability in accumulation, temperature, and snowmelt and their effect on firn compaction result in dH /dt corrections between −0.19 and +0.12 m a −1 on the grounded ice of the nAP, with generally positive (thickening) corrections on the western side and negative to the east. The small effect on the firn layer, and the high variability of accumulation both inter-annually and among the basin areas (Fig. 4a) make the correction relatively insignificant. We therefore report dH /dt as observed from the satellite data. From these observations,  Table 1), eastern basins with major ice shelf loss in the period 1986-2009 (c and d; basins 19, 21-25, and 27-30 in Table 1), and basins draining to the Scar Inlet ice shelf area (e and f; basins 31b, 32, and 33 in Table 1). Height is binned in 50 m intervals. Note that rates of elevation change trends at the highest elevations (> 2000 m a.s.l., right side of left column of panels) are based on less data and are not reliable.
we report mass change in Tables 1 and S2 as where (dH /dt) hyps is the elevation-band-weighted mean measured dH /dt, A is area of the glacier basin or island, and ρ is our assumed mean density of ice and firn lost by dynamics (900 kg m −3 ). We eliminated the nunatak areas from each of the basins, based on the Antarctic Digital Database mapping of rock outcroppings in the region similar to previous studies (e.g., Gardner et al., 2013). Our estimate of mass balance for the combined nAP region is −24.9 ± 7.8 Gt a −1 , with the great majority of the mass loss occurring at elevations below 1000 m a.s.l. (−21.9 ± 6.3 Gt a −1 , or 88 %; Table 1). Regionally, the eastern nAP basins dominate the mass loss at −17.7±3.7 Gt a −1 , or 72 % of the loss, and of this, −15.2 ± 3.2 Gt a −1 (60 %) is from 12 glacier basins flowing into embayments formerly occupied by the Prince Gustav, Larsen Inlet, Larsen A, and Larsen B ice shelves. For the 11 western nAP glacier basins and islands, the mass loss rate is similar at low and high elevations (−2.3 ± 0.7 Gt a −1 > 1000 m a.s.l., and −2.2 ± 1.0 Gt a −1 below). Overall, the nAP region accounts for ∼ 29 % of Antarc-tica mass imbalance during the study period (Shepherd et al., 2012).
We also examined the mass balance ratio of the basins and regional areas, based on mass input, primarily snow accumulation ; Tables 2 and S3). Surface mass balance (SMB) in the region has a very large gradient from west to east, with values of 1500 to 3000 km m −2 a −1 in RACMO2.1 grid cells for the western areas and high elevations dropping to ∼ 500 to 1500 kg m −2 a −1 in the lowelevation areas of the eastern nAP coast. A ratio of the mass balance divided by the mass accumulation input indicates the degree of imbalance in the glacier systems, and suggests the level of ice flux increase for glacier systems having recently accelerated due to ice front or ice shelf losses. We term this value the imbalance ratio. The imbalance ratio for the nAP as a whole is −0.46, implying that mass outflow is 45 % greater during the study period relative to a steady-state rate in the current climate. For the eastern nAP glaciers, the mean ratio is −0.81 and the major ISL glaciers in the Larsen A and Larsen B between −0.33 and −2.41. The upper areas of these glacier systems are essentially balanced (ratio ∼ −0.1). IFL glaciers along the western and northern coastlines have  imbalance ratios similar to the major ISL glaciers, ∼ −0.5 to −2.4. Our mass balance estimate for the nAP region agrees well with recently published gravimetric values. Recent GRACEbased estimates that can be most easily compared with our study yield values of −27.5 ± 10 Gt a −1 (summing the mascons encompassing and adjacent to our study area) (Luthcke et al., 2013), and 26 ± 3 Gt a −1 for a larger GRACE mascon extending to 70 • S (Sasgen et al., 2013). Due to the low spatial resolution of the gravimetric measurement, both these GRACE-derived results inherently include portions of the Larsen C ice shelf and adjacent ice-covered islands we did not measure (notably, King George Island) that lost elevation and mass during the ICESat period (Gardner et al., 2013). Similarly, the strong east-west gradient revealed in our study is not discernable by the GRACE system. Overall, however, the GRACE results provide a good summary confirmation for our study, and imply that nearly all of the net mass loss for the Antarctic Peninsula lies in the nAP region defined here.
For earlier ICESat-only studies of the mass balance in the area (Shepherd et al., 2012), the apparent agreement is likely fortuitous. Simple extrapolation methods that do not include information about individual basin dynamics (e.g., spatial/elevation extent, ice shelf loss, east-west variations), lead to very different values for total mass change. We conducted two experiments using only our cross-track-adjusted ICESat data to examine the scale of possible discrepancies. With an assumption of uniform elevation change for each elevation band throughout the nAP, the volume change from ICESat data into two sub-sets would be −36.6 km 3 a −1 . This overestimate derives from ISL glaciers forming too great a part of the net elevation change measurement data, especially for their lower elevations. This is, in part, due to more ICE-Sat data being acquired along the eastern nAP, likely a result of less cloud cover there. If one partially addresses this by separating ICESat data in two sub-sets (ISL basins vs. the rest), the volume change is still 10 % greater than our study, −30.6 km 3 a −1 .
The most recent assessment of the mass balance of the entire peninsula uses CryoSat-2 interferometric radar altimetry data to infer a mass balance of −13 ± 13 Gt a −1 in the most comparable basins of their study for a period following our evaluation, 2010(McMillan et al., 2014their basins 25 and 26). We suggest that the large difference between this study and ours is due to underrepresentation of the narrow deep fjord glaciers on both sides that are rapidly losing ice elevation. This is inherent to the use of a radar altimeter over rough terrain: the system unavoidably oversamples high areas within the beam footprint. However, detailed study of a set of glacier outlets that formerly fed the Larsen A and Prince Gustav ice shelves (which were the site of major ice shelf disintegrations in 1988 and 1995) suggests some parts of our study area have begun to see a slightly reduced level of ice mass loss in the 2011-2013 period (Rott et al., 2014).
We now examine the potential impact of further ice shelf loss in the Scar Inlet region, a remnant ice shelf section from the Larsen B ice shelf. Comparing high-resolution bathymetric mapping of the seabed exposed by nAP-wide ice shelf loss and glacier retreat with our data in Fig. 1 shows that it is the glaciers with deep (> 500 m) troughs and recent ice shelf loss that have the greatest elevation loss and mass imbalance (Zgur et al., 2007;Shuman et al., 2011;Rebesco et al., 2014). Recent ice-thickness maps of the tributary glaciers (Starbuck, Flask, and Leppard glaciers) of the still-intact Scar Inlet ice shelf (SIIS) indicate they have unusually deep glacier troughs just behind the grounding line, well in excess of 1000 m below sea level in the case of Flask Glacier, and 500 m below sea level for Starbuck Glacier (Farinotti et al., 2013(Farinotti et al., , 2014. From Table S3, the mean imbalance ratio of ISL glaciers with ice-front bathymetric troughs exceeding 500 m depth is −1.20, and −3.18 for those exceeding 1000 m depth (for comparison, it is +0.07 for trough areas less that 500 m depth). If we assume that the three primary tributary glaciers of SIIS will experience the same mean imbalance ratio following a collapse of their frontal ice shelf in Scar Inlet, we can anticipate increased mass imbalance in those basins, from the −1.36 Gt a −1 observed during our study period to ∼ −5.5 Gt a −1 .

Conclusions
Overall, our study suggests that the nAP mass imbalance pattern is a combination of several recent changes to the coastal glaciers and ice shelf systems, likely beginning several decades ago along the western coastal fjords and islands, with extensive inland propagation of mass loss to the ice divide area, and more recent ice shelf loss along the eastern flanks and islands with extensive and expanding inland propagation. Further, the large measured increase in snow accumulation over the past few decades has not created vast regions of positive mass balance suggesting that negative mass balances will continue into the future.
The Supplement related to this article is available online at doi:10.5194/tc-8-2135-2014-supplement.
Author contributions. T. Scambos led the writing and compilation of graphics and tables, and with T. Haran and J. Bohlander conducted the ICESat-based elevation change analysis. E. Berthier conducted the differential DEM analysis and integrated the ICESat data with the dDEM data. C. Shuman and A. Cook evaluated glacier front area changes and C. Shuman produced components of the spatial analysis. A. Cook provided the glacier basin outlines. S. Ligtenberg evaluated the firn compaction and accumulation variability estimates, and their impact on our results. All co-authors contributed to the writing of the paper.