Hostname: page-component-788cddb947-nxk7g Total loading time: 0 Render date: 2024-10-14T15:23:39.091Z Has data issue: false hasContentIssue false

Amundsen Sea Embayment ice-sheet mass-loss predictions to 2050 calibrated using observations of velocity and elevation change

Published online by Cambridge University Press:  14 August 2023

Suzanne Bevan*
Affiliation:
Department of Geography, Faculty of Science and Engineering, Swansea University, Swansea, SA2 8PP, UK
Stephen Cornford
Affiliation:
School of Geographical Sciences, University of Bristol, Bristol, BS8 1SS, UK
Lin Gilbert
Affiliation:
Mullard Space Science Laboratory, University College London, London, RH5 6NT, UK
Inés Otosaka
Affiliation:
Centre for Polar Observation and Modelling, School of Earth and Environment, University of Leeds, Leeds, LS2 9JT, UK
Daniel Martin
Affiliation:
Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA
Trystan Surawy-Stepney
Affiliation:
Centre for Polar Observation and Modelling, School of Earth and Environment, University of Leeds, Leeds, LS2 9JT, UK
*
Corresponding author: Suzanne Bevan; Email: s.l.bevan@swansea.ac.uk
Rights & Permissions [Opens in a new window]

Abstract

Mass loss from the Amundsen Sea Embayment of the West Antarctic Ice Sheet is a major contributor to global sea-level rise (SLR) and has been increasing over recent decades. Predictions of future SLR are increasingly modelled using ensembles of simulations within which model parameters and external forcings are varied within credible ranges. Accurately reporting the uncertainty associated with these predictions is crucial in enabling effective planning for, and construction of defences against, rising sea levels. Calibrating model simulations against current observations of ice-sheet behaviour enables the uncertainty to be reduced. Here we calibrate an ensemble of BISICLES ice-sheet model simulations of ice loss from the Amundsen Sea Embayment using remotely sensed observations of surface elevation and ice speed. Each calibration type is shown to be capable of reducing the 90% credibility bounds of predicted contributions to SLR by 34 and 43% respectively.

Type
Article
Creative Commons
Creative Common License - CCCreative Common License - BY
This is an Open Access article, distributed under the terms of the Creative Commons Attribution licence (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted re-use, distribution and reproduction, provided the original article is properly cited.
Copyright
Copyright © The Author(s), 2023. Published by Cambridge University Press on behalf of The International Glaciological Society

1. Introduction

Current Antarctic Ice Sheet thinning is concentrated in the Amundsen Sea Embayment (ASE) sector of the West Antarctic Ice Sheet (WAIS), and in Wilkes Basin, East Antarctica (McMillan and others, Reference McMillan2014; Shepherd and others, Reference Shepherd2019; Smith and others, Reference Smith2020). Between 2002 and 2013, the ASE sector contributed an estimated 3.2 ± 0.1 mm to total sea-level rise (SLR) (Sutterley and others, Reference Sutterley2014), with ice loss from the sector increasing by a factor of 5 from the decade ending in 2002 to the decade ending in 2017 (Shepherd and others, Reference Shepherd2019). The loss has been identified as being clearly dynamic in origin rather than being driven by surface mass-balance (SMB) variability (Shepherd and others, Reference Shepherd2019). The areas that are thinning and losing mass are spreading inland, and an estimated 59% of the ASE sector was out of balance by 2016 (Shepherd and others, Reference Shepherd2019). By this time, the two largest glaciers draining the ASE sector, Thwaites and Pine Island Glaciers, were each losing mass at rates of 55 ± 4 and 76 ± 6 Gt a−1 of ice, respectively. Much of the region's ice is grounded below sea level, with the ice–bed interface sloping inland away from the grounding line, making it potentially vulnerable to marine ice-sheet instability (MISI), whereby loss of buttressing at the periphery leads to rapid and runaway grounding line retreat (Hughes, Reference Hughes1973; Schoof, Reference Schoof2007; Joughin and others, Reference Joughin, Smith and Medley2014). Early warning indicators show that two out of three basal melt-rate tipping points for MISI may already have been crossed on Pine Island Glacier, with the third tipping point, one which would lead to a total retreat of the glacier, being an increase in basal melt rate corresponding to an increase in ocean temperatures of 1.2°C (Rosier and others, Reference Rosier2021). MISI in the ASE sector and subsequent collapse could raise global sea levels by 2.3 m (Martin and others, Reference Martin, Cornford and Payne2019).

Various studies have used stand-alone ice-sheet models to estimate contributions to future SLR from the ASE region. The modelled rates of 21st century SLR include: a modal rate value of 0.3 mm a−1 in an ensemble experiment calibrated with observations (Nias and others, Reference Nias, Cornford, Edwards, Gourmelen and Payne2019); total contributions by 2100 of 2.0–4.5 cm by 2100 under various ocean forcing scenarios and initialisation conditions (Alevropoulos-Borrill and others, Reference Alevropoulos-Borrill, Nias, Payne, Golledge and Bingham2020); 1.5–4 cm by 2100 in an experiment that varied future ice-shelf melt rates and meteoric accumulation (Cornford and others, Reference Cornford2015); a median rise by 2100 relative to 2000 of 2 cm using different ocean forcings (Levermann and others, Reference Levermann2020); and a median 50 year SLR of 18.9 mm (with 90% confidence limits of 13.9–24.8 mm) in a study testing ensemble calibration methods (Wernecke and others, Reference Wernecke, Edwards, Nias, Holden and Edwards2020). To put these ASE SLRs into context, current observed 2006–2018 rates of global SLR are 3.69 mm a−1 (3.21–4.18 mm a−1) (Fox-Kemper and others, Reference Fox-Kemper2021) meaning that estimates of ASE near-future contributions are of the order of 10% of global rates.

In considering global and local impacts of and adaptations to rising sea levels, Aschwanden and others (Reference Aschwanden, Bartholomaus, Brinkerhoff and Truffer2021), and Durand and others (Reference Durand2022) emphasise the need for a realistic quantification of the uncertainty associated with SLR predictions. Identifying the sources of uncertainty will, further, allow efforts to reduce them to be directed more effectively. By using a probability density function (pdf) to represent the outcome, for example SLR, of ensembles of ice-sheet model runs, the percentage probability bounds can be identified within which predicted outcomes may be expected to lie. The ensembles may represent different initial states, different parameters, different forcing or different boundary conditions. The resulting model pdfs, or priors, should ideally be conditioned or calibrated with observations to produce posterior pdfs (Nias and others, Reference Nias, Cornford, Edwards, Gourmelen and Payne2019; Aschwanden and others, Reference Aschwanden, Bartholomaus, Brinkerhoff and Truffer2021). Examples of using calibrated pdfs of SLR based on large ensembles of ice-sheet models calibrated with observations to quantify uncertainty include contributions to 2100 SLR from the ASE calibrated using observed rates of surface elevation change (Nias and others, Reference Nias, Cornford, Edwards, Gourmelen and Payne2019; Wernecke and others, Reference Wernecke, Edwards, Nias, Holden and Edwards2020), WAIS contributions to SLR based on a large ensemble calibrated with scores based on modern and reconstructed geological observations (Pollard and others, Reference Pollard, Chang, Haran, Applegate and DeConto2016), a Greenland ice-sheet model calibrated using profiles of surface elevation (Chang and others, Reference Chang, Applegate, Haran and Keller2014) and the use of little ice age constraints to condition an emulated response of the Antarctic Ice Sheet to future high emission scenarios (Gilford and others, Reference Gilford2020). Felikson and others (Reference Felikson2022) used observations of thickness change, mass change and velocity change to compare the impact of calibration type on the posterior distribution of an ensemble of model simulations for the Greenland ice sheet, and Aschwanden and Brinkerhoff (Reference Aschwanden and Brinkerhoff2022) explored the use of a two-stage calibration method for ensemble simulations, also of the Greenland ice sheet.

In this report, we consider the uncertainties in predicted sea-level equivalent (SLE) of ice loss from the ASE region associated with parameter specifications in the BISICLES ice-sheet model. In particular, we consider different parameter values in a regularised Coulomb friction law, in the specification of ice-shelf thinning rates, and in coefficients of viscosity and basal drag. We then use observations to calibrate a large ensemble of simulations that run to 2050. The work builds on the study by Nias and others (Reference Nias, Cornford, Edwards, Gourmelen and Payne2019) in that we now have available an improved bedrock and thickness dataset, implement a new sliding law in the model and as well as elevation change observations we now use velocity observations to calibrate the ensemble.

2. Methods

2.1. BISICLES model

We used the BISICLES ice-flow model (Cornford and others, Reference Cornford2013) to simulate the evolution of the ASE sector of the Antarctic ice sheet to the year 2050. BISICLES is a finite-volume model that incorporates longitudinal and lateral stresses and a simplified treatment of vertical shear based on Schoof and Hindmarsh (Reference Schoof and Hindmarsh2010). Adaptive mesh refinement allows the spatial resolution to be fine in the locality of grounding lines and coarser elsewhere, and to evolve with grounding line migration. Here, meshes are constructed from rectangular patches, each of which is a grid of square cells 4 km, 2 km, 1 km or 500 m across. (Fig. 1).

Figure 1. Model domain and nested mesh regions.

Initial-state datasets included bedrock and surface elevation, and ice thickness from MEaSUREs BedMachine Antarctica Version 2 (Morlighem and others, Reference Morlighem2020), and a three-dimensional temperature field based on the Jet Propulsion Laboratory Ice Sheet System Model submission to initMIP-Antarctica (Seroussi and others, Reference Seroussi2019). Surface accumulation rates were held constant throughout, at the mean of the monthly values from 1981 to 2018 from the Modéle Atmosphérique Régional (MAR) version 3.10 regional climate model (Agosta and others, Reference Agosta2019) (Fig. 2).

Figure 2. (a) Mean and (b) standard deviation of the monthly 1981–2018 Modéle Atmosphérique Régional (MAR) version 3.10 surface mass-balance datasets.

We denote the part of the domain containing floating ice at time t by Ωf(t). Ice-shelf melt rates were imposed indirectly by setting an ice-shelf thickness change rate ∂h/∂tf) < 0, allowing direct comparison with observations, and avoiding the need to invent depth-dependent melt-rate parameterisations. The imposition of ∂h/∂tf) < 0 could represent either ocean melt in excess of that needed to counter ice flow, or a weakening as a result of ice-shelf damage.

