TY - JOUR
AU - Gunter, Brian C.
AB - Abstract The measurement of ongoing ice-mass loss and associated melt water contribution to sea-level change from regions such as West Antarctica is dependent on a combination of remote sensing methods. A key method, the measurement of changes in Earth's gravity via the GRACE satellite mission, requires a potentially large correction to account for the isostatic response of the solid Earth to ice-load changes since the Last Glacial Maximum. In this study, we combine glacial isostatic adjustment modelling with a new GPS dataset of solid Earth deformation for the southern Antarctic Peninsula to test the current understanding of ice history in this region. A sufficiently complete history of past ice-load change is required for glacial isostatic adjustment models to accurately predict the spatial variation of ongoing solid Earth deformation, once the independently-constrained effects of present-day ice mass loss have been accounted for. Comparisons between the GPS data and glacial isostatic adjustment model predictions reveal a substantial misfit. The misfit is localized on the southwestern Weddell Sea, where current ice models under-predict uplift rates by approximately 2 mm yr−1. This under-prediction suggests that either the retreat of the ice sheet grounding line in this region occurred significantly later in the Holocene than currently assumed, or that the region previously hosted more ice than currently assumed. This finding demonstrates the need for further fieldwork to obtain direct constraints on the timing of Holocene grounding line retreat in the southwestern Weddell Sea and that GRACE estimates of ice sheet mass balance will be unreliable in this region until this is resolved. Sea level change, Space geodetic surveys, Global change from geodesy, Glaciology, Antarctica 1 INTRODUCTION Ice lost by the West Antarctic Ice Sheet, predominantly through the changing flux of ice streams into the ocean, makes a significant and increasing contribution to global sea level rise (King et al.2012; Shepherd et al.2012). Within West Antarctica, the Antarctic Peninsula (AP) is experiencing some of the largest climatic change found anywhere on Earth (Abram et al.2013) and it is therefore important to understand how the AP ice sheet is responding. Ice sheet response can be estimated using a range of remote sensing techniques, including measurements of changes in Earth's gravity field, such as from the Gravity Recovery and Climate Experiment (GRACE, e.g. Wahr et al.2000), and satellite altimetry (e.g. Wouters et al.2015). In addition to the signal from changes in ice mass, these measurements contain a signal due to changes in the solid Earth, which must be subtracted; this is especially important for GRACE (Bentley & Wahr 1998). The dominant cause of solid Earth motion on timescales relevant to estimates of ice-mass change is glacial isostatic adjustment (GIA), a product of the relatively long relaxation time of Earth's interior in response to past changes in surface loading. Estimates of ongoing deformation due to GIA (uplift rates) can be produced by calibrating numerical GIA models using data such as relative sea level (RSL) indicators or GPS measurements of Earth surface motion (e.g. Lambeck & Chappell 2001; Milne et al.2001). GIA models comprise an assumed Earth rheology and a reconstructed ice history. The ice history reconstruction must be reasonably accurate for the aforementioned calibration exercise to be effective. Palmer Land in the southern AP (Fig. 1) has undergone substantial ice loading changes since the Last Glacial Maximum (LGM), the time period of interest for modelling present-day GIA (Bentley et al.2006; Whitehouse et al.2012a). Despite the importance of the region, only low spatial and temporal resolution observations are available to constrain ice extent and retreat history in this region (Whitehouse et al.2012a). Regional high-resolution studies addressing the solid Earth response to ice-load changes during the last ∼150 yr have increased understanding of GIA in the AP (e.g. Nield et al.2012, 2014), but the longer-term ice history, particularly of Palmer Land, remains very poorly known. Further issues include the possibility for climate-driven ice thickness oscillations during the last 1–2 kyr, (e.g. Bertler et al.2011) and uncertainties in ice thickness change around major ice streams. The differences in reconstructed ice history express themselves as differences in modelled GIA. We discuss further details of the origin and uncertainties of various ice history reconstructions in Section 2. Figure 1. Open in new tabDownload slide Location of GPS stations in Palmer Land, red triangles are CAPGIA sites, open circles are other sites as described in the text. Locations of pre-existing data constraints are shown: the blue ellipse indicates the region in Marguerite Bay where the RSL data of Bentley et al. (2005) were collected. Yellow stars indicate the distribution of Holocene cosmogenic exposure ages (Bentley et al.2006; Hodgson et al.2009; Bentley et al.2011; Johnson et al.2012). Regions referred to by name in the text are indicated. Present-day ice shelf fronts are shown with a grey line. Figure 1. Open in new tabDownload slide Location of GPS stations in Palmer Land, red triangles are CAPGIA sites, open circles are other sites as described in the text. Locations of pre-existing data constraints are shown: the blue ellipse indicates the region in Marguerite Bay where the RSL data of Bentley et al. (2005) were collected. Yellow stars indicate the distribution of Holocene cosmogenic exposure ages (Bentley et al.2006; Hodgson et al.2009; Bentley et al.2011; Johnson et al.2012). Regions referred to by name in the text are indicated. Present-day ice shelf fronts are shown with a grey line. Data that can directly constrain GIA model predictions for Palmer Land are limited. The only RSL data set available is located in Marguerite Bay (Fig. 1; Bentley et al.2005) and rates from only four GPS stations have been published previously (Argus et al.2014). As such, GIA models for this region are uncertain, with different models producing substantially different patterns of uplift rates. For GRACE-based estimates of ice-mass change, this represents systematic uncertainty (King et al.2012). It has been observed that GIA models often overpredict GPS uplift rate measurements (e.g. Thomas et al.2011). This previous poor performance of GIA models in Palmer Land is most probably due to rheological uncertainty, poor knowledge of ice history and model resolution limitations when applied to such a small region. Recent model developments have produced a new set of predictions for this region (compare: Whitehouse et al.2012a; Ivins et al.2013; Argus et al.2014), although these are yet to be fully tested by independent data. In this work, we improve the GPS data coverage by analysing uplift rates measured at nine additional sites in Palmer Land (Fig. 1) and by modelling GIA at high resolution using a range of plausible Earth rheology parameters. 1.1 Geological background to GIA modelling The AP has a complex geological history as it has been on the Antarctic–Pacific Plate boundary since the late Precambrian (Barker 2001). There is a considerable literature on the tectonics of the last 100 Myr of the Antarctic Peninsula (e.g. Johnson & Swain 1995; McCarron & Larter 1998; Eagles et al.2009), but for our purposes it is sufficient to note that the margin of the AP is currently passive except for a small amount of subduction and associated backarc spreading to the north of our region, near the South Shetland trench. Subduction shut down in Palmer Land before 30 Ma (Johnson & Swain 1995). From this brief history, it is already clear that deformation on the AP region cannot be appropriately modelled with the same Earth rheology parameters as the generally much older and thicker Early Precambrian East Antarctic shield (Adie 1962; Morelli & Danesi 2004). When modelling viscoelastic deformation, the lithosphere thickness (LT) specified is the ‘effective elastic thickness’ of the lithosphere. This definition of the lithosphere means that the LT value that provides the best fit to a particular geodetic data set may not exactly match lithosphere thickness values derived by alternative methods. For example, the seismically-imaged thermal lithosphere of the East Antarctic Shield could be around 220 km thick (Morelli & Danesi 2004) but previous GIA modelling has found that elastic thicknesses between 65 and 120 km produce good fits to a range of data from across Antarctica (Whitehouse et al.2012a; Ivins et al.2013). The cause of this apparent discrepancy between geophysical observations and GIA model fit is unclear, but may be due to depth dependant viscosity within the lithosphere (Kuchar & Milne 2015). The range of effective elastic thickness we consider in this work is within the range of values widely considered in global GIA modelling (e.g. Lambeck et al.2014). The viscosity profile of Earth's mantle is uncertain and lateral variations are likely to be significant. The globally averaged upper mantle viscosity (UMV) is thought to lie in the range 1020–1021 Pa s (Mitrovica & Forte 2004), while the lower mantle viscosity (LMV) is even more uncertain. Lambeck et al. (2014) found that globally-averaged LMV values of 2 × 1021 and 1023 Pa s could both satisfy their global data constraints, although regional variation is likely. A fairly general study in the Antarctic region by Li & Yuen (2014) produced an estimate of 1022–1.6 × 1023 Pa s for the LMV. There is some evidence of higher upper mantle temperatures beneath the AP and it is often assumed that mantle viscosities are likely to be lower than average in this region (Morelli & Danesi 2004; Shapiro & Ritzwoller 2004; Maule et al.2005), although limitations on seismic tomography resolution in the region pose a significant problem. In addition to the uncertainties outlined above, implementing lithosphere thickness or mantle viscosity directly from the wider geophysical literature in a GIA model is not appropriate due to the necessarily simplified Maxwell rheology assumed within a GIA model. Instead, we use the information above to guide our choice for the range of parameters that can be considered plausible, given the poor state of knowledge of the rheology of the region. In this paper, we first describe the GPS data collection and processing, GIA modelling methods and our data-model fitting algorithm. We then report the GPS results with corrections for elastic deformation and tectonic plate rotation, followed by output from the GIA modelling. Finally, the residuals between GPS data and GIA model output are examined and their implications discussed. 2 METHODS 2.1 GPS During late 2009 and early 2010 nine continuous GPS sites were installed in southern Palmer Land within the framework of the Constraints on Antarctic Peninsula GIA (CAPGIA) project (triangles in Fig. 1). Aside from FOS1, all sites required the installation of new monumentation; the FOS1 monument was installed in 1995 as part of campaign-style GPS measurements in 1995, 1996 and 1998 (Dietrich et al.1996; Dietrich et al.2004). Part of the 2009/2010 CAPGIA FOS1 record was appended to those earlier data by Thomas et al. (2011). We supplement the CAPGIA sites with data from other sites in the region (circles, Fig. 1), namely BREN (Thomas et al.2011), SUGG and HAAG (Argus et al.2014). All sites sample the GPS signals every 30 s. Powering of remote equipment through the Antarctic winter remains a challenge and data gaps exist within each record, notably over April to September, although in many cases near year-round operation has been demonstrated. For all CAPGIA sites (including FOS1) we use data from late 2009 until the most recent site download in late 2013 or early 2014. Earlier campaign data exist for FOS1 but they are subject to different elastic corrections, which are less well known than for 2009–2014; they may also be subject to unknown offsets due to equipment changes. The BREN antenna was changed in 2012 and hence has an offset in the time-series at that point. As such we only make use of data up to a data gap in mid 2010 for BREN noting that the data used to calculate the elastic correction at this site were specific to the BREN epoch. Further details of the record from each site are provided in Table 1. Readers not interested in the specific details of the GPS analysis may wish to advance to Section 2.2. Table 1. Site names (4-character site IDs), locations, information for data analysed within this study after outlier removal, vertical GPS rates, errors and elastic correction. DD denotes decimal degrees. U is the vertical elastic correction, σU is the uncertainty on the vertical elastic correction, errors and uncertainties are 1σ. . . . . . . . . . . Elastic . . . . . . . . . . . . . . . . . . . . . . . Site . Lat . Lon . Ellipsoidal . Start . End . Data . RINEX data . GPS rate . Uncertainty . . . (DD) . (DD) . Height (m) . (YYYY DDD) . (YYYY DDD) . days . DOI . (mm yr−1) . (mm yr−1) . U . σU . BEAN − 75.96 290.70 907.00 2010 007 2012 015 397 10.7283/T55Q4T6R 7.54 1.24 1.24 0.82 FOS1 − 71.31 291.68 158.83 2009 352 2013 123 903 10.7283/T54T6GF7 3.89 1.07 1.87 0.75 GMEZ − 73.89 291.46 1417.00 2010 008 2013 314 731 10.7283/T58G8HT4 5.74 0.79 2.38 0.92 HTON − 74.08 298.27 949.00 2010 005 2013 328 413 10.7283/T5222RV6 6.36 0.87 1.08 0.61 JNSN − 73.08 293.90 1365.00 2009 364 2013 314 557 10.7283/T5SJ1HP1 5.33 0.86 2.26 0.86 LNTK − 74.84 286.10 1197.00 2009 365 2013 048 259 10.7283/T5J1017P 6.03 0.66 2.41 0.88 MKIB − 75.28 294.40 1155.00 2010 006 2013 320 479 10.7283/T5D798HD 6.92 0.50 1.37 0.77 TRVE − 69.99 292.45 1047.00 2009 356 2013 310 457 10.7283/T5NS0RZ9 4.67 0.63 2.66 0.78 WLCH − 70.73 296.18 1526.00 2010 010 2013 309 349 10.7283/T5×928DR 2.64 0.68 1.31 0.79 BREN − 72.67 296.97 1646.14 2006 363 2010 195 484 10.7283/T52V2D7X 3.22 0.77 0.59 0.81 HAAG − 77.04 281.71 1172.60 – – – – 7.49 2.15 1.75 0.76 SUGG − 75.28 287.82 1090.40 – – – – 4.5 1.85 1.76 0.88 . . . . . . . . . . Elastic . . . . . . . . . . . . . . . . . . . . . . . Site . Lat . Lon . Ellipsoidal . Start . End . Data . RINEX data . GPS rate . Uncertainty . . . (DD) . (DD) . Height (m) . (YYYY DDD) . (YYYY DDD) . days . DOI . (mm yr−1) . (mm yr−1) . U . σU . BEAN − 75.96 290.70 907.00 2010 007 2012 015 397 10.7283/T55Q4T6R 7.54 1.24 1.24 0.82 FOS1 − 71.31 291.68 158.83 2009 352 2013 123 903 10.7283/T54T6GF7 3.89 1.07 1.87 0.75 GMEZ − 73.89 291.46 1417.00 2010 008 2013 314 731 10.7283/T58G8HT4 5.74 0.79 2.38 0.92 HTON − 74.08 298.27 949.00 2010 005 2013 328 413 10.7283/T5222RV6 6.36 0.87 1.08 0.61 JNSN − 73.08 293.90 1365.00 2009 364 2013 314 557 10.7283/T5SJ1HP1 5.33 0.86 2.26 0.86 LNTK − 74.84 286.10 1197.00 2009 365 2013 048 259 10.7283/T5J1017P 6.03 0.66 2.41 0.88 MKIB − 75.28 294.40 1155.00 2010 006 2013 320 479 10.7283/T5D798HD 6.92 0.50 1.37 0.77 TRVE − 69.99 292.45 1047.00 2009 356 2013 310 457 10.7283/T5NS0RZ9 4.67 0.63 2.66 0.78 WLCH − 70.73 296.18 1526.00 2010 010 2013 309 349 10.7283/T5×928DR 2.64 0.68 1.31 0.79 BREN − 72.67 296.97 1646.14 2006 363 2010 195 484 10.7283/T52V2D7X 3.22 0.77 0.59 0.81 HAAG − 77.04 281.71 1172.60 – – – – 7.49 2.15 1.75 0.76 SUGG − 75.28 287.82 1090.40 – – – – 4.5 1.85 1.76 0.88 Open in new tab Table 1. Site names (4-character site IDs), locations, information for data analysed within this study after outlier removal, vertical GPS rates, errors and elastic correction. DD denotes decimal degrees. U is the vertical elastic correction, σU is the uncertainty on the vertical elastic correction, errors and uncertainties are 1σ. . . . . . . . . . . Elastic . . . . . . . . . . . . . . . . . . . . . . . Site . Lat . Lon . Ellipsoidal . Start . End . Data . RINEX data . GPS rate . Uncertainty . . . (DD) . (DD) . Height (m) . (YYYY DDD) . (YYYY DDD) . days . DOI . (mm yr−1) . (mm yr−1) . U . σU . BEAN − 75.96 290.70 907.00 2010 007 2012 015 397 10.7283/T55Q4T6R 7.54 1.24 1.24 0.82 FOS1 − 71.31 291.68 158.83 2009 352 2013 123 903 10.7283/T54T6GF7 3.89 1.07 1.87 0.75 GMEZ − 73.89 291.46 1417.00 2010 008 2013 314 731 10.7283/T58G8HT4 5.74 0.79 2.38 0.92 HTON − 74.08 298.27 949.00 2010 005 2013 328 413 10.7283/T5222RV6 6.36 0.87 1.08 0.61 JNSN − 73.08 293.90 1365.00 2009 364 2013 314 557 10.7283/T5SJ1HP1 5.33 0.86 2.26 0.86 LNTK − 74.84 286.10 1197.00 2009 365 2013 048 259 10.7283/T5J1017P 6.03 0.66 2.41 0.88 MKIB − 75.28 294.40 1155.00 2010 006 2013 320 479 10.7283/T5D798HD 6.92 0.50 1.37 0.77 TRVE − 69.99 292.45 1047.00 2009 356 2013 310 457 10.7283/T5NS0RZ9 4.67 0.63 2.66 0.78 WLCH − 70.73 296.18 1526.00 2010 010 2013 309 349 10.7283/T5×928DR 2.64 0.68 1.31 0.79 BREN − 72.67 296.97 1646.14 2006 363 2010 195 484 10.7283/T52V2D7X 3.22 0.77 0.59 0.81 HAAG − 77.04 281.71 1172.60 – – – – 7.49 2.15 1.75 0.76 SUGG − 75.28 287.82 1090.40 – – – – 4.5 1.85 1.76 0.88 . . . . . . . . . . Elastic . . . . . . . . . . . . . . . . . . . . . . . Site . Lat . Lon . Ellipsoidal . Start . End . Data . RINEX data . GPS rate . Uncertainty . . . (DD) . (DD) . Height (m) . (YYYY DDD) . (YYYY DDD) . days . DOI . (mm yr−1) . (mm yr−1) . U . σU . BEAN − 75.96 290.70 907.00 2010 007 2012 015 397 10.7283/T55Q4T6R 7.54 1.24 1.24 0.82 FOS1 − 71.31 291.68 158.83 2009 352 2013 123 903 10.7283/T54T6GF7 3.89 1.07 1.87 0.75 GMEZ − 73.89 291.46 1417.00 2010 008 2013 314 731 10.7283/T58G8HT4 5.74 0.79 2.38 0.92 HTON − 74.08 298.27 949.00 2010 005 2013 328 413 10.7283/T5222RV6 6.36 0.87 1.08 0.61 JNSN − 73.08 293.90 1365.00 2009 364 2013 314 557 10.7283/T5SJ1HP1 5.33 0.86 2.26 0.86 LNTK − 74.84 286.10 1197.00 2009 365 2013 048 259 10.7283/T5J1017P 6.03 0.66 2.41 0.88 MKIB − 75.28 294.40 1155.00 2010 006 2013 320 479 10.7283/T5D798HD 6.92 0.50 1.37 0.77 TRVE − 69.99 292.45 1047.00 2009 356 2013 310 457 10.7283/T5NS0RZ9 4.67 0.63 2.66 0.78 WLCH − 70.73 296.18 1526.00 2010 010 2013 309 349 10.7283/T5×928DR 2.64 0.68 1.31 0.79 BREN − 72.67 296.97 1646.14 2006 363 2010 195 484 10.7283/T52V2D7X 3.22 0.77 0.59 0.81 HAAG − 77.04 281.71 1172.60 – – – – 7.49 2.15 1.75 0.76 SUGG − 75.28 287.82 1090.40 – – – – 4.5 1.85 1.76 0.88 Open in new tab We processed undifferenced raw data for BREN and the CAPGIA sites in GIPSY/OASIS v6.2 using a largely standard approach holding fixed JPL's reprocessed (as of early 2014) non-fiducial satellite clock and orbit products. Raw GPS pseudo-range data were smoothed using the carrier phase data and, together with the phase data, were subsampled to 5-min intervals; the ionosphere-free linear combination was then formed. We did not model the higher order ionospheric terms as they are not included within the analysis that generated the satellite products. We mapped the slant tropospheric delays to the zenith using the gridded Vienna Mapping Functions v1 (VMF1; Boehm et al.2006) with a priori hydrostatic zenith delays based on ECMWF/VMF1 grids (Tregoning & Herring 2006). We applied antenna phase centre models as specified within the igs08_1781.atx file, modified to include the relatively minor effects of the POLENET monument (King et al.2011a) installed at all CAPGIA sites other than FOS1. We modelled elastic effects associated with the pole tide, Earth body tide and ocean tide loading according to the IERS2010 Conventions. The ocean tide loading displacements were computed in the centre of mass of the entire Earth system (CM) frame (Blewitt 2003) based on the TPXO7.2 ocean tide model (Egbert et al.2009) using the SPOTL software (Agnew 1997). Ocean tide loading displacements in our region are dominated by tides in the Weddell Sea. The TPXO7.2 model was chosen because a comparison to independent GPS measurements of ice shelf tidal motion showed it to be the most accurate of the models tested (King et al.2011b,c). We did not correct for time-variable elastic effects associated with recent polar motion as they are small at these locations (King & Watson 2014). We estimated site coordinates once per day, alongside tropospheric zenith delays and horizontal gradients every 5 min, receiver clocks at each measurement epoch, and real-valued ambiguity parameters. Ambiguities were fixed to integers using JPL's tabulated wide-lane phase biases using the approach of Bertiger et al. (2010). Typical GIPSY analysis involves weighting all data uniformly but this does not reflect the error distribution of typical GPS observations for which scatter tends to decrease with increasing satellite elevation angle. Instead, we performed an iterative analysis similar to that widely used in GAMIT software (Herring et al.2010), whereby an initial analysis with uniform weighting was performed and then an elevation-dependent model was fit to each of the phase and code residuals. We then repeated the analysis using the modified weighting models to produce the final parameter estimates. We transformed daily site coordinates to an International GNSS Service (IGS) realisation of the International Terrestrial Reference Frame 2008 (ITRF2008; Altamimi et al.2011) using daily 7-parameter transformations provided by JPL; before October 2012 the final frame was IGS08 and after this it was IGb08 but with negligible offset. We removed the effects of atmospheric loading displacements from the time-series using interpolated Global Geophysical Fluids Centre grids in the centre of Earth's figure (CF) frame, as is consistent with non-secular motions in frames such as ITRF2008 (Dong et al.1997). We identified time-series outliers by fitting a model comprised of constant, linear, annual, and semi-annual terms for each coordinate component (Herring 2003); there are no time-series offsets within the data periods we use here. Daily estimates more than three standard deviations from the model in any one component were regarded as outliers and removed. Some of our sites exhibit rapid transient behaviour during the wintertime that is most evident in the horizontal coordinate components. Similar transients have been identified in previous studies (e.g. Jaldehag et al.1996; Johansson et al.2002) and attributed to frost or rime ice, which may form on the side of the antenna radomes in autumn and winter. Following the processing approach of Argus et al. (2014), we remove the affected time intervals from our GPS data before analysis. The only sites where such transient behaviour is not evident are FOS1 and BREN—we note that FOS1 is much closer to sea level than the other sites and BREN does not have a radome—so the observations are consistent with the rime ice hypothesis but the source of this signal clearly requires further investigation. After removal of outliers, we computed three-component site velocities and uncertainties using the CATSv3.20 software (Williams 2008), estimating constant, annual, semi-annual and velocity terms and adopting a time-variable white plus power-law noise model. We estimated the power-law spectral index (Agnew 1992) of the sites for each coordinate component, with a mean value of −0.81 ± 0.22 for the vertical coordinates and −0.25 ± 0.09 for the horizontal coordinates. Velocities in ITRF2008 have their origin at an approximation to CM over the period 1983.0–2009.0. For SUGG and HAAG we did not analyse raw data but instead adopted the Argus et al. (2014) velocities and uncertainties, using their tabulated frame translation parameters and plate rotations to shift their horizontal and vertical components back into ITRF2008. 2.2 Elastic modelling We modelled the effects of elastic deformation of the solid Earth due to ice load changes occurring over the GPS data collection epoch. For sites that have GPS data in the 2010–2014 epoch (all sites including SUGG and HAAG but excluding BREN), we adopted 5 km gridded ice elevation trends for West Antarctica and Alexander Island based on a combination of CryoSat-2 (CS2) data (November 2010–September 2013 inclusive) as described by McMillan et al. (2014) and the detailed loading changes of Nield et al. (2014) for the region north of Marguerite Bay. We extrapolated the CS2 data to the grounding line, although this has negligible effect on our solutions. CS2 mass trends and uncertainties were derived using the method of McMillan et al. (2014), whereby CS2 elevation trends were converted to mass trends using a density model to separate regions changing mass at a density of snow or ice. In computing ice mass trends, snowfall fluctuations over the observational period complicate the volume to mass conversion, and affect the capacity to isolate long term trends in ice mass imbalance. To account for this source of uncertainty within the mass trends, an estimate of density model uncertainty is not explicitly made, but instead an additional uncertainty term is computed from the contemporary snowfall fluctuation and added to the elevation rate uncertainty (McMillan et al.2014). To calculate the elastic response to load changes we adopted the method and elastic Earth structure of Nield et al. (2014). Ensuring conservation of mass, we expanded the surface load to spherical harmonic degree 3700. The modelled magnitudes of the three elastic velocity components at the GPS sites are listed in Tables 1 and 2, with the spatial pattern of the elastic uplift rates shown in Fig. S1. We subtracted these modelled rates from the estimated GPS rates. Table 2. Horizontal GPS rates following the application of two different models of plate rotation, GPS horizontal uncertainty, and horizontal elastic corrections. N and E denote the north- and east-directed horizontal vectors. The uncertainties in the horizontal elastic correction are negligible and are not included. Euler poles for plate motion correction are: Altamimi et al. (2012) 58.5455, −129.843, 0.2090. Argus (personal communication Don Argus, 2015): 59.876, −127.277, 0.2178 (lat (deg), lon (deg), rate (deg/Ma)). Site . GPS rates with plate rotation removed . . . . . . . . . . . . . GPS uncertainty . Elastic correction . . Altamimi Model . Argus Model . (1σ) . . . N . E . N . E . Nerr . Eerr . N . E . BEAN 0.46 0.41 0.72 − 0.33 0.31 0.36 − 0.11 0.28 FOS1 − 0.29 1.51 − 0.04 0.69 0.21 0.21 − 0.02 0.04 GMEZ − 1.47 2.95 − 1.22 2.18 0.20 0.27 0.08 0.2 HTON − 0.1 1.13 0.09 0.33 0.21 0.22 − 0.19 0.42 JNSN − 0.67 2.45 − 0.44 1.65 0.48 0.32 − 0.04 0.24 LNTK − 0.77 0.62 − 0.47 − 0.12 0.21 0.21 − 0.01 0.4 MKIB 0.32 1.47 0.54 0.71 0.23 0.44 − 0.27 0.31 TRVE − 2.75 1.1 − 2.51 0.25 0.21 0.27 − 0.11 0.3 WLCH − 1.13 0.68 − 0.92 − 0.17 0.21 0.18 − 0.02 0.1 BREN 0.05 1.09 0.25 0.27 0.18 0.17 − 0.01 0.24 HAAG 0.49 − 0.27 0.81 − 0.93 0.40 0.85 − 0.22 0.39 SUGG − 0.66 0 − 0.39 − 0.71 0.55 0.65 − 0.26 0.33 Site . GPS rates with plate rotation removed . . . . . . . . . . . . . GPS uncertainty . Elastic correction . . Altamimi Model . Argus Model . (1σ) . . . N . E . N . E . Nerr . Eerr . N . E . BEAN 0.46 0.41 0.72 − 0.33 0.31 0.36 − 0.11 0.28 FOS1 − 0.29 1.51 − 0.04 0.69 0.21 0.21 − 0.02 0.04 GMEZ − 1.47 2.95 − 1.22 2.18 0.20 0.27 0.08 0.2 HTON − 0.1 1.13 0.09 0.33 0.21 0.22 − 0.19 0.42 JNSN − 0.67 2.45 − 0.44 1.65 0.48 0.32 − 0.04 0.24 LNTK − 0.77 0.62 − 0.47 − 0.12 0.21 0.21 − 0.01 0.4 MKIB 0.32 1.47 0.54 0.71 0.23 0.44 − 0.27 0.31 TRVE − 2.75 1.1 − 2.51 0.25 0.21 0.27 − 0.11 0.3 WLCH − 1.13 0.68 − 0.92 − 0.17 0.21 0.18 − 0.02 0.1 BREN 0.05 1.09 0.25 0.27 0.18 0.17 − 0.01 0.24 HAAG 0.49 − 0.27 0.81 − 0.93 0.40 0.85 − 0.22 0.39 SUGG − 0.66 0 − 0.39 − 0.71 0.55 0.65 − 0.26 0.33 Open in new tab Table 2. Horizontal GPS rates following the application of two different models of plate rotation, GPS horizontal uncertainty, and horizontal elastic corrections. N and E denote the north- and east-directed horizontal vectors. The uncertainties in the horizontal elastic correction are negligible and are not included. Euler poles for plate motion correction are: Altamimi et al. (2012) 58.5455, −129.843, 0.2090. Argus (personal communication Don Argus, 2015): 59.876, −127.277, 0.2178 (lat (deg), lon (deg), rate (deg/Ma)). Site . GPS rates with plate rotation removed . . . . . . . . . . . . . GPS uncertainty . Elastic correction . . Altamimi Model . Argus Model . (1σ) . . . N . E . N . E . Nerr . Eerr . N . E . BEAN 0.46 0.41 0.72 − 0.33 0.31 0.36 − 0.11 0.28 FOS1 − 0.29 1.51 − 0.04 0.69 0.21 0.21 − 0.02 0.04 GMEZ − 1.47 2.95 − 1.22 2.18 0.20 0.27 0.08 0.2 HTON − 0.1 1.13 0.09 0.33 0.21 0.22 − 0.19 0.42 JNSN − 0.67 2.45 − 0.44 1.65 0.48 0.32 − 0.04 0.24 LNTK − 0.77 0.62 − 0.47 − 0.12 0.21 0.21 − 0.01 0.4 MKIB 0.32 1.47 0.54 0.71 0.23 0.44 − 0.27 0.31 TRVE − 2.75 1.1 − 2.51 0.25 0.21 0.27 − 0.11 0.3 WLCH − 1.13 0.68 − 0.92 − 0.17 0.21 0.18 − 0.02 0.1 BREN 0.05 1.09 0.25 0.27 0.18 0.17 − 0.01 0.24 HAAG 0.49 − 0.27 0.81 − 0.93 0.40 0.85 − 0.22 0.39 SUGG − 0.66 0 − 0.39 − 0.71 0.55 0.65 − 0.26 0.33 Site . GPS rates with plate rotation removed . . . . . . . . . . . . . GPS uncertainty . Elastic correction . . Altamimi Model . Argus Model . (1σ) . . . N . E . N . E . Nerr . Eerr . N . E . BEAN 0.46 0.41 0.72 − 0.33 0.31 0.36 − 0.11 0.28 FOS1 − 0.29 1.51 − 0.04 0.69 0.21 0.21 − 0.02 0.04 GMEZ − 1.47 2.95 − 1.22 2.18 0.20 0.27 0.08 0.2 HTON − 0.1 1.13 0.09 0.33 0.21 0.22 − 0.19 0.42 JNSN − 0.67 2.45 − 0.44 1.65 0.48 0.32 − 0.04 0.24 LNTK − 0.77 0.62 − 0.47 − 0.12 0.21 0.21 − 0.01 0.4 MKIB 0.32 1.47 0.54 0.71 0.23 0.44 − 0.27 0.31 TRVE − 2.75 1.1 − 2.51 0.25 0.21 0.27 − 0.11 0.3 WLCH − 1.13 0.68 − 0.92 − 0.17 0.21 0.18 − 0.02 0.1 BREN 0.05 1.09 0.25 0.27 0.18 0.17 − 0.01 0.24 HAAG 0.49 − 0.27 0.81 − 0.93 0.40 0.85 − 0.22 0.39 SUGG − 0.66 0 − 0.39 − 0.71 0.55 0.65 − 0.26 0.33 Open in new tab To capture the impact of the uncertainty in the mass trend data we computed 1000 perturbations to the mass trend field, with each perturbation governed by the uncertainties in the data. To account for spatial correlation in mass trend uncertainties (due mainly to the spatial correlation in the density model) we adopted a correlation length of 100 km. Neglecting this spatial correlation produces negligible elastic uncertainties. 100 km is about half the width of Palmer Land so we take this to be a conservative correlation distance. For each perturbation we recomputed the elastic deformation at each site, from which the standard deviation was computed (Table 1) and propagated to the GPS uncertainties. This increased the uncertainty of the measured uplift rate at the GPS sites by about 0.4 mm yr−1 with negligible increases for the horizontal components. The different epoch of the GPS data used for BREN (2007–2010) required a different elastic correction for that site as CS2 data were not available before November 2010. The estimate we produce is based on ICESat data (2003–2009) and a surface mass balance and firn densification model derived from RACMO2 (Ligtenberg et al.2011; Lenaerts et al.2012), using the process described by Gunter et al. (2014). We filtered the input fields with a Gaussian filter with a 6σ width of 100 km and replaced the load changes north of 69.5°S with those used by Nield et al. (2014) for the same period. Our BREN vertical elastic correction (0.59 mm yr−1) is 1.26 mm yr−1 smaller than the elastic term adopted by Thomas et al. (2011) and 0.26 mm yr−1 smaller than their alternative correction (ICESat-based, their supplementary material). This updated (ICESat-based) elastic model results in larger estimates of GIA-related uplift at BREN than when using the elastic correction of Thomas et al. (2011). We validated our elastic computations by comparing them with solutions based on SPOTL (Agnew 1997) using three different models of shallow Earth structure (reflecting oceanic, continental and average structure). Comparing the SPOTL results using the average structure with our estimates based on PREM revealed differences in uplift rates up to 0.16 mm yr−1 but typically smaller than ±0.1 mm yr−1. The differences due to varying shallow Earth structure were no more than 0.08 mm yr−1, although we did not explore a fully 3-D elastic structure (Dill et al.2015). These effects are <∼10 per cent of the modelled signal and we regard them as negligible for our purposes. In addition to modelling the elastic effect of near-field ice mass changes, we also examined the potentially non-negligible signal due to far-field elastic effects (associated with non-Antarctic mass changes), including the relative motion between the centre of mass and centre of figure of the Earth. To do this, we estimated a trend in surface mass change over continental areas outside Antarctica from GRACE data, which provide a global coverage. The resulting vertical deformation over Antarctica was then computed in a gravitationally self-consistent way (Clarke et al.2005). We used GRACE RL05 monthly solutions from CSR for the period 2007–2013 (http://icgem.gfz-potsdam.de/ICGEM/shms/monthly/csr-rl05/), post-processed by means of the DDK5 filter (Kusche et al.2009), with the C20 coefficient from Satellite Laser Ranging data (Cheng et al.2011) and degree 1 coefficients from the approach by Swenson et al. (2008). The mass trend due to GIA was removed, based on the solution of A et al. (2013) (http://grace.jpl.nasa.gov). The estimated signal magnitude was ∼–0.5 mm yr−1 across Palmer Land. However, due to the degree of uncertainty associated with the driving processes and the spatial uniformity of the signal across the AP, we did not correct for these components of deformation but consider them further when comparing data with GIA model predictions. 2.3 Plate rotation correction In addition to elastic effects, the horizontal components of the processed GPS velocities contain the tectonic motion of the Antarctic Plate. In order to recover the horizontal GIA rates, a model of plate rotation must be removed. In common with previous work we assume that the Antarctic Plate is rigid, including this southern region of the Antarctic Peninsula (Dietrich et al.2001; Argus et al.2014), and adopt two plate rotation models to produce two variations of horizontal velocities that we interpret as horizontal rates associated with GIA (Altamimi et al.2012; Argus et al.2014). We note that these geodetically derived plate rotation estimates may have been biased by an unmodelled horizontal signal due to GIA (Klemann et al.2008). 2.4 GIA modelling After GPS processing and correction for elastic deformation and plate rotation, the residual GPS deformation rates (with error) are assumed to represent the predominantly viscous isostatic re-equilibration of Earth's interior to past ice load changes, which we refer to as GIA. The GIA code used in this study (Milne & Mitrovica 1998; Whitehouse et al.2012a) models deformation of a Maxwell (viscoelastic) Earth with spherically symmetric rheology, divided into three layers (after: Peltier 1974): (1) the lithosphere is assigned a viscosity of 1043 Pa s, approximating purely elastic behaviour over an effective thickness (LT) defined in kilometres; (2) the UMV value is assigned to the region extending from the base of the lithosphere to 660 km depth; (3) the LMV value is assigned to the region below 660 km and above the core–mantle boundary. A lithosphere and two layers in the mantle are considered adequate given the apparent resolving power of GIA measurements in most regions (Paulson et al.2007). Density and velocity (elastic) parameters were adopted from PREM (Dziewonski & Anderson 1981). The suite of Earth models used in this study is presented in Table 3. The Earth model abbreviation format is: LT (km)/UMV (×1021 Pa s)/LMV(×1021 Pa s) and is used throughout the remainder of the text. Table 3. Earth model parameters used within the GIA model. *We used the model: 90/VM2/VM2, after Peltier (2004), which consists of more mantle layers than our other Earth models, resulting in a more graded change in viscosity through the mantle. This was done in order to run ICE-5G with its preferred Earth model. Projected onto a three-layer structure, 90/VM2/VM2 is approximately equivalent to 90/0.5/3. †Only run with W12. ‡W12 standard Earth model. LT (km) . UMV (Pa s) . LMV (Pa s) . Abbreviation . 46 1.00E + 19 1.00E + 22 46/0.01/10 46 5.00E + 19 1.00E + 21 46/0.05/1 46 1.00E + 20 1.00E + 22 46/0.1/10 46 3.00E + 20 1.00E + 22 46/0.3/10 46 5.00E + 20 1.00E + 21 46/0.5/1 46 5.00E + 20 1.00E + 22 46/0.5/10 46 2.00E + 21 1.00E + 22 46/2/10 71 1.00E + 20 1.00E + 22 71/0.1/10 71 3.00E + 20 1.00E + 22 71/0.3/10 90 5.00E + 20 3.00E + 21 *90/VM2/VM2 96 5.00E + 20 1.00E + 22 96/0.5/10 120 3.00E + 19 1.00E + 22 120/0.03/10 120 1.00E + 20 1.00E + 22 120/0.1/10 120 3.00E + 20 1.00E + 21 120/0.3/1 120 3.00E + 20 1.00E + 22 120/0.3/10 120 5.00E + 20 1.00E + 21 †120/0.5/1 120 5.00E + 20 5.00E + 21 120/0.5/5 120 5.00E + 20 1.00E + 22 120/0.5/10 120 5.00E + 20 1.50E + 22 120/0.5/15 120 5.00E + 20 2.00E + 22 †120/0.5/20 120 5.00E + 20 1.00E + 23 †120/0.5/100 120 1.00E + 21 1.00E + 22 ‡120/1/10 120 2.00E + 21 1.00E + 22 120/2/10 LT (km) . UMV (Pa s) . LMV (Pa s) . Abbreviation . 46 1.00E + 19 1.00E + 22 46/0.01/10 46 5.00E + 19 1.00E + 21 46/0.05/1 46 1.00E + 20 1.00E + 22 46/0.1/10 46 3.00E + 20 1.00E + 22 46/0.3/10 46 5.00E + 20 1.00E + 21 46/0.5/1 46 5.00E + 20 1.00E + 22 46/0.5/10 46 2.00E + 21 1.00E + 22 46/2/10 71 1.00E + 20 1.00E + 22 71/0.1/10 71 3.00E + 20 1.00E + 22 71/0.3/10 90 5.00E + 20 3.00E + 21 *90/VM2/VM2 96 5.00E + 20 1.00E + 22 96/0.5/10 120 3.00E + 19 1.00E + 22 120/0.03/10 120 1.00E + 20 1.00E + 22 120/0.1/10 120 3.00E + 20 1.00E + 21 120/0.3/1 120 3.00E + 20 1.00E + 22 120/0.3/10 120 5.00E + 20 1.00E + 21 †120/0.5/1 120 5.00E + 20 5.00E + 21 120/0.5/5 120 5.00E + 20 1.00E + 22 120/0.5/10 120 5.00E + 20 1.50E + 22 120/0.5/15 120 5.00E + 20 2.00E + 22 †120/0.5/20 120 5.00E + 20 1.00E + 23 †120/0.5/100 120 1.00E + 21 1.00E + 22 ‡120/1/10 120 2.00E + 21 1.00E + 22 120/2/10 Open in new tab Table 3. Earth model parameters used within the GIA model. *We used the model: 90/VM2/VM2, after Peltier (2004), which consists of more mantle layers than our other Earth models, resulting in a more graded change in viscosity through the mantle. This was done in order to run ICE-5G with its preferred Earth model. Projected onto a three-layer structure, 90/VM2/VM2 is approximately equivalent to 90/0.5/3. †Only run with W12. ‡W12 standard Earth model. LT (km) . UMV (Pa s) . LMV (Pa s) . Abbreviation . 46 1.00E + 19 1.00E + 22 46/0.01/10 46 5.00E + 19 1.00E + 21 46/0.05/1 46 1.00E + 20 1.00E + 22 46/0.1/10 46 3.00E + 20 1.00E + 22 46/0.3/10 46 5.00E + 20 1.00E + 21 46/0.5/1 46 5.00E + 20 1.00E + 22 46/0.5/10 46 2.00E + 21 1.00E + 22 46/2/10 71 1.00E + 20 1.00E + 22 71/0.1/10 71 3.00E + 20 1.00E + 22 71/0.3/10 90 5.00E + 20 3.00E + 21 *90/VM2/VM2 96 5.00E + 20 1.00E + 22 96/0.5/10 120 3.00E + 19 1.00E + 22 120/0.03/10 120 1.00E + 20 1.00E + 22 120/0.1/10 120 3.00E + 20 1.00E + 21 120/0.3/1 120 3.00E + 20 1.00E + 22 120/0.3/10 120 5.00E + 20 1.00E + 21 †120/0.5/1 120 5.00E + 20 5.00E + 21 120/0.5/5 120 5.00E + 20 1.00E + 22 120/0.5/10 120 5.00E + 20 1.50E + 22 120/0.5/15 120 5.00E + 20 2.00E + 22 †120/0.5/20 120 5.00E + 20 1.00E + 23 †120/0.5/100 120 1.00E + 21 1.00E + 22 ‡120/1/10 120 2.00E + 21 1.00E + 22 120/2/10 LT (km) . UMV (Pa s) . LMV (Pa s) . Abbreviation . 46 1.00E + 19 1.00E + 22 46/0.01/10 46 5.00E + 19 1.00E + 21 46/0.05/1 46 1.00E + 20 1.00E + 22 46/0.1/10 46 3.00E + 20 1.00E + 22 46/0.3/10 46 5.00E + 20 1.00E + 21 46/0.5/1 46 5.00E + 20 1.00E + 22 46/0.5/10 46 2.00E + 21 1.00E + 22 46/2/10 71 1.00E + 20 1.00E + 22 71/0.1/10 71 3.00E + 20 1.00E + 22 71/0.3/10 90 5.00E + 20 3.00E + 21 *90/VM2/VM2 96 5.00E + 20 1.00E + 22 96/0.5/10 120 3.00E + 19 1.00E + 22 120/0.03/10 120 1.00E + 20 1.00E + 22 120/0.1/10 120 3.00E + 20 1.00E + 21 120/0.3/1 120 3.00E + 20 1.00E + 22 120/0.3/10 120 5.00E + 20 1.00E + 21 †120/0.5/1 120 5.00E + 20 5.00E + 21 120/0.5/5 120 5.00E + 20 1.00E + 22 120/0.5/10 120 5.00E + 20 1.50E + 22 120/0.5/15 120 5.00E + 20 2.00E + 22 †120/0.5/20 120 5.00E + 20 1.00E + 23 †120/0.5/100 120 1.00E + 21 1.00E + 22 ‡120/1/10 120 2.00E + 21 1.00E + 22 120/2/10 Open in new tab The model is implemented in a spherical harmonic formalism where the load changes (ice thickness and ocean depth) are applied as an upper surface boundary condition. The modelling code has more typically been truncated to spherical harmonic degree 256 for global to continental-scale predictions (e.g. Whitehouse et al.2012a). Here we extend the maximum truncation to degree 512, equivalent to a 21 arc minute (∼35 km at equator) spatial resolution. This resolution increase allows the load changes within the narrow geography of Palmer Land to be better resolved, particularly for the higher resolution local ice load histories we implement here (described below). Resolution is particularly important when exploring relatively low viscosity Earth models where spatial gradients in deformation can be large. The extension of spherical harmonic truncation to higher degrees (512 and above) has previously been shown capable of capturing fine-grained isostatic motions for similar scenarios (Ivins et al.2011; Nield et al.2014; Wolstencroft et al.2014). The GPS sites with the closest spacing (LNTK-SUGG and BEAN-SUGG) are 42 and 60 arcmin apart, respectively, the other sites (including LNTK and BEAN) are all >60 arcmin apart. Spatial model resolution should therefore be sufficient to allow comparison between model output and data for the majority of sites. However, gradients in uplift rate between LNTK-SUGG and BEAN-SUGG sites may not be fully resolvable by the GIA modelling. The ocean loading is generated from the ice history by solving an advanced formulation of the sea level equation (e.g. Mitrovica & Milne 2003). The effects of the ocean-ice sheet mass exchange on Earth's moment of inertia and the migration of coastlines are also modelled (Milne & Mitrovica 1998). To generate uplift rate predictions comparable to the processed GPS data, modelled deformations were averaged over a period ±6 months from the present, defined as 2010 AD. As we correct the GPS data for elastic deformation (Section 2.2), no ice load changes were applied after 2000 AD, preventing any ‘double accounting’ of elastic effects for those ice models with load changes up to the present-day. This approach will also exclude any viscous deformation effects of load change since 2000 AD; however, as the minimum Maxwell time in this study is around 100 years, the impact is negligible. The additional rotational feedback on the linear-momentum balance (Martinec & Hagedoorn 2014) was not modelled. The degree-2 effect from this mechanism is predicted to be essentially uniform over Palmer Land; hence, in this study the signal is inseparable from LMV and reference frame uncertainty and will be filtered out by the offset correction applied (see Section 2.5). As the impact on present-day uplift rates is currently uncertain, we do not consider it further. Our GIA predictions are made in the CE (centre-of-solid-Earth excluding atmosphere and ocean) reference frame. 2.4.1 Ice history reconstructions Ice history models used in this study are the product of field data constraints combined with glaciological insights and, in some cases, physical ice sheet models of varying complexity. The main aspect of a given ice history, which effects GIA, is the timing and magnitude of ice mass change. These drivers interact with the rheology of the modelled mantle to produce the modelled GIA rate. A total of 9 different ice history reconstructions were used in this study. With the exception of ICE-5G (Peltier 2004), all were modifications on the W12 ice model (Whitehouse et al.2012a), where changes were applied in order to determine whether the GPS data are able to constrain Holocene regional-scale ice load changes across Palmer Land. Details of the origin and construction of each ice history are summarized in Table 4 and described in detail below. Table 4. Ice models used in this study, all run for a full glacial cycle: 122 ka to present (2010 AD). Name . Description . ICE-5G v1.2 after Peltier (2004), expanded to degree 512 W12 Whitehouse et al. (2012b), expanded from degree 256 to degree 512 W12_N12 Modelled ice history since 1850 (N12, Nield et al.2012) was added to W12 W12_AP1Ka Maximum MWP-LIA-like scenario, added to the last 1350 years of W12 W12_AP1Kb Minimum MWP-LIA-like scenario; 20 per cent of AP1Ka W12_Nev No ice loss in the Evans Ice Stream region after 5 ka W12_05Ev 50 per cent of W12 ice loss in the Evans Ice Stream region after 5 ka W12_07Ev 70 per cent of W12 ice loss in the Evans Ice Stream region after 5 ka, consistent with Bentley et al. (2006) W12_b5Ev Post 5 ka W12 ice loss in Evans Ice Stream region instead removed 10–5 ka, no change after 5 ka Name . Description . ICE-5G v1.2 after Peltier (2004), expanded to degree 512 W12 Whitehouse et al. (2012b), expanded from degree 256 to degree 512 W12_N12 Modelled ice history since 1850 (N12, Nield et al.2012) was added to W12 W12_AP1Ka Maximum MWP-LIA-like scenario, added to the last 1350 years of W12 W12_AP1Kb Minimum MWP-LIA-like scenario; 20 per cent of AP1Ka W12_Nev No ice loss in the Evans Ice Stream region after 5 ka W12_05Ev 50 per cent of W12 ice loss in the Evans Ice Stream region after 5 ka W12_07Ev 70 per cent of W12 ice loss in the Evans Ice Stream region after 5 ka, consistent with Bentley et al. (2006) W12_b5Ev Post 5 ka W12 ice loss in Evans Ice Stream region instead removed 10–5 ka, no change after 5 ka Open in new tab Table 4. Ice models used in this study, all run for a full glacial cycle: 122 ka to present (2010 AD). Name . Description . ICE-5G v1.2 after Peltier (2004), expanded to degree 512 W12 Whitehouse et al. (2012b), expanded from degree 256 to degree 512 W12_N12 Modelled ice history since 1850 (N12, Nield et al.2012) was added to W12 W12_AP1Ka Maximum MWP-LIA-like scenario, added to the last 1350 years of W12 W12_AP1Kb Minimum MWP-LIA-like scenario; 20 per cent of AP1Ka W12_Nev No ice loss in the Evans Ice Stream region after 5 ka W12_05Ev 50 per cent of W12 ice loss in the Evans Ice Stream region after 5 ka W12_07Ev 70 per cent of W12 ice loss in the Evans Ice Stream region after 5 ka, consistent with Bentley et al. (2006) W12_b5Ev Post 5 ka W12 ice loss in Evans Ice Stream region instead removed 10–5 ka, no change after 5 ka Name . Description . ICE-5G v1.2 after Peltier (2004), expanded to degree 512 W12 Whitehouse et al. (2012b), expanded from degree 256 to degree 512 W12_N12 Modelled ice history since 1850 (N12, Nield et al.2012) was added to W12 W12_AP1Ka Maximum MWP-LIA-like scenario, added to the last 1350 years of W12 W12_AP1Kb Minimum MWP-LIA-like scenario; 20 per cent of AP1Ka W12_Nev No ice loss in the Evans Ice Stream region after 5 ka W12_05Ev 50 per cent of W12 ice loss in the Evans Ice Stream region after 5 ka W12_07Ev 70 per cent of W12 ice loss in the Evans Ice Stream region after 5 ka, consistent with Bentley et al. (2006) W12_b5Ev Post 5 ka W12 ice loss in Evans Ice Stream region instead removed 10–5 ka, no change after 5 ka Open in new tab ICE-5G and W12. Two ice loading histories from the literature, ICE-5G v1.2 (Peltier 2004) and W12 (Whitehouse et al.2012a), were expanded to spherical harmonic degree 512 for this study. Care was taken to preserve both total global ice volume and ice margin location during resampling. The global ice volume in the rescaled models at 20 ka (LGM in W12) is within 0.5 per cent of the original models. This resampling process does not; however, alter the original spatial accuracy of the ice models. ICE-5G has been superseded by ICE-6G_C (Argus et al.2014; Peltier et al.2015), but at the time of writing, the ICE-6G_C loading history has not been made available (although a comparison to published model output from ICE-6G is made below). Consideration of ICE-5G is still relevant as this model is used in GRACE ice mass balance studies (e.g. King et al.2012; Shepherd et al.2012). The W12a ice model (Whitehouse et al.2012a) is not used in this study as it includes a somewhat simplistic modification to W12, which was applied in order to address the apparent overprediction of previously published GPS rates in the AP when using an Earth model more suited to an Antarctic average. W12 with N12. The W12_N12 ice history is a modification to W12, which incorporates ice load changes in the AP over the last 150 years as detailed by Nield et al. (2012), referred to here as N12. W12_N12 allows us to test the GPS data for sensitivity to relatively recent load changes, although it is important to note that due to the upper mantle viscosities adopted within our modelling (described above), a large effect is not to be expected. The N12 model was produced from 5 km resolution ice sheet model output (details in: Nield et al.2012) and has a resolution defined by spherical harmonic degree 512, thus the accuracy matches that of the underlying model resolution and is higher than that of W12, even when the latter has been resampled as described above. The combined W12_N12 model (Fig. S2) was produced by adding ice load changes from N12 to the last 150 yr of W12 at 10-yr intervals. Present-day ice thickness was preserved by subtracting the additional load of N12 from W12 for all time intervals after 20 ka but prior to the addition of N12. AP1K. With N12, Nield et al. (2012) demonstrated how changes in temperature and precipitation can lead to changes in ice thickness on the AP in the last 150 years. Extending this concept further back into the Holocene, it is plausible that climate oscillations similar to the northern hemisphere Medieval Warm Period (MWP) and Little Ice Age (LIA) could result in snow/ice loading changes sufficient to influence local viscous deformation rates in the AP. We refer to this class of ice model as: ‘Antarctic Peninsula (last) 1 kyr’, abbreviated to AP1K. Evidence of MWP-like and LIA-like climate history in the AP is equivocal; there is support for a global long-term cooling trend over the last 2 kyr, which ended in the 19th century, but within this there is significant regional variation (e.g. Neukom et al.2014). Analysis of the available Antarctic record suggests a MWP-like period (+0.35°C relative to present) ended at +1288.0 with the onset of a LIA-like period (−2.0 °C relative to present), which persisted until +1808.0 (Bertler et al.2011). In West Antarctica, analysis of the WAIS-divide ice core suggests colder-than-average temperatures between +1300.0 and +1800.0, with warmer temperatures (+0.5 °C) before this time (Orsi et al.2012; Abram et al.2014). The James Ross Island ice core shows no MWP/LIA-like signal (Mulvaney et al.2012; Abram et al.2014). However, the glacier on James Ross Island, which is located 700 km north of Palmer Land, is somewhat unique and unlikely to fully represent the wider AP or Palmer Land (Sime et al.2009). The modelled scenarios described below assume that a MWP/LIA-like climatic oscillation occurred in Palmer Land and that ice thickness variations are precipitation-dominated (Thomas & Bracegirdle 2009). The ice history was generated using the Glimmer ice sheet model (Rutt et al.2009) in a similar manner to N12 (Nield et al.2012). The ice sheet model was initially run for 30 kyr to reach steady state using accumulation patterns for +1850.0 with present day surface temperatures. We used this long spin-up time as we were only interested in changes in ice thickness caused by the imposed climatic changes. A MWP-like climate was then prescribed using accumulation rates from +2000.0 and temperatures +0.5 °C relative to present (Orsi et al.2012). These conditions were maintained for 600 yr, after which a LIA-like transition was simulated by prescribing +1850.0 accumulation rates and temperatures −2.0 °C relative to present (Bertler et al.2011). These conditions were also maintained for 600 yr. In this scenario, the warm period starts in +650.0, ends in +1250.0 and is followed by the cool period which ends in +1850.0. Ice thicknesses were output at 50 yr intervals (see Fig. S3 for more details). The accumulation history used to force the ice-sheet model has a strong dependence on the Gomez ice core (Nield et al.2012), which shows generally larger accumulation rates than other nearby ice cores (Thomas et al.2008). Therefore, we assume that the AP1K modelled changes in ice load are maximum plausible values and call this model AP1Ka. Ice thickness changes in AP1Ka were then uniformly reduced by 80 per cent and called AP1Kb, our assumed minimum model for the scenario. AP1Ka and AP1Kb were each added to the last 1360 years of W12 in the same manner as N12. As with W12_N12, the AP1K portion of W12_AP1K* is produced from a 5 km Glimmer model grid but has a resolution that reflects its expansion to maximum spherical harmonic degree 512. We accept that this approach may not reflect actual climate variations in the AP over the last ∼1 kyr. Other approaches are possible; Ivins et al. (2000) adopted the assumption that ice thickness change is temperature-dominated. As a result they modelled an oscillatory load change in the northern and central AP which has the opposite sense to our AP1K scenario. There may also be further complexities in the spatial and temporal distribution of precipitation, accumulation and ice sheet flow, which are yet to be constrained. Evans Ice Stream. Recent ice loading changes within the Evans Ice Stream region (Fig. 1) could influence modelled uplift rates in Palmer Land as it lies close to several GPS sites. There is a lack of information on the timing of grounding line retreat in the western Weddell Sea, which buttresses this ice stream. As a result, the Evans Ice Stream is relatively unconstrained in W12 (Whitehouse et al.2012b). The single nunatak (rock outcrop) constraint for this region lies in the nearby Behrendt Mountains (Fig. 1) and indicates <300 m of ice thinning since 7.2 ka (Bentley et al.2006). Within W12, there is ∼450 m of ice thinning in this region between 7 and 2 ka, suggesting that W12 overpredicts unloading by at least 30 per cent (Whitehouse et al.2012b). To constrain the possible impact of this overprediction, we produced a series of scenarios where the magnitude of ice load change in the Evans region between 5 and 2 ka was reduced (Table 4; Fig. S4). These models are named using the format W12_*Ev, where * describes the modification applied (individual modifications described in Table 4). Present-day ice thickness was maintained for W12_NEv, W12_05Ev and W12_07Ev by subtracting the total magnitude of the ice thickness modifications between 5 ka and the present from the W12 ice history between 20 ka (LGM) and 5 ka. 2.5 Model fit analysis approach To compare the model predictions of vertical uplift rates with the vertical GPS rates we adopted a weighted root mean square (WRMS) fit analysis after the weighted mean difference (termed offset) between data and predictions for a given model has been subtracted from the model predictions. A different uniform offset is applied to each model; the magnitude of this offset may be attributed to uncertainty in the LMV, uncertainty in the far-field elastic signal or possible discrepancies between the reference frames of the GPS data and the GIA model (e.g. Argus 2012). Our offset correction does not affect the regional-scale spatial patterns of GIA, which we are investigating here. Limits to the magnitude of the offset that can reasonably be applied are presented in Section 4.3. The fit to the GPS data for each model prediction was determined by first calculating the data-model prediction offset \begin{equation} D = \frac{{\sum\nolimits_1^n {w_i (o_i - p_i )} }}{{\sum\nolimits_1^n {w_i } }}, \end{equation} (1) where the weight |$w_i = \frac{1}{{\sigma _i^2 }}$|⁠, oi is the observation, σi is the observation error, pi is the prediction and n is the number of GPS sites being compared. The offset-corrected WRMS is \begin{equation} DWRMS = \sqrt {\frac{{\sum\nolimits_1^n {w_i (o_i - (p_i + D))^2 } }}{{\sum\nolimits_1^n {w_i } }}} , \end{equation} (2)D is negative when the mean of the predictions is greater than the mean of the observations and positive when it is smaller. Addition of D aligns the means of the data and predictions by shifting the predictions downwards if they were originally too high or upwards if they were originally too low. Applying this offset makes our data-model comparisons insensitive to the chosen LMV. 3 GPS RESULTS The rates of surface deformation measured at the 12 Palmer Land GPS sites (corrected for elastic deformation and plate rotation) are presented in Tables 1, 2, Figs 2 and 3. To first order, vertical rates show a reducing trend from south to north (Fig. 2B). On closer inspection, the HTON, MKIB and BEAN sites (the ‘SE sites’) located on the southeast coast of Palmer Land (red squares in Fig. 2b) clearly yield significantly higher rates than the rest of the GPS sites. The elevated rates of the SE sites are spatially coherent and suggest a source of uplift in the region of the southwestern Ronne Ice Shelf. The agreement between HTON, MKIB and BEAN and their varying distances from local glaciers, makes it unlikely that the high rates observed are caused by unmodelled site-specific elastic effects. Rates at SUGG and HAAG (Argus et al.2014; in grey in Fig. 2b) have relatively larger error but are consistent with the observation of large uplift in this region. Figure 2. Open in new tabDownload slide Vertical GPS rates after elastic deformation correction. (a) Map view. The symbol size is scaled by the inverse of the 1σ uncertainties of the vertical GPS rates. (b) Vertical rates by latitude with red squares highlighting the SE sites (see text for detail). The empirically-derived uplift rates of Gunter et al. (2014; solution CSR RL05 DDK5, corrected for 1.9 mm yr−1 GIA bias) are also shown, extracted along the blue transect line in panel (a). Figure 2. Open in new tabDownload slide Vertical GPS rates after elastic deformation correction. (a) Map view. The symbol size is scaled by the inverse of the 1σ uncertainties of the vertical GPS rates. (b) Vertical rates by latitude with red squares highlighting the SE sites (see text for detail). The empirically-derived uplift rates of Gunter et al. (2014; solution CSR RL05 DDK5, corrected for 1.9 mm yr−1 GIA bias) are also shown, extracted along the blue transect line in panel (a). Figure 3. Open in new tabDownload slide Plate rotation corrected horizontal rates and elastic correction. The horizontal elastic correction is shown in blue and is of negligible magnitude. The brown and black arrows represent GPS rates first corrected for elastic effects and then rotated using either the Argus et al. (2014) or ITRF2008 (Altamimi et al.2011) plate rotation models to illustrate the significant differences. Ellipses show 1σ errors of the GPS rates; the large errors for SUGG and HAAG adopted from Argus et al. (2014) are omitted for clarity. Figure 3. Open in new tabDownload slide Plate rotation corrected horizontal rates and elastic correction. The horizontal elastic correction is shown in blue and is of negligible magnitude. The brown and black arrows represent GPS rates first corrected for elastic effects and then rotated using either the Argus et al. (2014) or ITRF2008 (Altamimi et al.2011) plate rotation models to illustrate the significant differences. Ellipses show 1σ errors of the GPS rates; the large errors for SUGG and HAAG adopted from Argus et al. (2014) are omitted for clarity. In Fig. 3, we show horizontal GPS velocities after subtracting the elastic model and the two plate rotation models. Differences between the two plate rotation models are much larger than the horizontal elastic correction. Noting that horizontal vectors will tend to point away from regions of past ice mass loss, three centres of deformation are indicated by the horizontal velocities: the southern Ronne Ice Shelf (consistent with the vertical rates of the SE sites in Fig. 2), offshore western Palmer Land, and Marguerite Bay. This result is independent of the choice of plate rotation model. The horizontal velocities at TRVE (Fig. 3) are interesting in that the application of the Altamimi et al. (2012) plate rotation model produces a vector that points away from the centre of Marguerite Bay, whereas the application of the Argus et al. (2014) model produces a vector pointing directly away from the present-day grounding line. The former would suggest a large centre of uplift within Marguerite Bay, whereas the latter would suggest that more recent mass loss, associated with the 1970s collapse of Wordie Ice Shelf (Doake & Vaughan 1991), dominates the deformation pattern at TRVE. We are not able to separate these two competing hypotheses using the horizontal velocities. The interpretation of horizontal GPS measurements may also be complicated by lateral viscosity variations across the wider region, which have the potential to perturb predictions of horizontal surface deformation in comparison to those calculated assuming a purely radially varying Earth model (Kaufmann et al.2005). Indeed, previous work in the northern AP identified a low viscosity mantle in that region (Nield et al.2014), suggesting that lateral variations may indeed be significant along the length of the AP. We focus solely on the vertical velocities for the remainder of the paper. 4 GIA MODELLING OUTPUT The GIA modelling comprised a wide range of individual Earth and ice model combinations and focussed on vertical deformation rates. Predictions were output at the GPS site locations and along the AP transect indicated in Fig. 2(a). The spatial deformation field is presented at GIA model resolution for map plots. 4.1 Comparison to published models Fig. 4 presents the elastic-corrected GPS rates plotted over modelled uplift rates from: ICE-5G, with its preferred Earth model (VM2), ICE-6G_C VM5a, W12 and two IJ05_R2 uplift predictions provided by E. Ivins (ICE-5G and W12 outputs were produced with the GIA code described in this study). It is clear that none of the models presented in Fig. 4 can fully satisfy the spatial pattern of the measured uplift rates, suggesting that the misfit cannot be explained by a simple vertical offset of these published models. For this reason the misfit must lie in the detailed Earth model, the ice model, or some combination of both. We address this question in the following two sections. Figure 4. Open in new tabDownload slide Comparison of predicted uplift rates from published models (background colour) with GPS data (symbols scaled by inverse of 1σ error): (a) ICE-5G v1.2 with 90/VM2/VM2 (Peltier 2004); (b) ICE-6G_C VM5a (Argus et al.2014; Peltier et al.2015); (c, d) IJ05_R2 with 115/0.2/4 and 65/0.2/1.5, respectively (Ivins et al.2013); (e) W12 with 120/1/10 (Whitehouse et al.2012a). Figure 4. Open in new tabDownload slide Comparison of predicted uplift rates from published models (background colour) with GPS data (symbols scaled by inverse of 1σ error): (a) ICE-5G v1.2 with 90/VM2/VM2 (Peltier 2004); (b) ICE-6G_C VM5a (Argus et al.2014; Peltier et al.2015); (c, d) IJ05_R2 with 115/0.2/4 and 65/0.2/1.5, respectively (Ivins et al.2013); (e) W12 with 120/1/10 (Whitehouse et al.2012a). 4.2 Sensitivity of uplift rate predictions to Earth model parameters Using W12 as a reference, Fig. 5 illustrates the impact of the uncertainty in Earth model parameters on possible deformation rates and spatial patterns in Palmer Land. Deformation rates in Fig. 5 are plotted along the transect shown in Fig. 2(a). The choice of LT and UMV influence the shorter wavelength spatial pattern. The lithosphere acts as a low pass filter, with greater LT leading to a longer wavelength pattern when all other parameters are held constant (Fig. 5a). The modelled UMV interacts strongly with the timing of the ice model load changes. Adopting a low UMV results in low present-day uplift rates due to the short relaxation time; moderate UMVs produce the maximum peak uplift signal; higher UMVs respond more slowly to load changes, producing lower peak rates at present that preserve more of the long term regional uplift signal (−69°S to −74°S Fig. 5b). Variations in the LMV have a very minor impact on the pattern of spatial deformation in a small region like Palmer Land but produce an essentially uniform vertical shift in the overall magnitude of deformation rates (Peltier 2004; Fig. 5c). Figure 5. Open in new tabDownload slide Sensitivity of predicted uplift rates to different rheology parameters. Predictions are given for points along the transect shown in Fig. 2(a) using the W12 ice model. (a–c) Model sensitivity to the indicated rheology parameter. Earth model name format: LT (km)/UMV (1021 Pa s)/LMV (1021 Pa s). Figure 5. Open in new tabDownload slide Sensitivity of predicted uplift rates to different rheology parameters. Predictions are given for points along the transect shown in Fig. 2(a) using the W12 ice model. (a–c) Model sensitivity to the indicated rheology parameter. Earth model name format: LT (km)/UMV (1021 Pa s)/LMV (1021 Pa s). 4.3 Acceptable range of D When comparing the GPS data and GIA model predictions we corrected for the weighted mean difference between data and model prediction D (eq. 1). Several factors contribute to the range of acceptable D values. Relative motions between the reference frames of the GPS data and GIA modelling could contribute an offset, primarily along the Z-axis of the reference frame. Recent estimates of the motion of CM (GPS) relative to CE (GIA model) suggest the two differ by less than ±0.5 mm yr−1 (Argus et al.2014). Far-field elastic effects, due predominantly to ongoing changes in global water distribution, could also influence deformation in the AP. As outlined in Section 2.2, we estimate that the combined effect of these two processes is uniformly ≈–0.5 mm yr−1 in Palmer Land over the GPS data collection period (Fig. S5). Uncertainty as to the most appropriate modelled value for LMV can introduce significant variation to our modelled uplift rates. Starting with a reference value of 1022 Pa s (Whitehouse et al.2012a), we modelled 2 orders of magnitude uncertainty in LMV (Fig. 5c). The model output indicates that only negative offsets are acceptable in Palmer Land. This asymmetric distribution is a product of the interaction of the LMV with the applied ice history. Only the longest time and length-scale variations in global ice and ocean loading excite deformation in the lower mantle. If the LMV is low, this signal decays more rapidly; if LMV is high, there is higher resistance to deformation and the overall signal is smaller. The result is a ‘sweet spot’ value of LMV, which produces a maximum lower mantle response (in terms of uplift rates) to the applied ice loading history. In Palmer Land this value is close to our reference value of 1022 Pa s. We obtained a range of 0 to −3 mm yr−1 for the LMV effect by taking the average difference between uplift rates along our AP transect for the 120/0.5/1 and 120/0.5/10 (format: LT km/UMV 1021 Pa s/LMV 1021 Pa s) Earth models from −77°S to −69°S (Fig. 5c). By combining the potential offset due to far field/reference frame effects with the LMV effect range we defined the reasonable range for D as 0 to −3.5 mm yr−1. A given model combination is considered unreasonable if D, rounded to 1 decimal place, falls outside this range. The GIA bias corrections applied by Gunter et al. (2014) to their empirical solutions is of a similar magnitude and direction (1.6–2.1 mm yr−1 subtracted from their empirical GIA estimate depending on the GRACE solution), but our method accounts for different sources of offset and is not comparable to the GIA bias correction of Gunter et al. (2014). 4.4 Sensitivity of uplift rate predictions to ice model variations Varying the 1-D Earth model alters the magnitude and time-sensitivity of the solid Earth response to the imposed load, varying ice load history impacts the prediction in a different manner. Fig. 6 shows the range of uplift rates produced when different ice models are combined with a fixed Earth model. ICE-5G predicts higher uplift than all W12 cases in northern Palmer Land and rates similar to W12_NEv (no recent ice loss in the Evans region over the last 5 ka) between −75°S and −76°S. W12 produces the greatest peak uplift rates while W12_N12 and W12_AP1Ka/b have a lower and flatter profile because the modelled late Holocene load increases counter the longer term signal. The modifications to the W12 Evans Ice Stream retreat influence the peak uplift, producing a ∼2 mm yr−1 variation in local uplift rates. These results suggest that local or sharp differences between data and model prediction are more likely to represent errors in the assumed ice load history than uncertainties in the Earth model. Figure 6. Open in new tabDownload slide Predicted uplift rates for the ice models used in this study: ICE-5G and a range of perturbations to the W12 ice model. All calculations use the same Earth model: 120/0.3/10. Values plotted along the AP transect from Fig. 2(a). Figure 6. Open in new tabDownload slide Predicted uplift rates for the ice models used in this study: ICE-5G and a range of perturbations to the W12 ice model. All calculations use the same Earth model: 120/0.3/10. Values plotted along the AP transect from Fig. 2(a). 4.5 Model fit to the data In this section, we move from comparing regional trends along the AP transect to analysing the fit of the GIA rate predictions to the GPS rates. We calculated D for each Earth and ice model combination. Models with a D value outside the acceptable range were rejected. The remaining model combinations were then ranked by lowest DWRMS (Table 5; Table S1). Table 5. The 10 highest-ranking (by DWRMS) model combinations for each of the three analysis scenarios described in the main text. Model combinations which have an offset (D) value outside the acceptable range (0 to −3.5 mm yr−1) are shaded grey. *DWRMS is not directly comparable between scenarios as the fit is calculated only with respect to the GPS sites under consideration. A full list of data-model comparisons can be found in Table S1. All sites . Ice . Earth . D . *DWRMS . ICE-5G 120/0.3/10 − 1.5 1.28 W12_b5Ev 120/0.3/10 − 0.6 1.38 W12_N12 120/0.3/10 − 0.5 1.40 W12_05Ev 120/0.3/10 − 0.7 1.42 W12_07Ev 120/0.3/10 − 0.8 1.43 W12_NEv 120/0.3/10 − 0.4 1.45 W12_AP1Ka 120/0.1/10 − 0.5 1.47 W12 120/0.3/10 − 1.0 1.49 W12_AP1Kb 120/0.3/10 − 1.0 1.49 W12_AP1Ka 120/0.3/10 − 0.9 1.52 HTON, MKIB and BEAN only Ice Earth D *DWRMS ICE-5G 120/0.5/10 − 0.1 0.55 W12_07Ev 90/VM2/VM2 4.3 0.04 W12_05Ev 71/0.1/10 3.4 0.05 W12_05Ev 71/0.3/10 3.1 0.05 W12_AP1Ka 90/VM2/VM2 4.3 0.07 W12_b5Ev 71/0.3/10 3.1 0.09 W12_05Ev 46/2/10 3.3 0.10 W12_07Ev 46/2/10 3.3 0.10 W12_AP1Ka 46/2/10 3.2 0.10 W12_b5Ev 46/0.5/10 4.7 0.10 All sites excluding HTON, MKIB and BEAN Ice Earth D *DWRMS ICE-5G 120/0.1/10 − 0.5 0.45 W12 120/0.1/10 − 0.1 0.51 W12_AP1Kb 120/0.1/10 − 0.4 0.54 W12_07Ev 120/0.1/10 0.0 0.56 ICE-5G 120/0.3/10 − 2.3 0.57 W12_b5Ev 120/0.3/10 − 1.5 0.58 W12_05Ev 120/0.3/10 − 1.6 0.63 ICE-5G 71/0.1/10 − 0.3 0.65 W12_NEv 120/0.3/10 − 1.2 0.67 W12_07Ev 120/0.3/10 − 1.6 0.70 All sites . Ice . Earth . D . *DWRMS . ICE-5G 120/0.3/10 − 1.5 1.28 W12_b5Ev 120/0.3/10 − 0.6 1.38 W12_N12 120/0.3/10 − 0.5 1.40 W12_05Ev 120/0.3/10 − 0.7 1.42 W12_07Ev 120/0.3/10 − 0.8 1.43 W12_NEv 120/0.3/10 − 0.4 1.45 W12_AP1Ka 120/0.1/10 − 0.5 1.47 W12 120/0.3/10 − 1.0 1.49 W12_AP1Kb 120/0.3/10 − 1.0 1.49 W12_AP1Ka 120/0.3/10 − 0.9 1.52 HTON, MKIB and BEAN only Ice Earth D *DWRMS ICE-5G 120/0.5/10 − 0.1 0.55 W12_07Ev 90/VM2/VM2 4.3 0.04 W12_05Ev 71/0.1/10 3.4 0.05 W12_05Ev 71/0.3/10 3.1 0.05 W12_AP1Ka 90/VM2/VM2 4.3 0.07 W12_b5Ev 71/0.3/10 3.1 0.09 W12_05Ev 46/2/10 3.3 0.10 W12_07Ev 46/2/10 3.3 0.10 W12_AP1Ka 46/2/10 3.2 0.10 W12_b5Ev 46/0.5/10 4.7 0.10 All sites excluding HTON, MKIB and BEAN Ice Earth D *DWRMS ICE-5G 120/0.1/10 − 0.5 0.45 W12 120/0.1/10 − 0.1 0.51 W12_AP1Kb 120/0.1/10 − 0.4 0.54 W12_07Ev 120/0.1/10 0.0 0.56 ICE-5G 120/0.3/10 − 2.3 0.57 W12_b5Ev 120/0.3/10 − 1.5 0.58 W12_05Ev 120/0.3/10 − 1.6 0.63 ICE-5G 71/0.1/10 − 0.3 0.65 W12_NEv 120/0.3/10 − 1.2 0.67 W12_07Ev 120/0.3/10 − 1.6 0.70 Open in new tab Table 5. The 10 highest-ranking (by DWRMS) model combinations for each of the three analysis scenarios described in the main text. Model combinations which have an offset (D) value outside the acceptable range (0 to −3.5 mm yr−1) are shaded grey. *DWRMS is not directly comparable between scenarios as the fit is calculated only with respect to the GPS sites under consideration. A full list of data-model comparisons can be found in Table S1. All sites . Ice . Earth . D . *DWRMS . ICE-5G 120/0.3/10 − 1.5 1.28 W12_b5Ev 120/0.3/10 − 0.6 1.38 W12_N12 120/0.3/10 − 0.5 1.40 W12_05Ev 120/0.3/10 − 0.7 1.42 W12_07Ev 120/0.3/10 − 0.8 1.43 W12_NEv 120/0.3/10 − 0.4 1.45 W12_AP1Ka 120/0.1/10 − 0.5 1.47 W12 120/0.3/10 − 1.0 1.49 W12_AP1Kb 120/0.3/10 − 1.0 1.49 W12_AP1Ka 120/0.3/10 − 0.9 1.52 HTON, MKIB and BEAN only Ice Earth D *DWRMS ICE-5G 120/0.5/10 − 0.1 0.55 W12_07Ev 90/VM2/VM2 4.3 0.04 W12_05Ev 71/0.1/10 3.4 0.05 W12_05Ev 71/0.3/10 3.1 0.05 W12_AP1Ka 90/VM2/VM2 4.3 0.07 W12_b5Ev 71/0.3/10 3.1 0.09 W12_05Ev 46/2/10 3.3 0.10 W12_07Ev 46/2/10 3.3 0.10 W12_AP1Ka 46/2/10 3.2 0.10 W12_b5Ev 46/0.5/10 4.7 0.10 All sites excluding HTON, MKIB and BEAN Ice Earth D *DWRMS ICE-5G 120/0.1/10 − 0.5 0.45 W12 120/0.1/10 − 0.1 0.51 W12_AP1Kb 120/0.1/10 − 0.4 0.54 W12_07Ev 120/0.1/10 0.0 0.56 ICE-5G 120/0.3/10 − 2.3 0.57 W12_b5Ev 120/0.3/10 − 1.5 0.58 W12_05Ev 120/0.3/10 − 1.6 0.63 ICE-5G 71/0.1/10 − 0.3 0.65 W12_NEv 120/0.3/10 − 1.2 0.67 W12_07Ev 120/0.3/10 − 1.6 0.70 All sites . Ice . Earth . D . *DWRMS . ICE-5G 120/0.3/10 − 1.5 1.28 W12_b5Ev 120/0.3/10 − 0.6 1.38 W12_N12 120/0.3/10 − 0.5 1.40 W12_05Ev 120/0.3/10 − 0.7 1.42 W12_07Ev 120/0.3/10 − 0.8 1.43 W12_NEv 120/0.3/10 − 0.4 1.45 W12_AP1Ka 120/0.1/10 − 0.5 1.47 W12 120/0.3/10 − 1.0 1.49 W12_AP1Kb 120/0.3/10 − 1.0 1.49 W12_AP1Ka 120/0.3/10 − 0.9 1.52 HTON, MKIB and BEAN only Ice Earth D *DWRMS ICE-5G 120/0.5/10 − 0.1 0.55 W12_07Ev 90/VM2/VM2 4.3 0.04 W12_05Ev 71/0.1/10 3.4 0.05 W12_05Ev 71/0.3/10 3.1 0.05 W12_AP1Ka 90/VM2/VM2 4.3 0.07 W12_b5Ev 71/0.3/10 3.1 0.09 W12_05Ev 46/2/10 3.3 0.10 W12_07Ev 46/2/10 3.3 0.10 W12_AP1Ka 46/2/10 3.2 0.10 W12_b5Ev 46/0.5/10 4.7 0.10 All sites excluding HTON, MKIB and BEAN Ice Earth D *DWRMS ICE-5G 120/0.1/10 − 0.5 0.45 W12 120/0.1/10 − 0.1 0.51 W12_AP1Kb 120/0.1/10 − 0.4 0.54 W12_07Ev 120/0.1/10 0.0 0.56 ICE-5G 120/0.3/10 − 2.3 0.57 W12_b5Ev 120/0.3/10 − 1.5 0.58 W12_05Ev 120/0.3/10 − 1.6 0.63 ICE-5G 71/0.1/10 − 0.3 0.65 W12_NEv 120/0.3/10 − 1.2 0.67 W12_07Ev 120/0.3/10 − 1.6 0.70 Open in new tab When all 12 GPS sites are considered, model combinations with a reasonable D value display poor DWRMS fits (Table 5, Figs 7a and b). To investigate the source of the overall misfit, we analysed the SE sites (HTON, MKIB and BEAN) and other GPS sites separately. We divided the sites in this manner because the SE sites display a spatially coherent pattern of much higher uplift rates than the other sites (Fig. 2). Figs 7(c)–(f) show the best-fitting model combinations for just the SE sites and for the remaining sites. Figure 7. Open in new tabDownload slide (a and b) Best model fit to all GPS sites; ICE-5G 120/0.3/10. (c and d) Best model fit to just the SE sites (highlighted in red in b, d and f); ICE-5G 120/0.5/10. (e and f) Best model fit to all sites except the SE sites; ICE-5G 120/0.1/10. The underprediction at the SE sites by this model is given in brackets above site names in (f) in mm yr−1. Figure 7. Open in new tabDownload slide (a and b) Best model fit to all GPS sites; ICE-5G 120/0.3/10. (c and d) Best model fit to just the SE sites (highlighted in red in b, d and f); ICE-5G 120/0.5/10. (e and f) Best model fit to all sites except the SE sites; ICE-5G 120/0.1/10. The underprediction at the SE sites by this model is given in brackets above site names in (f) in mm yr−1. Taken alone, rates at the SE sites can be fit (Fig. 7d) but only by a model that does not fit the other sites. When we excluded the SE sites from the fit analysis (Figs 7e and f), a good fit to the remaining sites is achieved with a reasonable offset. The best fit model results in an underprediction at HTON, MKIB and BEAN of ∼2 mm yr−1 (indicated in Fig. 7f). 4.6 Comparison between GPS uplift rates and independent uplift estimates The GPS uplift rates are compared with the empirical estimates of GIA-related uplift produced by Gunter et al. (2014) in Fig. 8; the spatial pattern and magnitude of our GPS and the empirical rates are in good agreement. The empirical rates from Gunter et al. (2014) used for comparison were those derived from the variant using the CSR RL05 DDK GRACE solutions. Figure 8. Open in new tabDownload slide (a) Comparison of the elastic-corrected GPS rates (symbols) with the Gunter et al. (2014) empirical uplift rates (background); (b) site comparison with D offset correction and fit values shown. Figure 8. Open in new tabDownload slide (a) Comparison of the elastic-corrected GPS rates (symbols) with the Gunter et al. (2014) empirical uplift rates (background); (b) site comparison with D offset correction and fit values shown. The spatial detail of the Gunter et al. (2014) model is limited by the spherical harmonic degree 60 resolution of the GRACE data and the 400-km Gaussian smoothing during processing. We note that the method used by Gunter et al. (2014) relies on coverage by ICESat, GRACE and surface mass balance and firn densification information derived from the RACMO2 regional atmospheric climate model (Ligtenberg et al.2011; Lenaerts et al.2012). In particular, since ICESat data coverage is reduced in the northern Peninsula the accuracy of their estimate may be affected in the northern part of our region (e.g. near TRVE). The estimates of Gunter et al. (2014) and those presented here are also based on data covering two different periods and hence any errors in removing short-term effects will alter this comparison. Nevertheless, the agreement is impressive and supports the accuracy of our GPS velocities and the Gunter et al. (2014) results and indirectly corroborates the ice mass balance estimate of this region from Gunter et al. (2014). The original Gunter et al. (2014) empirical uplift rates lie below the rates we measure in Palmer Land (Fig. 2b). This offset is probably due to the GIA bias correction that Gunter et al. (2014) applied to their empirical uplift rates. These rates were corrected by assuming that a low precipitation zone in East Antarctica (LPZ; their term) is free of GIA. Following this assumption they calibrated the whole of Antarctic GIA to zero in the LPZ. As Gunter et al. (2014) acknowledge, this approach implicitly removes any far-field GIA signal and assumes that such a signal is uniform across Antarctica. Therefore, the difference between the Gunter et al. (2014) rates corrected for ‘GIA bias’ and our GPS rates may be due to differences in this far field GIA signal in the AP relative to their LPZ. The D value required to align the means of the Gunter et al. (2014) rates and the GPS rates (Fig. 8) is 1.3 mm yr−1 upwards, effectively undoing a significant proportion of the 1.9 mm yr−1 downward GIA bias correction applied to the original solution by Gunter et al. (2014). Once this correction is applied, the Gunter et al. (2014) empirical rates have a better DWRMS fit than the best fitting GIA model to all sites (0.74 mm yr−1 versus 1.28 mm yr−1). Our observations also agree with the spatial pattern of the GIA estimate produced by Sasgen et al. (2013; fig. 4a), although the magnitude is different. The origin of this difference is not clear but could be due to the limited GPS data available for Sasgen et al. (2013) to constrain their correction in the AP sector or the smoothing of regional offsets due to lateral variations in rheology. 5 DISCUSSION The results presented above comprise the most detailed study of GIA in Palmer Land to date. We present a new set of GPS uplift rates and have carried out high resolution GIA modelling of the region using a broad range of rheological parameters. The use of a bounded offset correction (D) in the fit analysis between data and predictions eliminates long wavelength signals related to uncertainty in LMV, far-field elastic signals and possible reference frame offsets between GPS data and model predictions. The failure of any of our models to closely replicate the GPS-observed uplift pattern suggests that the ice histories adopted here are the primary source of error in the modelling. The agreement in spatial pattern and magnitude with the Gunter et al. (2014) estimates of vertical uplift (after removing an offset) adds weight to the robustness of our new GPS-observed uplift estimates and the consequent need for improved GIA modelling, and particularly ice history, in this region. This strongly implies that there is an important data gap in the ice history reconstruction in or around Palmer Land, most probably in the region of the southwestern Ronne Ice Shelf. Whitehouse et al. (2012b) noted that there are currently few data that constrain our understanding of ice sheet behaviour in the Weddell Sea since the LGM. This has not improved significantly at the time of writing (see review by: Hillenbrand et al.2013) and ice models such as IJ05_R2 (Ivins et al.2013) and ICE-6G_C (Argus et al.2014) published since W12 have not drawn on additional ice margin constraints. Model output from ICE-6G_C (Argus et al.2014) shows a similar spatial pattern to the empirical estimates of Gunter et al. (2014) but has much higher uplift rates (∼13 mm yr−1) over the Weddell Sea than Gunter et al. (5–8 mm yr−1). Clearly this is an area where it is very difficult to constrain past ice loading and as such it represents one of the main areas of uncertainty in current GIA models. With the improved accuracy and coverage of GPS data, there is the risk of overfitting to uplift rates, which may contain poorly constrained elastic or recent viscous effects (e.g. Nield et al.2012, 2014). This issue is best avoided by ensuring that field evidence of ice thickness and extent is honoured by ice history reconstructions. Improved ice history constraints are therefore essential if we are to understand the GPS data in southeast Palmer Land. Specifically, more constraints are required on the timing of grounding line retreat and on former ice thickness, such as data from nunataks and ice cores (e.g. Bentley et al.2010). With these data it will be possible to determine whether grounding line retreat in the southwestern Weddell Sea was significantly later than assumed in models such as W12 or whether more ice than currently assumed was lost from that region during the late Holocene, or some combination of the two. Reduced uplift rate uncertainties at all sites and especially SUGG and HAAG, together with new GPS estimates from sites further east, will further assist in testing competing models. The relatively sharp east–west gradient in the GPS uplift rates suggests a low viscosity Earth rheology in this region. Such rheology is also favoured by our modelling (despite inaccuracies in our adopted ice histories), with the set of preferred Earth rheology models sharing a thick LT (120 km), a relatively low UMV (1–3 × 1020 Pa s) and an unconstrained LMV. The absence of significant uplift on the western side of Palmer Land, despite evidence of substantial post-LGM ice unloading in this region (Bentley et al.2006; Hodgson et al.2009; Johnson et al.2012) also suggests that the related signal has decayed or been supressed by recent accumulation increases (Nield et al.2012), as would occur with a low UMV. We note that very low viscosities are required to explain recent rapid uplift in the northern AP (Nield et al.2014). It is conceivable that east–west variations in LT or mantle viscosity could explain some of the higher rates at the SE sites. However, in such a geographically small region, the filtering effect of the thick LT favoured by the data-model fits (Table 5) makes it unlikely that the entire residual can be explained in this way. Ultimately, an understanding of the impact of such lateral variations of Earth rheology must await the application of the next generation of GIA modelling codes and the data required to constrain them (Latychev et al.2005; Whitehouse et al.2006; A et al.2013; Austermann et al.2013; van der Wal et al.2013, 2015). Thus the difference between the SE sites and the rest of the GPS sites may contain a component due to lateral rheology variations, but it is highly unlikely that the pattern of uplift rates can be fully explained without modifications to the assumed ice history. 6 CONCLUSION We present data from a spatially dense network of GPS sites in Palmer Land in the southern Antarctic Peninsula. After removing a model of elastic effects and tectonic plate motion, horizontal GPS velocities suggest centres of former ice loss in Marguerite Bay, west of Palmer Land and in the region of the southwest Ronne Ice Shelf. Differences between competing models of plate tectonic motion represent an important limitation to the further interpretation of the horizontal velocities. We compared our new GPS uplift data with the spatially continuous, but independent, estimates of Gunter et al. (2014) and show good agreement, supporting the empirical GIA methodology of Gunter et al. (2014). We attempted to reproduce the GPS uplift rates using a model of GIA considering a range of ice histories, including modifications to a published ice history (W12). Vertical GPS rates in the central and western Antarctic Peninsula can be fit with a range of ice histories by adjusting for the likely weaker rheological structure of the mantle beneath the AP (relative to Whitehouse et al.2012a). The relatively sharp east–west gradient of uplift in southern Palmer Land suggests a low viscosity upper mantle rheology, and hence the observed uplift rates are more sensitive to ice loading changes in the late Holocene than would be the case with a higher viscosity upper mantle. Uplift rates at our sites on the southeast coast of Palmer Land cannot be fit with existing ice models within a reasonable range of uncertainty. The failure to reproduce the observed spatial pattern of uplift suggests that existing Holocene deglaciation models for the southwestern Weddell Sea are wrong. Either grounding line retreat occurred significantly later in the Holocene than previously assumed or the region hosted more ice at the LGM. Given evidence for low upper mantle viscosity in this region, which would cause more recent ice-load changes to dominate the present-day uplift, we consider later retreat as being most likely. This outcome joins other recent observations (e.g. Bradley et al.2015) in highlighting the need for further glaciological, geological and geodetic work in the southern Weddell Sea to constrain the timing of grounding line retreat in the Holocene. The current GPS data set cannot reliably discriminate between different post-5 ka ice history scenarios. GPS uncertainties will reduce as time-series grow, especially at SUGG and HAAG and, together with unpublished GPS sites further east, will allow the spatial pattern of present-day uplift in this region to be well characterised. Most importantly, improved field constraints on the local late Holocene ice history are required to significantly enhance our understanding of present-day uplift rates across Palmer Land and the southern Weddell Sea regions. This work was funded by NERC grant: NE/F01452X/1. MAK is a recipient of an Australian Research Council Future Fellowship (project number FT110100207). PLW is a recipient of a NERC Independent Research Fellowship (NE/K009958/1). REMR is recipient of a Netherlands Organisation for Scientific Research VIDI grant (project number 864.12.012). OD is supported by the New Netherlands Polar Programme of the Netherlands Organization for Scientific Research. We thank Matt Burke for deploying the CAPGIA sites and British Antarctic Survey engineers for maintaining them and downloading data. We acknowledge the hard work of the ANET team in deploying and maintaining SUGG and HAAG. This work is based on data services provided by the UNAVCO Facility with support from the National Science Foundation (NSF) and National Aeronautics and Space Administration (NASA) under NSF Cooperative Agreement No. EAR-0735156. We thank NASA JPL for making GIPSY available, Glenn Milne for providing the GIA code and Erik Ivins and Richard Peltier for making their model outputs available. Erik Ivins, Zdeněk Martinec and Editor Duncan Agnew provided comments that helped improve the clarity of the manuscript. REFERENCES A G. Wahr J. Zhong S. Computations of the viscoelastic response of a 3-D compressible Earth to surface loading: an application to Glacial Isostatic Adjustment in Antarctica and Canada Geophys. J. Int. 2013 192 557 572 Google Scholar Crossref Search ADS WorldCat Abram N.J. et al. Acceleration of snow melt in an Antarctic Peninsula ice core during the twentieth century Nat. Geosci. 2013 6 404 411 Google Scholar Crossref Search ADS WorldCat Abram N.J. Mulvaney R. Vimeux F. Phipps S.J. Turner J. England M.H. Evolution of the Southern Annular Mode during the past millennium Nat. Climate Change 2014 4 564 569 Google Scholar Crossref Search ADS WorldCat Adie R.J. The geology of Antarctica Antarctic Research: The Matthew Fontaine Maury Memorial Symposium 1962 American Geophysical Union 26 39 pp. Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Agnew D.C. The time-domain behavior of power-law noises Geophys. Res. Lett. 1992 19 333 336 Google Scholar Crossref Search ADS WorldCat Agnew D.C. NLOADF: a program for computing ocean-tide loading J. geophys. Res. 1997 102 5109 5110 Google Scholar Crossref Search ADS WorldCat Altamimi Z. Collilieux X. Metivier L. ITRF2008: an improved solution of the international terrestrial reference frame J. Geod. 2011 85 457 473 Google Scholar Crossref Search ADS WorldCat Altamimi Z. Metivier L. Collilieux X. ITRF2008 plate motion model J. geophys. Res. 2012 117 doi:10.1029/2011jb008930 Google Scholar OpenURL Placeholder Text WorldCat Argus D.F. Uncertainty in the velocity between the mass center and surface of Earth J. geophys. Res.: Solid Earth 2012 117 doi:10.1029/2012JB009196 Google Scholar OpenURL Placeholder Text WorldCat Argus D.F. Peltier W.R. Drummond R. Moore A.W. The Antarctica component of postglacial rebound model ICE-6G_C (VM5a) based on GPS positioning, exposure age dating of ice thicknesses, and relative sea level histories Geophys. J. Int. 2014 198 1 537 563 Google Scholar Crossref Search ADS WorldCat Austermann J. Mitrovica J.X. Latychev K. Milne G.A. Barbados-based estimate of ice volume at Last Glacial Maximum affected by subducted plate Nat. Geosci. 2013 6 553 557 Google Scholar Crossref Search ADS WorldCat Barker P.F. Scotia Sea regional tectonic evolution: implications for mantle flow and palaeocirculation Earth-Sci. Rev. 2001 55 1 39 Google Scholar Crossref Search ADS WorldCat Bentley C.R. Wahr J.M. Satellite gravity and the mass balance of the Antarctic ice sheet J. Glaciol. 1998 44 207 213 Google Scholar Crossref Search ADS WorldCat Bentley M.J. Hodgson D.A. Smith J.A. Cox N.J. Relative sea level curves for the South Shetlands and Marguerite Bay regions, Antarctic Peninsula Quatern. Sci. Rev. 2005 24 1203 1216 Google Scholar Crossref Search ADS WorldCat Bentley M.J. Fogwill C.J. Kubik P.W. Sugden D.E. Geomorphological evidence and cosmogenic 10Be/26Al exposure ages for the Last Glacial Maximum and deglaciation of the Antarctic Peninsula Ice Sheet Bull. geol. Soc. Am. 2006 118 9–10 1149 1159 Google Scholar Crossref Search ADS WorldCat Bentley M.J. Fogwill C.J. Le Brocq A.M. Hubbard A.L. Sugden D.E. Dunai T.J. Freeman S.P. Deglacial history of the West Antarctic Ice Sheet in the Weddell Sea embayment: constraints on past ice volume change Geology 2010 38 411 414 Google Scholar Crossref Search ADS WorldCat Bentley M.J. Johnson J.S. Hodgson D.A. Dunai T. Freeman S.P.H.T. Ó Cofaigh C. Rapid deglaciation of Marguerite Bay, western Antarctic Peninsula in the early Holocene Quatern. Sci. Rev. 2011 30 3338 3349 Google Scholar Crossref Search ADS WorldCat Bertiger W. Desai S. Haines B. Harvey N. Moore A. Owen S. Weiss J. Single receiver phase ambiguity resolution with GPS data J. Geod. 2010 84 327 337 Google Scholar Crossref Search ADS WorldCat Bertler N.A.N. Mayewski P.A. Carter L. Cold conditions in Antarctica during the Little Ice Age—Implications for abrupt climate change mechanisms Earth planet. Sci. Lett. 2011 308 1 41 51 Google Scholar Crossref Search ADS WorldCat Blewitt G. Self-consistency in reference frames, geocenter definition, and surface loading of the solid Earth J. geophys. Res. 2003 108 B2 2103 doi:10.1029/2002JB002082 Google Scholar Crossref Search ADS WorldCat Boehm J. Werl B. Schuh H. Troposphere mapping functions for GPS and very long baseline interferometry from European Centre for Medium-Range Weather Forecasts operational analysis data J. geophys. Res. 2006 111 B02406 doi:10.1029/2005JB003629 Google Scholar Crossref Search ADS WorldCat Bradley S.L. Hindmarsh R.C. Whitehouse P.L. Bentley M.J. King M.A. Low post-glacial rebound rates in the Weddell Sea due to Late Holocene ice-sheet readvance Earth planet. Sci. Lett. 2015 413 79 89 Google Scholar Crossref Search ADS WorldCat Cheng M. Ries J.C. Tapley B.D. Variations of the Earth's figure axis from satellite laser ranging and GRACE J. geophys. Res. 2011 116 B01409 doi:10.1029/2010JB000850 Google Scholar OpenURL Placeholder Text WorldCat Clarke P.J. Lavallee D.A. Blewitt G. van Dam T.M. Wahr J.M. Effect of gravitational consistency and mass conservation on seasonal surface mass loading models Geophys. Res. Lett. 2005 32 L08306 doi:10.1029/2005GL022441 Google Scholar OpenURL Placeholder Text WorldCat Dietrich R. et al. The SCAR 95 GPS campaign: objectives, data analysis and final solution The Geodetic Antarctic Project GAP95 - German Contributions to the SCAR 95 Epoch Campaign 1996 Deutsche Geodätische Kommission 9 14 pp. Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Dietrich R. et al. ITRF coordinates and plate velocities from repeated GPS campaigns in Antarctica—an analysis based on different individual solutions J. Geod. 2001 74 756 766 Google Scholar Crossref Search ADS WorldCat Dietrich R. Rulke A. Ihde J. Lindner K. Miller H. Niemeier W. Schenke H.W. Seeber G. Plate kinematics and deformation status of the Antarctic Peninsula based on GPS Global Planet. Change 2004 42 313 321 Google Scholar Crossref Search ADS WorldCat Dill A.R. Klemann A.V. Martinec A.Z. Tesauro A.M. Applying local Green's functions to study the influence of the crustal structure on hydrological loading displacements J. Geodyn. 2015 88 14 22 Google Scholar Crossref Search ADS WorldCat Doake C.S.M. Vaughan D.G. Rapid disintegration of the Wordie Ice Shelf in response to atmospheric warming Nature 1991 350 328 330 Google Scholar Crossref Search ADS WorldCat Dong D. Dickey J.O. Chao Y. Cheng M.K. Geocenter variations caused by atmosphere, ocean and surface ground water Geophys. Res. Lett. 1997 24 1867 1870 Google Scholar Crossref Search ADS WorldCat Dziewonski A.M. Anderson D.L. Preliminary reference Earth model Phys. Earth planet. Inter. 1981 25 4 297 356 Google Scholar Crossref Search ADS WorldCat Eagles G. Gohl K. Larter R.D. Animated tectonic reconstruction of the Southern Pacific and alkaline volcanism at its convergent margins since Eocene times Tectonophysics 2009 464 1 21 29 Google Scholar Crossref Search ADS WorldCat Egbert G.D. Erofeeva S.Y. Han S.C. Luthcke S.B. Ray R.D. Assimilation of GRACE tide solutions into a numerical hydrodynamic inverse model Geophys. Res. Lett. 2009 36 doi:10.1029/2009gl040376 Google Scholar OpenURL Placeholder Text WorldCat Gunter B.C. Didova O. Riva R.E.M. Ligtenberg S.R.M. Lanaerts J.T.M. King M. van den Broeke M.R. Urban T. Empirical estimation of present-day Antarctic glacial isostatic adjustment and ice mass change Cryosphere 2014 8 743 760 Google Scholar Crossref Search ADS WorldCat Herring T. MATLAB Tools for viewing GPS velocities and time series GPS Sol. 2003 7 194 199 Google Scholar Crossref Search ADS WorldCat Herring T.A. King R.W. McClusky S.C. 2010 Cambridge Massachusetts Institute of Technology Documentation for the GAMIT GPS Analysis Software version 10.40 Hillenbrand C.D. et al. Reconstruction of changes in the Weddell Sea sector of the Antarctic Ice Sheet since the Last Glacial Maximum Quatern. Sci. Rev. 2013 100 111 136 Google Scholar Crossref Search ADS WorldCat Hodgson D.A. et al. Exploring former subglacial Hodgson Lake, Antarctica Paper I: site description, geomorphology and limnology Quatern. Sci. Rev. 2009 28 2295 2309 Google Scholar Crossref Search ADS WorldCat Ivins E. Raymond C. James T. The influence of 5000 year-old and younger glacial mass variability on present-day crustal rebound in the Antarctic Peninsula Earth Planets Space 2000 52 11 1023 1029 Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Ivins E.R. Watkins M.M. Yuan D.N. Dietrich R. Casassa G. Rülke A. On-land ice loss and glacial isostatic adjustment at the Drake Passage: 2003–2009 J. geophys. Res.: Solid Earth 1978–2012 2011 116 doi:10.1029/2010JB007607 Google Scholar OpenURL Placeholder Text WorldCat Ivins E.R. James T.S. Wahr J. Schrama O. Ernst J. Landerer F.W. Simon K.M. Antarctic contribution to sea level rise observed by GRACE with improved GIA correction J. geophys. Res.: Solid Earth 2013 118 6 3126 3141 Google Scholar Crossref Search ADS WorldCat Jaldehag R.T. Johansson J.M. Davis J.L. Elósegui P. Geodesy using the Swedish Permanent GPS Network? Effects of snow accumulation on estimates of site positions Geophys. Res. Lett. 1996 23 13 1601 1604 Google Scholar Crossref Search ADS WorldCat Johnson A.C. Swain C.J. Further evidence of fracture-zone induced tectonic segmentation of the Antarctic Peninsula from detailed aeromagnetic anomalies Geophys. Res. Lett. 1995 22 1917 1920 Google Scholar Crossref Search ADS WorldCat Johnson J.S. Everest J.D. Leat P.T. Golledge N.R. Rood D.H. Stuart F.M. The deglacial history of NW Alexander Island, Antarctica, from surface exposure dating Quatern. Res. 2012 77 273 280 Google Scholar Crossref Search ADS WorldCat Johansson J.M. et al. Continuous GPS measurements of postglacial adjustment in Fennoscandia 1. Geodetic results J. geophys. Res.: Solid Earth (1978–2012) 2002 107 B8 do i:10.1029/2001JB000400 Google Scholar OpenURL Placeholder Text WorldCat Kaufmann G. Wu P. Ivins E.R. Lateral viscosity variations beneath Antarctica and their implications on regional rebound motions and seismotectonics J. Geodyn. 2005 39 2 165 181 Google Scholar Crossref Search ADS WorldCat King M.A. Watson C.S. Geodetic vertical velocities affected by recent rapid changes in polar motion Geophys. J. Int. 2014 199 1161 1165 Google Scholar Crossref Search ADS WorldCat King M.A. Bevis M. Wilson T. Johns B. Blume F. Monument-antenna effects on GPS coordinate time series with application to vertical rates in Antarctica J. Geod. 2011a 86 53 63 Google Scholar Crossref Search ADS WorldCat King M.A. Padman L. Nicholls K. Clarke P.J. Gudmundsson G.H. Kulessa B. Shepherd A. Ocean tides in the Weddell Sea: new observations on the Filchner-Ronne and Larsen C ice shelves and model validation J. geophys. Res. 2011b 116 C06006 doi:10.1029/2011JC006949 Google Scholar OpenURL Placeholder Text WorldCat King M.A. Padman L. Nicholls K. Clarke P.J. Gudmundsson G.H. Kulessa B. Shepherd A. Gourmelen N. Correction to: ocean tides in the Weddell Sea: new observations on the Filchner-Ronne and Larsen C ice shelves and model validation J. geophys. Res. 2011c 116 C08026 doi:10.1029/2011JC006949 Google Scholar OpenURL Placeholder Text WorldCat King M.A. Bingham R.J. Moore P. Whitehouse P.L. Bentley M.J. Milne G.A. Lower satellite-gravimetry estimates of Antarctic sea-level contribution Nature 2012 491 586 589 Google Scholar Crossref Search ADS PubMed WorldCat Klemann V. Martinec Z. Ivins E.R. Glacial isostasy and plate motion J. Geodyn. 2008 46 3–5 95 103 Google Scholar Crossref Search ADS WorldCat Kuchar J. Milne G.A. The influence of viscosity structure in the lithosphere on predictions from models of glacial isostatic adjustment J. Geodyn. 2015 86 1 9 Google Scholar Crossref Search ADS WorldCat Kusche J. Schmidt R. Petrovic S. Rietbroek R. Decorrelated GRACE time-variable gravity solutions by GFZ, and their validation using a hydrological model J. Geod. 2009 83 10 903 913 Google Scholar Crossref Search ADS WorldCat Lambeck K. Chappell J. Sea level change through the last glacial cycle Science 2001 292 679 686 Google Scholar Crossref Search ADS PubMed WorldCat Lambeck K. Rouby H. Purcell A. Sun Y. Sambridge M. Sea level and global ice volumes from the Last Glacial Maximum to the Holocene Proc. Natl. Acad. Sci. U.S.A. 2014 111 43 15 296 15 303 Google Scholar Crossref Search ADS WorldCat Latychev K. Mitrovica J.X. Tromp J. Tamisiea M.E. Komatitsch D. Christara C.C. Glacial isostatic adjustment on 3-D Earth models: a finite-volume formulation Geophys. J. Int. 2005 161 421 444 Google Scholar Crossref Search ADS WorldCat Lenaerts J.T.M. van den Broeke M.R. van de Berg W.J. van Meijgaard E. Kuipers Munneke P. A new, high-resolution surface mass balance map of Antarctica (1979–2010) based on regional atmospheric climate modeling Geophys. Res. Lett. 2012 39 L04501 doi:10.1029/2011GL050713 Google Scholar Crossref Search ADS WorldCat Li J. Yuen D.A. Mid-mantle heterogeneities associated with Izanagi plate: implications for regional mantle viscosity Earth planet. Sci. Lett. 2014 385 137 144 Google Scholar Crossref Search ADS WorldCat Ligtenberg S.R.M. Helsen M.M. van den Broeke M.R. An improved semi-empirical model for the densification of Antarctic firn The Cryosphere 2011 5 809 819 Google Scholar Crossref Search ADS WorldCat Martinec Z. Hagedoorn J. The rotational feedback on linear-momentum balance in glacial isostatic adjustment Geophys. J. Int. 2014 199 3 1823 1846 Google Scholar Crossref Search ADS WorldCat Maule C.F. Purucker M.E. Olsen N. Mosegaard K. Heat flux anomalies in Antarctica revealed by satellite magnetic data Science 2005 309 5733 464 467 Google Scholar Crossref Search ADS PubMed WorldCat McCarron J.J. Larter R.D. Late Cretaceous to early Tertiary subduction history of the Antarctic Peninsula J. geol. Soc. 1998 155 2 255 268 Google Scholar Crossref Search ADS WorldCat McMillan M. Shepherd A. Sundal A. Briggs K. Muir A. Ridout A. Hogg A. Wingham D. Increased ice losses from Antarctica detected by CryoSat-2 Geophys. Res. Lett. 2014 41 3899 3905 Google Scholar Crossref Search ADS WorldCat Milne G.A. Mitrovica J.X. Postglacial sea-level change on a rotating Earth Geophys. J. Int. 1998 133 1 1 19 Google Scholar Crossref Search ADS WorldCat Milne G.A. Davis J.L. Mitrovica J.X. Scherneck H.G. Johansson J.M. Vermeer M. Koivula H. Space-geodetic constraints on glacial isostatic adjustment in Fennoscandia Science 2001 291 2381 2385 Google Scholar Crossref Search ADS PubMed WorldCat Mitrovica J.X. Forte A.M. A new inference of mantle viscosity based upon joint inversion of convection and glacial isostatic adjustment data Earth planet. Sci. Lett. 2004 225 1 177 189 Google Scholar Crossref Search ADS WorldCat Mitrovica J.X. Milne G.A. On post-glacial sea level: I. General theory Geophys. J. Int. 2003 154 2 253 267 Google Scholar Crossref Search ADS WorldCat Morelli A. Danesi S. Seismological imaging of the Antarctic continental lithosphere: a review Global Planet. Change. 2004 42 155 165 Google Scholar Crossref Search ADS WorldCat Mulvaney R. et al. Recent Antarctic Peninsula warming relative to Holocene climate and ice-shelf history Nature 2012 489 7414 141 144 Google Scholar Crossref Search ADS PubMed WorldCat Neukom R. et al. Inter-hemispheric temperature variability over the past millennium Nat. Climate Change 2014 4 5 362 367 Google Scholar Crossref Search ADS WorldCat Nield G.A. Whitehouse P.L. King M.A. Clarke P.J. Bentley M.J. Increased ice loading in the Antarctic Peninsula since the 1850s and its effect on glacial isostatic adjustment Geophys. Res. Lett. 2012 39 17 doi:10.1029/2012GL052559 Google Scholar OpenURL Placeholder Text WorldCat Nield G.A. et al. Rapid bedrock uplift in the Antarctic Peninsula explained by viscoelastic response to recent ice unloading Earth planet. Sci. Lett. 2014 397 32 41 Google Scholar Crossref Search ADS WorldCat Orsi A.J. Cornuelle B.D. Severinghaus J.P. Little Ice Age cold interval in West Antarctica: evidence from borehole temperature at the West Antarctic Ice Sheet (WAIS) divide Geophys. Res. Lett. 2012 39 9 doi:10.1029/2012GL051260 Google Scholar OpenURL Placeholder Text WorldCat Paulson A. Zhong S. Wahr J. Inference of mantle viscosity from GRACE and relative sea level data Geophys. J. Int. 2007 171 497 508 Google Scholar Crossref Search ADS WorldCat Peltier W.R. The impulse response of a Maxwell Earth Rev. Geophys. 1974 12 4 649 669 Google Scholar Crossref Search ADS WorldCat Peltier W.R. Global glacial isostasy and the surface of the ice-age Earth: the ICE-5G (VM2) model and GRACE Annu. Rev. Earth planet. Sci. 2004 32 111 149 Google Scholar Crossref Search ADS WorldCat Peltier W.R. Argus D.F. Drummond R. Space geodesy constrains ice-age terminal deglaciation: the global ICE-6G_C (VM5a) model J. geophys. Res.: Solid Earth 2015 120 450 487 Google Scholar Crossref Search ADS WorldCat Rutt I.C. Hagdorn M. Hulton N.R.J. Payne A.J. The Glimmer community ice sheet model J. geophys. Res.: Earth Surface 2003–2012 2009 114 doi:10.1029/2008JF001015 Google Scholar OpenURL Placeholder Text WorldCat Sasgen I. Konrad H. Ivins E.R. Van den Broeke M.R. Bamber J.L. Martinec Z. Klemann V. Antarctic ice-mass balance 2003 to 2012: regional reanalysis of GRACE satellite gravimetry measurements with improved estimate of glacial-isostatic adjustment based on GPS uplift rates Cryosphere 2013 7 1499 1512 Google Scholar Crossref Search ADS WorldCat Shapiro N.M. Ritzwoller M.H. Inferring surface heat flux distributions guided by a global seismic model: particular application to Antarctica Earth planet. Sci. Lett. 2004 223 1 213 224 Google Scholar Crossref Search ADS WorldCat Shepherd A. et al. A Reconciled Estimate of Ice-Sheet Mass Balance Science 2012 338 1183 1189 Google Scholar Crossref Search ADS PubMed WorldCat Sime L.C. Marshall G.J. Mulvaney R. Thomas E.R. Interpreting temperature information from ice cores along the Antarctic Peninsula: ERA40 analysis Geophys. Res. Lett. 2009 36 18 doi:10.1029/2009GL038982 Google Scholar OpenURL Placeholder Text WorldCat Swenson S. Chambers D. Wahr J. Estimating geocenter variations from a combination of GRACE and ocean model output J. geophys. Res. 2008 113 8410 Google Scholar Crossref Search ADS WorldCat Thomas E.R. Bracegirdle T.J. Improving ice core interpretation using in situ and reanalysis data J. geophys. Res., Atmos., 1984–2012 2009 114 doi:10.1029/2009JD012263 Google Scholar OpenURL Placeholder Text WorldCat Thomas E.R. Marshall G.J. McConnell J.R. A doubling in snow accumulation in the western Antarctic Peninsula since 1850 Geophys. Res. Lett. 2008 35 1 doi:10.1029/2007GL032529 Google Scholar OpenURL Placeholder Text WorldCat Thomas I.D. et al. Widespread low rates of Antarctic glacial isostatic adjustment revealed by GPS observations Geophys. Res. Lett. 2011 38 L22302 doi:10.1029/2011GL049277 Google Scholar OpenURL Placeholder Text WorldCat Tregoning P. Herring T.A. Impact of a priori zenith hydrostatic delay errors on GPS estimates of station heights and zenith total delays Geophys. Res. Lett. 2006 33 L23303 doi:10.1029/2006GL027706 Google Scholar Crossref Search ADS WorldCat van der Wal W. Barnhoorn A. Stocchi P. Gradmann S. Wu P. Drury M. Vermeersen B. Glacial isostatic adjustment model with composite 3-D Earth rheology for Fennoscandia Geophys. J. Int. 2013 194 61 77 Google Scholar Crossref Search ADS WorldCat van der Wal W. Whitehouse P.L. Schrama E.J.O. Effect of GIA models with 3D mantle viscosity on GRACE mass balance estimates for Antarctica Earth planet. Sci. Lett. 2015 414 134 143 Google Scholar Crossref Search ADS WorldCat Wahr J. Wingham D. Bentley C. A method of combining ICESat and GRACE satellite data to constrain Antarctic mass balance J. geophys. Res.: Solid Earth (1978–2012) 2000 105 B7 16 279 16 294 Google Scholar Crossref Search ADS WorldCat Whitehouse P. Latychev K. Milne G.A. Mitrovica J.X. Kendall R. Impact of 3-D Earth structure on Fennoscandian glacial isostatic adjustment: implications for space-geodetic estimates of present-day crustal deformations Geophys. Res. Lett. 2006 33 L13502 doi:10.1029/2006GL026568 Google Scholar Crossref Search ADS WorldCat Whitehouse P.L. Bentley M.J. Milne G.A. King M.A. Thomas I.D. A new glacial isostatic adjustment model for Antarctica: calibrated and tested using observations of relative sea-level change and present-day uplift rates Geophys. J. Int. 2012a 190 3 1464 1482 Google Scholar Crossref Search ADS WorldCat Whitehouse P.L. Bentley M.J. Le Brocq A.M. A deglacial model for Antarctica: geological constraints and glaciological modelling as a basis for a new model of Antarctic glacial isostatic adjustment Quatern. Sci. Rev. 2012b 32 1 24 Google Scholar Crossref Search ADS WorldCat Williams S. CATS: GPS coordinate time series analysis software GPS Sol. 2008 12 147 153 Google Scholar Crossref Search ADS WorldCat Wolstencroft M. Shen Z. Törnqvist T.E. Milne G.A. Kulp M. Understanding subsidence in the Mississippi Delta region due to sediment, ice, and water loading: insights from geophysical modeling J. geophys. Res.: Solid Earth 2014 119 3838 3856 Google Scholar Crossref Search ADS WorldCat Wouters B. Martin-Español A. Helm V. Flament T. van Wessem J.M. Ligtenberg S.R.M. van den Broeke M.R. Bamber J.L. Dynamic thinning of glaciers on the Southern Antarctic Peninsula Science 2015 348 6237 899 903 Google Scholar Crossref Search ADS PubMed WorldCat SUPPORTING INFORMATION Additional Supporting Information may be found in the online version of this paper: Figure S1. Modelled vertical elastic rate correction applied to GPS data. Ice load changes from McMillan et al. (2014). Figure S2. Assembly of composite ice model (left to right) W12 + N12 = W12_N12; example showing load change for period 1905–1915. W12 has no load change in Antarctica over the last 150 yr, although the global model was interpolated to account for changes in the northern hemisphere. Scale indicates ice thickness change, not total ice thickness. Figure S3. Top panel: time slices at 50-yr intervals of ice thickness change for the AP1Ka ice history. Note the opposite sign of response at high altitude and low altitude/coastal regions. Bottom panel: total ice-volume change on the AP over 1200 yr of the AP1Ka ice sheet reconstruction; prescribed climatic conditions indicated. AP1Kb has the same spatial and temporal pattern but with 20 per cent of the magnitude. Figure S4. Region masked (dark grey) to produce the W12_*Ev ice histories (Table 3). The jaggedness of the mask is a product of the model resolution. Figure S5. Modelled far field elastic deformation rates and reference frame offset effects for the period 2007–2014. Note variation between Palmer Land GPS sites is <0.05 mm yr–1, almost two orders of magnitude smaller than GPS error. Table S1. Full data-model comparison from which Table 5 (main text) was extracted, sorted by DWRMS, the grey cells indicate an unreasonable D value (outside 0 to −3.5 mm yr−1). Model combination naming: [Ice model]_LT (km)_UMV (1021 Pa s)_LMV (always 1022 Pa s), ‘p’ denotes decimal point (Supplementary Data). Please note: Oxford University Press is not responsible for the content or functionality of any supporting materials supplied by the authors. Any queries (other than missing material) should be directed to the corresponding author for the paper. © The Authors 2015. Published by Oxford University Press on behalf of The Royal Astronomical Society. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted reuse, distribution, and reproduction in any medium, provided the original work is properly cited. © The Authors 2015. Published by Oxford University Press on behalf of The Royal Astronomical Society.
TI - Uplift rates from a new high-density GPS network in Palmer Land indicate significant late Holocene ice loss in the southwestern Weddell Sea
JF - Geophysical Journal International
DO - 10.1093/gji/ggv327
DA - 2015-10-01
UR - https://www.deepdyve.com/lp/oxford-university-press/uplift-rates-from-a-new-high-density-gps-network-in-palmer-land-rPyuQ3gUIo
SP - 737
EP - 754
VL - 203
IS - 1
DP - DeepDyve
ER -