In addition, we specified ice fronts based on observations (Fig. 3). The complete ice front corresponds to the MEaSUREs V2 coastline (Mouginot and others, Reference Mouginot, Scheuchl and Rignot2017c), Pine Island Glacier fronts were adjusted to the observed locations in September 2017, October 2018 and February 2020 coinciding with major calving events; and Thwaites Glacier fronts were defined at monthly intervals based on a convolutional neural network applied to Sentinel-1 SAR backscatter images (Surawy-Stepney and others, Reference Surawy-Stepney, Hogg, Cornford and Davison2023).

Figure 3. Ice front positions based on satellite observations and the MEaSUREs V2 coastline.

The formulation of the model requires that a stiffening coefficient (φ, equivalent to an enhancement cofficient E = φ−3) and a basal traction coefficient (β 2) are defined by solving an inverse problem (Cornford and others, Reference Cornford2015). The stiffening coefficient is applied to the vertically integrated effective viscosity to compensate for uncertainties in ice temperature, uncertainties in the rate factor in Glen's flow law and large-scale damage. The goal of the inverse problem is to create smooth continuous fields of φ and effective drag (τ b = β 2|u|), that produce a good match between modelled and observed ice speeds (Cornford and others, Reference Cornford2015).

The optimal combination of β 2 and φ, subject to 0 < φ < 1, is given by minimising the function

(1)$$J = \Vert ( u_{\rm mi} - u_{\rm ob}) \Vert ^2 + {\sigma^2\over \lambda_{\beta^2}^2}\Vert \nabla \beta^2\Vert ^2 + {\sigma^2\over \lambda_\varphi^2}\Vert \nabla \varphi\Vert ^2,\; $$

where u mi and u o are the modelled and observed ice speeds, $\Vert .\Vert$2 means $\int \vert.\vert ^2{\rm d}\Omega$, and σ 2 represents the combined error of modelled and observed velocity. $\lambda _{\beta ^2}^2$ and $\lambda _\varphi ^2$ are regularisation parameters. Their values are chosen through L-curve analysis (Hansen, Reference Hansen2007). Note that simultaneous estimation of β 2 and φ over grounded ice is underdetermined, but is necessary if φ is to be continuous across the grounding line. The regularisation terms, together with the constraint 0 < φ < 1, mitigate against this to some extent (e.g. by prohibiting regions of φ > >1 that would arise in regions of near-zero rate-of-strain) but ultimately mean that an initial guess for both is required and important. Our initial guess assumes that viscous stresses can be neglected over the majority of the ice sheet, so that basal stress τ b equals gravitational driving stress τ d and that the relationship between temperature and rate factor (Cuffey and Paterson, Reference Cuffey and Paterson2010) holds where viscous stresses are unimportant.

The observed ice speeds used for the inversion were the MEaSUREs InSAR-Based Antarctica Ice Velocity Map, Version 2 (Rignot and others, Reference Rignot, Mouginot and Scheuchl2011; Mouginot and others, Reference Mouginot, Scheuchl and Rignot2012, Reference Mouginot, Rignot, Scheuchl and Millan2017a; Rignot aand others, Reference Rignot, Mouginot and Scheuchl2017); those for the ASE region are based mostly on data acquired in 2007, hence we assign 2007 as the start date of our simulations.

The complete model initial state comprises ice thickness h, β 2 and φ, which should be consistent in the sense that the initial u and ∂h/∂t fields are consistent with observations. Simply taking h to be the BedMachine V2 ice thickness results in (subtle) changes to the ice geometry that lead over a few years to substantial reduction in ice flow. To address this, the grounded ice thickness was relaxed over 10 years (see Fig. 4) with the ice-shelf thickness held constant and β 2 and φ periodically recomputed.

Figure 4. (a) Surface elevation from BedMachine V2. (b) Difference in surface elevation from BedMachine V2 after 10 years of model relaxation.

Using the optimal φ and β 2 coefficients and the relaxed geometry, we ran the model to the year 2050, defining the start year to be 2007 and imposing the calving-front retreats described above between 2016 and 2021. For these simulations, we switched to a regularised Coulomb friction law (Joughin and others, Reference Joughin, Smith and Schoof2019):

(2)$$\tau^b = -C_{u_0}\left({\left\vert \bf\,{u}\right\vert \over \left\vert \bf\,{u}\right\vert + u_0}\right)^{1/3},\; $$

where

(3)$$C_{u_0} = \beta^2 u_{\rm mi}^{2/3}\left(u_{\rm mi} + u_0 \right)^{1/3},\; $$

and u mi is the initial model speed, i.e the speed computed when solving the inverse problem.

This formulation of the friction law allows friction to transition towards a plastic Coulomb law at sliding velocities  ≫ u 0, i.e. for fast-flowing ice streams and glaciers with soft deformable sediments at the bed, where a plastic bed is required to realistically model recent observed increases in flow, for example on Pine Island Glacier (De Rydt and others, Reference De Rydt, Reese, Paolo and Gudmundsson2021).

2.2. Ensemble design

An ensemble of model simulations was designed to explore the range in contributions to SLR originating from uncertainties in the basal traction and ice viscosity fields, and rate of ice shelf thinning.

The first set combined the optimal fields of $C_{u_0}$ and φ derived from the inversion, five values of ∂h/∂tf) (0, −2, −5, −10 and −15  m a−1), and five values of u 0 in the friction law (20, 75, 150, 300 and 600 m a−1) – 25 runs. The set of u 0 values are centred on and expand upon the value of 300 m a−1 found to perform well for Pine Island Glacier (Joughin and others, Reference Joughin, Smith and Schoof2019). We then ran model simulations in which we scaled each of the optimal $C_{u_0}$ and φ coefficients for each combination of u 0 and ∂h/∂tf) by factors of 0.9 and 1.1, generating a further 100 parameter sets and outcomes. Using a maximin Latin hypercube sampling, a further set of 100 simulations was drawn, subsampling the four-dimensional parameter space defined by all of the above values of ∂h/∂tf) = (0, − 2, − 5, − 10 and −15 m a−1), u 0 = (20,  75,  150,  300 and 600 m a−1), and $C_{u_0}$ and φ coefficient factors between 0.9 and 1.1. The complete ensemble size consisted of 225 simulations, minus 12 that did not complete.

Finally, we ran a set of simulations with a lower overall SMB than MAR in order to reveal the effect of SMB forcing on ice-sheet dynamics. For these simulations, the SMB was set to 0.3 m a−1 everywhere, we used the optimum $C_{u_0}$ and φ fields, five values of ∂h/∂tf)(0, − 2, − 5, − 10 and −15 m a−1), and u 0 values of 20, 75, 150, 300, 600, 1200 and 2400 m a−1.

This final set of simulations with lower SMB was not included in the ensemble used in model calibration, but is presented simply to investigate the impact of varying SMB. The 0.3 m a−1 SMB, integrated over the model domain, is about 25% lower than the MAR SMB. If accumulating this difference over the full model run time and removing it from the final SLE generated by the MAR SMB simulations produced the same SLE as the SMB = 0.3 m a−1 runs, then we could conclude that an increased SMB has no effect on ice-sheet dynamics.

2.3. Model calibration

Our ensemble of simulations were used to calculate the SLE of ice loss by 2050 relative to 2007, from all ice above floatation within the Thwaites and Pine Island Glacier drainage basins. The ensemble SLEs were used to generate a smoothed pdf for SLE using a Gaussian kernel density estimator with a Silverman rule-of-thumb bandwidth (Silverman, Reference Silverman1986). This pdf, or prior, is denoted P(SLE).

According to Bayes’ theorem a posterior probability distribution P(SLE|O) can be found by weighting or conditioning the prior with a likelihood function P(O|SLE).

(4)$$P\left(SLE \vert O \right)\propto P\left(SLE\right)P\left(O \vert SLE \right).$$

By using a kernel density estimator, we are able to estimate the modes and percentiles of the distributions of prior and posterior SLEs.

Following Nias and others (Reference Nias, Cornford, Edwards, Gourmelen and Payne2019), we substituted for P(O|SLE) a set of weights based on the difference between the model output and observations (see also Edwards and others, Reference Edwards2014; Pollard and others, Reference Pollard, Chang, Haran, Applegate and DeConto2016). First we scored the model simulations against an observed spatial distribution of surface elevation change rate, ∂h/∂t, averaged over 2007–2017 (section 2.4). For each model output $m^{j}_{i}$ at location i from simulation j, we used the corresponding observation o i to calculate a simulation misfit M j:

(5)$$M_j = {1\over n}\sum_{i}{\left(m^{\,j}_{i} - o_{i}\right)^2\over \sigma^{2}_{i}},\; $$

where n is the total number of locations or pixels involved in the summation.

A likelihood score for each simulation, S j, was calculated by:

(6)$$S_j = \exp\left[-{1\over 2} M_j \right].$$

We also scored the simulations against observed surface velocity magnitude for 2007 and for 2020 (section 2.4). In this case, the summation in Eqn (5) is over the two years as well as over space. The area used for both calibrations was restricted to the Thwaites and Pine Island Glacier drainage basins.

Equation (6) assumes that the discrepancies are identically distributed in space, and independent (Edwards and others, Reference Edwards2014). In order to maximise the independence, we down-scaled the spatial resolution of model–observation discrepancies to reduce autocorrelation. Using the R (R Core Team, 2022) function variogram from the gstat package, we plotted the spatial autocorrelation of model–observation discrepancies in ∂h/∂t and u (Fig. 5).

Figure 5. Semi-variograms for model–observation discrepancy in directions 0, 45, 90 and 135 degrees. (a) Mean ∂h/∂t between 2007 and 2017, (b) year 2020 u. The plots are based on the simulation with u 0 = 20 m a−1, and ∂h/∂tf) = −2 m a−1.

For ∂h/∂t the range, the point beyond which the semi-variance no longer increases was ~100 km. The sill was less well determined for u but the semi-variance increases most steeply over the first 25 km (Fig. 5b). Therefore, we resampled the ∂h/∂t and u discrepancies to 100 and 25 km, respectively (Fig. 6).

Figure 6. Example simulation-minus-observation discrepancy maps used for scoring. (a) Mean ∂h/∂t discrepancy at 100 km resolution, (b) 2020 u discrepancy at 25 km resolution. The plots are based on the simulation with u 0 = 20 m a−1, and ∂h/∂tf) = −5 m a−1.

As described in Nias and others (Reference Nias, Cornford, Edwards, Gourmelen and Payne2019), $\sigma ^2_{i}$ in Eqn (5) is the discrepancy variance and it represents how well we expect the model to agree with an observation. $\sigma ^2_{i}$ needs to include both model and observational uncertainty. To set the uncertainty in observations of the 11-year mean ∂h/∂t we relied on the estimates of uncertainty accompanying the supplied data, these estimates are space and time dependent and were summed in quadrature for the period 2007–17. For surface speed, we also used the errors supplied with the data which are dependent on the sensor and the degree of data stacking within the measurement (Mouginot and others, Reference Mouginot, Rignot, Scheuchl and Millan2017a).

Model uncertainty is harder to put a number on, and it is common to choose a model error that is some multiple of observation error (e.g Nias and others, Reference Nias, Cornford, Edwards, Gourmelen and Payne2019; Wernecke and others, Reference Wernecke, Edwards, Nias, Holden and Edwards2020; Felikson and others, Reference Felikson2022). Here, we decided to base the model error on the magnitude of the observation on the basis that observation error is mostly a function of observing technique and available observations, and not necessarily relevant to the magnitude of error we can tolerate in the simulations. However, we can expect and tolerate a larger model error where velocities are high and ∂h/∂t is large. The goal is to avoid concentrating the weights on a small number of models, and yet to provide some discrimination between simulations: a very large discrepancy variance would result in a posterior distribution no narrower than the prior. For each calibration type, we tested model errors equivalent to 1, 5, 10 and 50% of the observed values.

Having calculated the set of likelihood scores, the weight for each ensemble member j is given by the normalised likelihood score:

(7)$$W_j = {S_j\over \sum_{\,j}{S_j}},\; $$

which is P(O|SLE) in Eqn (4).

2.4. Observations of surface elevation change and velocity

Observed surface elevation change rates used in the calibration consist of annual ∂h/∂t rates computed over 3-year periods using data from ERS-1/2, ENVISAT and CryoSat-2. We used data from 2007 to 2017. The method used for the derivation of these elevation change rates is described in Shepherd and others (Reference Shepherd2019). The 2007–17 mean detected ∂h/∂t is shown in Figure 7a.

Figure 7. (a) Observed 2007–2017 mean surface elevation change rates. (b) Example simulated 2007–2017 mean elevation change rates, for u 0 = 20 m a−1 and ∂h/∂tf) = −5 m a−1. Values are masked by the drainage basin boundaries for Pine Island and Thwaites Glaciers.

As well as the MEaSUREs Version 2 ice velocity, we also used the 2019–2020 velocity from the MEaSUREs Version 1, 1 km annual velocities dataset (Mouginot and others, Reference Mouginot, Scheuchl and Rignot2017b) (Fig. 8a), to calibrate the model response. The MEaSUREs Version 2 velocities were used to infer the initial optimum $C_{u_0}$ and φ fields; we would therefore expect these simulations to have a low mismatch with these early velocities meaning we are mainly testing the non-optimal $C_{u_0}$ and φ simulations. By including a second, later in time velocity, we are also testing the ability of all simulations to match the acceleration of the ice.

Figure 8. (a) MEaSUREs V1 2019–2020 surface velocities. (b) Example of 2020 simulated velocities, for u 0 = 20 ma−1 and ∂h/∂tf) = −5 m a−1.

For the purposes of calculating modelled and observed losses of volume of ice above floatation, the region was restricted to the combined Thwaites and Pine Island Glacier drainage basins (basins 21 and 22 in Zwally and others (Reference Zwally, Giovinetto, Beckley and Saba2012), e.g. Fig. 7). Both observed and modelled volume losses were converted to SLE by equating 103 km3 of ice to 2.5 mm of sea level.

3. Results

Examples of modelled ∂h/∂t (mean 2007–17 rate) and u (2020) are shown in Figures 7b and 8b, respectively, where they can be compared with observed values.

The ensemble of simulations described in section 2.2 resulted in losses/gains of ice ranging from 63.7 to −4.4 mm SLE by 2050 (Fig. 9). The histogram of SLEs and the derived pdf is shown in all panels of Figure 11. The median SLE of the ensemble of runs was 16.2 mm. The 5th and 95th percentiles of the pdf are −0.2 and 39.1 mm SLE, respectively, and the modal value is 10.7 mm. The observation-based SLE of 14.5 mm by 2050 is based on a linear extrapolation of the mean rate of observed volume change between 2007 and 2017 (black dashed line in Fig. 11).

Figure 9. Time evolution of SLE ice loss for the large ensemble described in section 2.2. The observations are based on observed dh/dt.

3.1. Calibration using ∂h/∂t and u

We used pair plots of SLE, and S j for ∂h/∂t and u to trial the different scaling factors for model error. Using a 10% scaling produced reasonable results in that scores are non-zero over a large number of simulations (Figs 10a,b) and the two scoring methods are strongly correlated with a Pearson's correlation coefficient of 0.90, significant at the 99% level (Figs 10c). Figure 10 also highlights that the simulations with optimum $C_{u_0}$ and φ score highly. The latter being true indicates that these simulations are capturing well both the rate of mass loss and the acceleration of ice during the calibration period.

Figure 10. Pair plots of the full ensemble, axis labels are given by the text in the diagonal panels. (a) S j (vel) distribution as a function of SLE. (b) S j (dh/dt) distribution as a function of SLE. (c) Relationship between S j (vel) and S j (dh/dt). SLE is the 2050 SLE of mass loss. S j (vel) is the score for each simulation calculated using velocity and a model error equivalent to 10% of the observed value of velocity. S j (dh/dt) is the score using ∂h/∂t and a model error also equivalent to 10% of the observed value of ∂h/∂t. CC = 0.9 is the Pearson's correlation coefficient between S j (dh/dt) and S j (vel).

All ten of the highest scoring simulations scored using ∂h/∂t are with optimum $C_{u_0}$ and φ values, of these five have a shelf thinning rate of 10 m a−1, the other five a rate of 15 m a−1. Of the top ten highest scoring simulations scored by u, six have optimum $C_{u_0}$ and φ values but have shelf thinning rates ranging from 0 to 10 m a−1.

Weighting the prior pdfs using likelihood scores based on observations of ∂h/∂t and u and with model errors set equal to 10% of the observed values narrows the 90% credibility bands of likely 2050 SLE of ice loss to 1.7–27.5 and 2.8–25.3 mm (Fig. 11a); with median (mode) values of 13.2 and 13.0 mm (12.1 and 12.1 mm) (Table 1). Increasing (decreasing) the allowance for model error broadens (narrows) the credibility range (Figs 11b,c). Note that in Figure 11b, the calibrated curve with model error equal to 1% of observed value is barely distinguishable from the 10% curve. The calibration process has the effect of suppressing many of the simulations that resulted in extreme high SLEs.

Figure 11. Normalised histogram and prior probability density functions of the full ensemble of simulations. The vertical black dashed line shows the mean 2007–17 rate of observed SLE volume loss extrapolated to 2050. (a) With ∂h/∂t and u calibrations for model error equal to 10% of the observational error. (b) Calibrated SLE pdfs for the ∂h/∂t calibration method with varying model error. (c) As for (b) but for the u calibration.

Table 1. Ensemble statistics for prior and posterior (calibrated) SLE of ice loss shown in Figure 11a

3.2. Sensitivity of total SLE to the combined effect of varying u 0, ∂h/∂tf) and SMB

Figure 12 shows the 2050 SLE of ice loss according to u 0 and ∂h/∂tf) for the MAR SMB and SMB =0.3 m a−1 fields. Removing the excess SMB that MAR provides compared with 0.3 m a−1 everywhere, an excess equivalent to 6.7 mm SLE by 2050, shows that the dynamic forcing (shown in black in Fig. 12) provided by this excess SMB is of the order of 10% of total discharge.

Figure 12. Total sea-level equivalent of loss of ice above floatation for varying u 0 values, grouped by the imposed rate of thinning of floating ice (∂h/∂tf)). The grey bars show the SLE using the mean MAR SMB, the grey-plus-coloured bars show the SLE for SMB = 0.3 m a−1, and the black sections show the additional SLE after the excess (MAR - 0.3 m a−1) accumulation has been removed from the MAR SMB simulations. (Note that the simulations for ∂h/∂tf) = −15 m a−1 with MAR SMB and ∂h/∂tf) = −10 m a−1 with SMB = 0.3 m a−1, failed to complete.)

Figure 12 also reveals that for ∂h/∂tf) ≥ −5 m a−1, increasing u o results in an increased total SLE by 2050. At higher thinning rates the reverse is true.

4. Discussion

4.1. Sliding law parameters, shelf thinning rates and SMB

Whilst the aim of these experiments was to demonstrate the utility of using observations to reduce uncertainty in modelling results, we can also consider which model parameter sets lead to simulations that agree best with extrapolated observations of current ASE ice loss. With the optimum $C_{u_0}$ and φ fields and the MAR SMB, the simulations that best agree with the extrapolated observed SLE of 14.5 mm by 2050 are those with ∂h/∂tf), of between −5 and −10 m a−1 (Fig. 12). This shelf thinning rate is greater than observed thickness changes of −19.4 m/decade for the region (Paolo and others, Reference Paolo, Fricker and Padman2015). However, as indicated in section 2.1, the imposed thinning rates also act as a proxy for ice-shelf damage, weakening the back stress buttressing the inland grounded ice. Damage such as fractures and crevasses in shear zones is observed to be increasing on both Thwaites and Pine Island Glacier ice shelves (Alley and others, Reference Alley, Scambos, Alley and Holschuh2019; Lhermitte and others, Reference Lhermitte2020) and may explain why we need to impose quite high thinning rates to keep pace with extrapolated observed ice loss rates. A model that evolves damage could circumvent this requirement for artificially high melt rates.

Within the range of values tested, the specification of shelf-thinning rate has a much greater effect on ice loss than the u 0 value in the sliding law (Fig. 12). Recent projections for 2100 under RCP8.5 indicate ice-shelf melt rates increasing by factors of 1.4–2.2 compared with present day (Jourdain and others, Reference Jourdain, Mathiot, Burgard, Caillet and Kittel2022). A doubling of shelf-thinning rates from 5 to 10 m a−1 in these simulations increases 2050 SLE ice loss by about 30%.

The correlation between the magnitude of u 0 in the regularised Coulomb friction law and the total SLE switches sign between ∂h/∂tf) = −5 m a−1 and ∂h/∂tf) = −10 m a−1 (Fig. 12). This relationship indicates that greater bed plasticity leads to greater sensitivity to loss of buttressing. At higher shelf thinning rates, greater overall volume loss is seen when u 0 is lower; that is, when a larger portion of the ice experiences plastic, Coulomb-like drag. At lower thinning rates, the converse is true. In the higher melt-rate simulations, the ice tends to flow faster, with the speed increase limited more by the increasing basal traction in Weertman-like cases (large u 0). At lower thinning rates, the ice tends to slow down, and in these cases the speed decrease is limited more by the decreasing basal traction in Weertman-like models. In everyday terms, the higher u 0, the more the bed will both increase its grip on the ice as it accelerates, and decrease its grip when the ice slows.

Simulations run with an SMB of 0.3 m a−1 everywhere, corresponding to a reduction in total accumulated mass of 61 Gt a−1 by 2050, or about 30% compared with the MAR SMB, resulted in around a 10% decrease in the dynamic discharge of ice to the ocean. These experiments with reduced SMB suggest that the increases in SMB predicted for the sector by 2080–2100 of 30–40% (Donat-Magnin and others, Reference Donat-Magnin2021) may be slightly offset by increases in discharge forced by the additional accumulation.

4.2. Large ensemble

The results of the large ensemble (Fig. 9) exhibit a wide range of 2050 SLE ice loss and encompass current observations, indicating that we have adequately sampled the parameter space; a conclusion confirmed by the calibrated curves (Fig. 11). Confidence in the calibration method is strengthened by the fact that the two independently calibrated posterior distributions are in good agreement with each other in terms of the maximum likelihood SLE. This agreement between calibration types is in contrast to an experiment using an ensemble of model simulations for the Greenland ice sheet (Felikson and others, Reference Felikson2022). In that experiment using velocity change to calibrate the ensemble resulted in a distinctly higher SLE by 2100 compared with using thickness change, which actually indicated that the ice sheet was gaining mass. There are many differences between that study and ours, but the authors concluded that using a firn densification model to convert observed thickness change to dynamic thickness may have had an effect on the calibration. In our methodology, as the majority of the region is in dynamic imbalance, we make the assumption that observed elevation changes are dominated by dynamics and thus avoid the need for any firn modelling.

The credibility bands of the calibrated curves are affected by the chosen σ values. We have little knowledge of the model structural error but using spatially variable estimates of observational uncertainty concentrates the contribution of grid cells to the likelihood values on regions where observational uncertainties are lowest. Increasing model error to be a greater proportion of observed values broadens the confidence intervals, and vice versa. Parameter pair plots (Fig. 10) confirm that for both ∂h/∂t and u, including a model σ value equal to 10% of the observed ∂h/∂t and u results in scores that avoid heavily weighting only a small number of simulations and yet do provide a useful constraint on the ensemble. Assuming a very large permissible model error would effectively make the calibration process without value.

The full 213-member ensemble included varying $C_{u_0}$ and φ by factors ranging between 0.9 and 1.1 to explore underdetermination in the initialisation. The calibration process gives low credibility to the simulations that result in extreme SLEs of ice loss. At the low end, six of the seven simulations resulting in a negative SLE were a result of setting both $C_{u_0}$ and φ to 1.1 times the optimum value obtained by the inverse problem. Similarly, at the high end, six of the seven highest SLE simulations were a result of setting both $C_{u_0}$ and φ to 0.9 times the optimum value. Thus, together with the result that the highest scoring simulations have optimum C and φ, we conclude that the calibration process adds confidence to the solution obtained to the inverse problem.

The calibrated curves assign a higher likelihood to 2050 SLE values lower than those indicated by an extrapolation of observations. This discrepancy could be a result of imposing calving front retreats for years that overlap with the calibration periods: during these periods, the simulations can match the observed rates of loss through a combination of shelf thinning and loss of shelf area. Beyond 2021, we impose no calving-front retreat and back-stress reduction is then achieved only by shelf thinning, resulting in a lower 2050 SLE compared with extrapolating observations. Alternatively, a simple extrapolation of observed mass loss cannot take into account any changes in grounding line retreat rates owing to upstream variations in bed geometry. In addition, the model is relaxing towards a stable equilibrium in which losses decay over time, whereas extrapolating observations assumes that the current imbalance is preserved.

5. Conclusions

An ensemble of decadal-scale simulations of ice-sheet evolution in the ASE region of WAIS was created using the BISICLES model, with a regularised Coulomb friction law, and with sliding and rheology parameters defined via the solution of the inverse problem.

Reproducing current observed mass loss rates requires the inclusion of some form of reduction of ice-shelf buttressing. In the experiments presented in this study, that loss of buttressing was supplied by a combination of imposed ice-shelf thinning and, in the early years, a short period of calving front retreats based on observations. Simulations with shelf thinning rates varying from 0 to 15 m a−1 highlighted that ice loss by 2050 was more sensitive to ice-shelf thinning or damage, than to the chosen value of u 0 in the sliding law.

Simulations also revealed that SMB changes have an impact on ice flow dynamics that may partially offset any predicted future increases in SMB.

An expanded ensemble including up to 10% variation in sliding and viscosity parameters was then calibrated using spatial fields of misfits between model output and remotely sensed observations of changes in surface elevation and ice speed. Calibration by u and by ∂h/∂t produced consistent posterior likelihood distributions, and reduced uncertainty compared with the prior distribution by 43 and 34%, respectively.

Data

The surface elevation change observational datasets are available at https://zenodo.org/record/8117577 (doi: 10.5281/zenodo.8117577). A summary file of the ensemble simulated annual SLE values, and annual NetCDF files of land-ice thickness, and the u and v components of velocity at 1 km spatial resolution are available at https://zenodo.org/record/8120460 (doi: 10.5281/zenodo.8120460).

Acknowledgements

This publication was supported by PROTECT. This project has received funding from the European Union's Horizon 2020 research and innovation programme under grant agreement No 869304, PROTECT contribution number 71. Support for this work was provided through the Scientific Discovery through Advanced Computing (SciDAC) program funded by the U.S. Department of Energy (DOE), Office of Science, Biological and Environmental Research and Advanced Scientific Computing Research programs, as a part of the ProSPect SciDAC Partnership. Work at Berkeley Lab was supported by the Director, Office of Science, of the U.S. Department of Energy under Contract No. DE-AC02-05CH11231. This research used resources of the National Energy Research Scientific Computing Center (NERSC), a U.S. Department of Energy Office of Science User Facility located at Lawrence Berkeley National Laboratory, operated under Contract No. DE-AC02-05CH11231 using NERSC award ASCR-ERCAPm1041.

Author contributions

Suzanne Bevan designed and ran the model simulations, analysed the results and wrote most of the text. Stephen Cornford contributed to data analysis and writing of the text. Lin Gilbert, Inés Otosaka and Trystan Surawy-Stepney contributed data. Stephen Cornford and Daniel Martin supported Suzanne Bevan in use of the BISICLES model and Berkeley Lab supercomputing facilities.

References

Agosta, C and 10 others (2019) Estimation of the Antarctic surface mass balance using the regional climate model MAR (1979–2015) and identification of dominant processes. The Cryosphere 13(1), 281296. doi:10.5194/tc-13-281-2019CrossRefGoogle Scholar
Alevropoulos-Borrill, AV, Nias, IJ, Payne, AJ, Golledge, NR and Bingham, RJ (2020) Ocean-forced evolution of the Amundsen sea catchment, West Antarctica, by 2100. The Cryosphere 14(4), 12451258. doi:10.5194/tc-14-1245-2020CrossRefGoogle Scholar
Alley, KE, Scambos, TA, Alley, RB and Holschuh, N (2019) Troughs developed in ice-stream shear margins precondition ice shelves for ocean-driven breakup. Science Advances 5(10), eaax2215. doi:10.1126/sciadv.aax2215CrossRefGoogle ScholarPubMed
Aschwanden, A, Bartholomaus, TC, Brinkerhoff, DJ and Truffer, M (2021) Brief communication: a roadmap towards credible projections of ice sheet contribution to sea level. The Cryosphere 15(12), 57055715. doi:10.5194/tc-15-5705-2021CrossRefGoogle Scholar
Aschwanden, A and Brinkerhoff, DJ (2022) Calibrated mass loss predictions for the Greenland Ice Sheet. Geophysical Research Letters 49(19), e2022GL099058. doi:10.1029/2022GL099058CrossRefGoogle Scholar
Chang, W, Applegate, PJ, Haran, M and Keller, K (2014) Probabilistic calibration of a Greenland Ice Sheet model using spatially resolved synthetic observations: toward projections of ice mass loss with uncertainties. Geoscientific Model Development 7(5), 19331943. doi:10.5194/gmd-7-1933-2014CrossRefGoogle Scholar
Cornford, SL, and 8 others (2013) Adaptive mesh, finite volume modeling of marine ice sheets. Journal of Computational Physics 232(1), 529549. doi:10.1016/j.jcp.2012.08.037CrossRefGoogle Scholar
Cornford, SL and 14 others (2015) Century-scale simulations of the response of the West Antarctic Ice Sheet to a warming climate. The Cryosphere 9(4), 15791600. doi:10.5194/tc-9-1579-2015CrossRefGoogle Scholar
Cuffey, K and Paterson, WSB (2010) The Physics of Glaciers. Amsterdam, the Netherlands: Academic Press.Google Scholar
De Rydt, J, Reese, R, Paolo, FS and Gudmundsson, GH (2021) Drivers of Pine Island Glacier speed-up between 1996 and 2016. The Cryosphere 15(1), 113132. doi:10.5194/tc-15-113-2021CrossRefGoogle Scholar
Donat-Magnin, M and 7 others (2021) Future surface mass balance and surface melt in the Amundsen sector of the West Antarctic Ice Sheet. The Cryosphere 15(2), 571593. doi:10.5194/tc-15-571-2021CrossRefGoogle Scholar
Durand, G and 17 others (2022) Sea-level rise: from global perspectives to local services. Frontiers in Marine Science 8, 709595. doi:10.3389/fmars.2021.709595CrossRefGoogle Scholar
Edwards, TL and 12 others (2014) Probabilistic parameterisation of the surface mass balance–elevation feedback in regional climate model simulations of the Greenland ice sheet. The Cryosphere 8(1), 181194. doi:10.5194/tc-8-181-2014CrossRefGoogle Scholar
Felikson, D and 6 others (2022) Choice of observation type affects Bayesian calibration of ice sheet model projections. EGUsphere 118. publisher: Copernicus GmbH. doi:10.5194/egusphere-2022-1213Google Scholar
Fox-Kemper, B and 17 others (2021) Ocean, cryosphere and sea level change. In Masson-Delmotte V, Zhai P, Pirani A, Connors S, Péan C, Berger S, Caud N, Chen Y, Goldfarb L, Gomis M, Huang M, Leitzell K, Lonnoy E, Matthews J, Maycock T, Waterfield T, Yelekçi O, Yu R and Zhou B (eds.), Climate Change 2021: The Physical Science Basis. Contribution of Working Group I to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change, Cambridge, United Kingdom and New York, NY, USA: Cambridge University Press, pp. 1211–1362.Google Scholar
Gilford, DM and 5 others (2020) Could the last interglacial constrain projections of future Antarctic ice mass loss and sea-level rise?. Journal of Geophysical Research: Earth Surface 125(10), 005418. doi:10.1029/2019JF005418Google Scholar
Hansen, PC (2007) Regularization tools version 4.0 for Matlab 7.3. Numerical Algorithms 46(2), 189194. doi:10.1007/s11075-007-9136-9CrossRefGoogle Scholar
Hughes, T (1973) Is the West Antarctic ice sheet disintegrating?. Journal of Geophysical Research 78(33), 78847910. doi:10.1029/JC078i033p07884CrossRefGoogle Scholar
Joughin, I, Smith, BE and Medley, B (2014) Marine ice sheet collapse potentially under way for the Thwaites glacier basin, West Antarctica. Science 344(6185), 735738. doi:10.1126/science.1249055CrossRefGoogle ScholarPubMed
Joughin, I, Smith, BE and Schoof, CG (2019) Regularized Coulomb friction laws for ice sheet sliding: application to pine Island Glacier, Antarctica. Geophysical Research Letters 46(9), 47644771. doi:10.1029/2019GL082526CrossRefGoogle ScholarPubMed
Jourdain, NC, Mathiot, P, Burgard, C, Caillet, J and Kittel, C (2022) Ice shelf basal melt rates in the Amundsen sea at the end of the 21st century. Geophysical Research Letters 49(22), e2022GL100629. doi:10.1029/2022GL100629CrossRefGoogle Scholar
Levermann, A and 36 others (2020) Projecting Antarctica's contribution to future sea level rise from basal ice shelf melt using linear response functions of 16 ice sheet models (LARMIP-2). Earth System Dynamics 11(1), 3576. doi:10.5194/esd-11-35-2020CrossRefGoogle Scholar
Lhermitte, S and 7 others (2020) Damage accelerates ice shelf instability and mass loss in Amundsen Sea Embayment. Proceedings of the National Academy of Sciences 117(40), 2473524741. doi:10.1073/pnas.1912890117CrossRefGoogle ScholarPubMed
Martin, DF, Cornford, SL and Payne, AJ (2019) Millennial-scale vulnerability of the Antarctic ice sheet to regional ice shelf collapse. Geophysical Research Letters 46(3), 14671475. doi:10.1029/2018gl081229CrossRefGoogle Scholar
McMillan, M and 7 others (2014) Increased ice losses from Antarctica detected by CryoSat-2. Geophysical Research Letters 41(11), 38993905. doi:10.1002/2014GL060111CrossRefGoogle Scholar
Morlighem, M and 36 others (2020) Deep glacial troughs and stabilizing ridges unveiled beneath the margins of the Antarctic ice sheet. Nature Geoscience 13(2), 132137. bandiera_abtest: a Cg_type: Nature Research Journals Number: 2 Primary_atype: Research Publisher: Nature Publishing Group Subject_term: Cryospheric science;Geomorphology;Projection and prediction Subject_term_id: cryospheric-science;geomorphology;projection-and-prediction. doi:10.1038/s41561-019-0510-8CrossRefGoogle Scholar
Mouginot, J, Rignot, E, Scheuchl, B and Millan, R (2017a) Comprehensive annual ice sheet velocity mapping using Landsat-8, Sentinel-1, and RADARSAT-2 data. Remote Sensing 9(4), 364. doi:10.3390/rs9040364CrossRefGoogle Scholar
Mouginot, J, Scheuchl, B and Rignot, E (2012) Mapping of ice motion in Antarctica using synthetic-aperture radar data. Remote Sensing 4(9), 27532767. doi:10.3390/rs4092753CrossRefGoogle Scholar
Mouginot, J, Scheuchl, B and Rignot, E (2017b) MEaSUREs Annual Antarctic Ice Velocity Maps, Version 1 [Data Set]. Boulder, Colorado USA. NASA National Snow and Ice Data Center Distributed Active Archive Center (doi: 10.5067/9T4EPQXTJYW9), Date Accessed: 12-12-2021.CrossRefGoogle Scholar
Mouginot, J, Scheuchl, B and Rignot, E (2017c) MEaSUREs Antarctic Boundaries for IPY 2007–2009 from Satellite Radar, Version 2. Boulder, Colorado, USA: NASA National Snow and Ice Data Center Distributed Active Archive Center (doi: 10.5067/AXE4121732AD), Date Accessed: 06-06-2022.CrossRefGoogle Scholar
Nias, IJ, Cornford, SL, Edwards, TL, Gourmelen, N and Payne, AJ (2019) Assessing uncertainty in the dynamical ice response to ocean warming in the Amundsen Sea Embayment, West Antarctica. Geophysical Research Letters 46(20), 1125311260. doi:10.1029/2019GL084941CrossRefGoogle Scholar
Paolo, FS, Fricker, HA and Padman, L (2015) Volume loss from Antarctic ice shelves is accelerating. Science 348(6232), 327331. doi:10.1126/science.aaa0940CrossRefGoogle ScholarPubMed
Pollard, D, Chang, W, Haran, M, Applegate, P and DeConto, R (2016) Large ensemble modeling of the last deglacial retreat of the West Antarctic Ice Sheet: comparison of simple and advanced statistical techniques. Geoscientific Model Development 9(5), 16971723. doi:10.5194/gmd-9-1697-2016CrossRefGoogle Scholar
R Core Team (2022) R: a language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria, (https://www.R-project.org/).Google Scholar
Rignot, E, Mouginot, J and Scheuchl, B (2011) Ice flow of the Antarctic ice sheet. Science 333(6048), 14271430. doi:10.1126/science.1208336CrossRefGoogle ScholarPubMed
Rignot, E, Mouginot, J and Scheuchl, B (2017) MEaSUREs InSAR-Based Antarctica Ice Velocity Map, Version 2. Boulder, Colorado, USA: NASA DAAC at the National Snow and Ice Data Center (doi: 10.5067/D7GK8F5J8M8R), Date Accessed: 01-02-2022.CrossRefGoogle Scholar
Rosier, SHR and 5 others (2021) The tipping points and early warning indicators for Pine Island Glacier, West Antarctica. The Cryosphere 15(3), 15011516. doi:10.5194/tc-15-1501-2021CrossRefGoogle Scholar
Schoof, C (2007) Ice sheet grounding line dynamics: steady states, stability, and hysteresis. Journal of Geophysical Research: Earth Surface 112(F3), F03S28. doi:10.1029/2006jf000664CrossRefGoogle Scholar
Schoof, C and Hindmarsh, RCA (2010) Thin-film flows with wall slip: an asymptotic analysis of higher order glacier flow models. The Quarterly Journal of Mechanics and Applied Mathematics 63(1), 73114. doi:10.1093/qjmam/hbp025CrossRefGoogle Scholar
Seroussi, H and 38 others (2019) initMIP-Antarctica: an ice sheet model initialization experiment of ISMIP6. The Cryosphere 13(5), 14411471. doi:10.5194/tc-13-1441-2019CrossRefGoogle Scholar
Shepherd, A and 9 others (2019) Trends in Antarctic ice sheet elevation and mass. Geophysical Research Letters 46(14), 81748183. doi:10.1029/2019GL082182CrossRefGoogle ScholarPubMed
Silverman, BW (1986) Density Estimation for Statistics and Data Analysis. London, UK: CRC Press.Google Scholar
Smith, B and 14 others (2020) Pervasive ice sheet mass loss reflects competing ocean and atmosphere processes. Science 368(6496), 12391242. doi:10.1126/science.aaz5845CrossRefGoogle ScholarPubMed
Surawy-Stepney, T, Hogg, AE, Cornford, SL and Davison, BJ (2023) Episodic dynamic change linked to damage on the Thwaites Glacier Ice Tongue. Nature Geoscience 16(1), 3743. doi:10.1038/s41561-022-01097-9CrossRefGoogle Scholar
Sutterley, TC and 7 others (2014) Mass loss of the Amundsen Sea Embayment of West Antarctica from four independent techniques. Geophysical Research Letters 41(23), 84218428. doi:10.1002/2014GL061940CrossRefGoogle Scholar
Wernecke, A, Edwards, TL, Nias, IJ, Holden, PB and Edwards, NR (2020) Spatial probabilistic calibration of a high-resolution Amundsen Sea Embayment ice sheet model with satellite altimeter data. The Cryosphere 14(5), 14591474. doi:10.5194/tc-14-1459-2020CrossRefGoogle Scholar
Zwally, HJ, Giovinetto, MB, Beckley, MA and Saba, JL (2012) Antarctic and Greenland Drainage Systems. Greenbelt, Maryland, US: GSFC Cryospheric Sciences Laboratory. https://earth.gsfc.nasa.gov/cryo/data/polar-altimetry/antarctic-and-greenland-drainage-systems.Google Scholar
Figure 0

Figure 1. Model domain and nested mesh regions.

Figure 1

Figure 2. (a) Mean and (b) standard deviation of the monthly 1981–2018 Modéle Atmosphérique Régional (MAR) version 3.10 surface mass-balance datasets.

Figure 2

Figure 3. Ice front positions based on satellite observations and the MEaSUREs V2 coastline.

Figure 3

Figure 4. (a) Surface elevation from BedMachine V2. (b) Difference in surface elevation from BedMachine V2 after 10 years of model relaxation.

Figure 4

Figure 5. Semi-variograms for model–observation discrepancy in directions 0, 45, 90 and 135 degrees. (a) Mean ∂h/∂t between 2007 and 2017, (b) year 2020 u. The plots are based on the simulation with u0 = 20 m a−1, and ∂h/∂tf) = −2 m a−1.

Figure 5

Figure 6. Example simulation-minus-observation discrepancy maps used for scoring. (a) Mean ∂h/∂t discrepancy at 100 km resolution, (b) 2020 u discrepancy at 25 km resolution. The plots are based on the simulation with u0 = 20 m a−1, and ∂h/∂tf) = −5 m a−1.

Figure 6

Figure 7. (a) Observed 2007–2017 mean surface elevation change rates. (b) Example simulated 2007–2017 mean elevation change rates, for u0 = 20 m a−1 and ∂h/∂tf) = −5 m a−1. Values are masked by the drainage basin boundaries for Pine Island and Thwaites Glaciers.

Figure 7

Figure 8. (a) MEaSUREs V1 2019–2020 surface velocities. (b) Example of 2020 simulated velocities, for u0 = 20 ma−1 and ∂h/∂tf) = −5 m a−1.

Figure 8

Figure 9. Time evolution of SLE ice loss for the large ensemble described in section 2.2. The observations are based on observed dh/dt.

Figure 9

Figure 10. Pair plots of the full ensemble, axis labels are given by the text in the diagonal panels. (a) Sj (vel) distribution as a function of SLE. (b) Sj (dh/dt) distribution as a function of SLE. (c) Relationship between Sj (vel) and Sj (dh/dt). SLE is the 2050 SLE of mass loss. Sj (vel) is the score for each simulation calculated using velocity and a model error equivalent to 10% of the observed value of velocity. Sj (dh/dt) is the score using ∂h/∂t and a model error also equivalent to 10% of the observed value of ∂h/∂t. CC = 0.9 is the Pearson's correlation coefficient between Sj (dh/dt) and Sj (vel).

Figure 10

Figure 11. Normalised histogram and prior probability density functions of the full ensemble of simulations. The vertical black dashed line shows the mean 2007–17 rate of observed SLE volume loss extrapolated to 2050. (a) With ∂h/∂t and u calibrations for model error equal to 10% of the observational error. (b) Calibrated SLE pdfs for the ∂h/∂t calibration method with varying model error. (c) As for (b) but for the u calibration.

Figure 11

Table 1. Ensemble statistics for prior and posterior (calibrated) SLE of ice loss shown in Figure 11a

Figure 12

Figure 12. Total sea-level equivalent of loss of ice above floatation for varying u0 values, grouped by the imposed rate of thinning of floating ice (∂h/∂tf)). The grey bars show the SLE using the mean MAR SMB, the grey-plus-coloured bars show the SLE for SMB = 0.3 m a−1, and the black sections show the additional SLE after the excess (MAR - 0.3 m a−1) accumulation has been removed from the MAR SMB simulations. (Note that the simulations for ∂h/∂tf) = −15 m a−1 with MAR SMB and ∂h/∂tf) = −10 m a−1 with SMB = 0.3 m a−1, failed to complete.)