Bathymetry critically influences the intrusion of warm Circumpolar Deep Water onto the continental shelf and ice shelf cavities in Antarctica, thereby forcing ice melting, grounding line retreat, ice discharge and sea level rise. We present a novel and comprehensive bathymetry of Antarctica that includes all its ice shelf cavities and previously uncovered continental shelf areas. The new map is based on a 3D inversion of a circumpolar compilation of gravity anomalies constrained by measurements from the International Bathymetry Chart of the Southern Ocean, BedMachine Antarctica, and discrete seafloor measurements from seismic and ocean robotic probes. Previously unknown troughs with deeper ice shelf cavities are revealed in many parts of Antarctica, especially East Antarctica. The greater depths of troughs on the continental shelf and ice shelf cavities imply that many glaciers are more vulnerable to ocean subsurface warming than previously thought, which will impacts the projections of sea level rise from Antarctica.
AntGG gravity data. Ground-based and airborne Antarctic gravity measurements have been compiled in the frame of the International Association of Geodesy (IAG) Sub-Commission 2.4f “Gravity and Geoid in Antarctica” (AntGG) led by one of the co-authors (MS). A first Antarctic-wide compilation of consistently combined and gridded gravity anomalies (AntGG2016) was published by Scheinert et al 2016a. Since then, a number of new gravity data sets have been made available, especially over regions where in-situ data were sparse or missed at all, e.g. in the South Pole region (Forsberg et al., 2017), in Dronning Maud Land, East Antarctica (Eagles et al., 2018a), over the Pensacola Basin (Paxman et al., 2019a) or in Marie Byrd Land, West Antarctica (Young et al., 2017). Applying an improved processing and taking all available in-situ gravity data into account –- both data already incorporated in AntGG2016 and new data –- the new AntGG compilation features several advantages: the grid resolution is refined from 10 to 5 km, now entire Antarctica (south of 60oS) is covered, and different data sets (in terms of functionals of the disturbing potential, T) are provided.
The improved processing scheme incorporates the so-called remove-compute-restore (RCR) technique, e.g. (Sanso et al., 2013). In the remove step, the observations – mostly gravity anomalies in the case of ground-based measurements, and gravity disturbances in the case of airborne gravimetry – were reduced using a-priori information. This step was done to count for biases and long-wavelengths signal parts, which cannot be resolved over a limited area like Antarctica, as well as to count for short-wavelengths signal parts, which are mainly generated by topography. The latter reduction also leads to a smoothing of the data, which facilitates interpolation and prediction). The reduced signal parts are later restored in terms of the disturbing potential, T, and its desired functionals.
For these two (remove and restore) steps, a high-resolution background model was constructed based on the satellite-only model GOCO-05S (Mayer and Guerr, 2015) and the EARTH2014 topographic model (Hirt et al., 2015). The method to combine the information of these two models is based on a regularization of the satellite normal equation using the topographic model as a-priori information. The combination itself was performed in the spheroidal-harmonic domain, resulting in the model SATOP-1 of up to harmonic degree and order 5480 (see Zingerle et al., 2019a). Then, the necessary evaluation of this background model was done on a 3D grid that covered a latitude range from 60o to 90oS (at 1’ resolution), a longitude range from 0o to 360o (at 2’ (60 to 75.4oS) and 4’ (75.4 to 90oS) resolution) as well as a vertical range from -1 to 6 km (at 200 m resolution). At this grid, the disturbing potential, T, and its functionals (like gravity anomaly) were synthesized. To infer the respective quantities for the remove step (e.g. at specific data points) and the restore step (at points of the final 2D grid at different specific height levels), respectively, a suitable interpolation scheme was applied using this pre-calculated 3D grid. Zingerle et al. (2019a) also discussed the strategy to evaluate and validate ground-based and airborne gravity anomaly data sets, starting with the first AntGG compilation (Scheinert et al., 2016a), using SATOP-1, thereby taking different spectral bands into account.
For the compute step, the method of least-squares collocation (LSC) (elaborated by Moritz et al., 1980, see also Sanso et al., 2013) was applied to the residual observations. LSC provides a prediction method solely based on the stochastic properties of the input data. An advantage is that all needed functionals of the disturbing potential, T, can be inferred at any given point in space (hence, at a regular grid with desired resolution as well as at different height levels), and that also accuracy measures are provided. To circumvent a major problem of LSC, namely the set-up and inversion of the covariance matrix of the (residual) observations, which would have the dimension of the number of observations squared, a localization of the covariance function was introduced. This step is possible since the disturbing potential, T, as harmonic function is defined at or close to a sphere (Zingerle et al., 2021), leading to locally dominating homogeneous-isotropic covariance functions. Similar to the technique described above, the covariance values were pre-calculated at an intermediate regular grid to facilitate simple and fast interpolation algorithms. For the details of this procedure and further approximations (Zingerle et al., 2021). Furthermore, an optimal partitioning of the entire collocation region (hence, Antarctica as a whole) into sub-regions was realized, which significantly reduced computation time and memory requirements (Zingerle et al. 2021). This step was accomplished taking a certain regional overlapping into account as well as minimizing fringe effects. The entire method, called partition-enhanced LSC (PE-LSC), was elaborated in detail by Zingerle et al. (2021) and validated, using also data sets from Antarctica.
As a result, the following data sets are provided at 5 km grid resolution: gravity anomaly, gravity disturbance, Bouguer gravity anomaly, height anomaly, 2nd radial derivative of the disturbing potential and the accuracy measure inferred from LSC. Details on the processing and resulting data sets will be subject to a separate publication (Scheinert et al., in prep.). For the study presented here the gravity anomaly data set was used (shortly: AntGG2022, Scheinert et al., 2021a).
Model setup. For each ice shelf, the model domain comprises an ice layer (density 0.9167 g/cm3), a water layer (density 1.028 g/cm3), and a bedrock layer (density of 2.67 g/cm3). For bedrock elevation beneath the glaciers on land, we use BMv3.7. BMv3.7 has been harmonized with Version 2 of the International Bathymetric Chart of the Southern Ocean (Dorschel et al., 2022). Together with seismic records for individual ice shelves, these data provide observational constrain for the inversion of AntGG2022.
Seismic records include Pine Island, Larsen C, Fimbul, Amery, George VI, Ross, Ronne and Totten ice shelves (Muto et al., 2016; Brisbourne et al., 2021; Smith et al., 2020; Galton-fenzi et al., 2008; Rosier et al., 2018; Vankova et al., 2023). Autonomous underwater vehicle data include Pine Island (Jenkins et al., 2011). We use seafloor depth measurements from the Marine mammals Exploring the Ocean Pole to Pole project (MEOP) (Charrassin et al., 2008). For each available mammal dive, we calculate the maximum depth reached by the mammals during their dive and assemble an Antarctic-wide minimum depth map (Fig. 1C, 1D) on a regular grid by calculating the minimum depth from all dives within a window size of 2500 x 2500 m, or half the resolution of AntGG2022. To capture the general trend in seafloor bathymetry, the grid is filtered by calculating the maximum depth, median, and standard deviation inside a rolling window of 10x10 km. If the 2-sigma value is greater than 100-m, or the accuracy of gravity inversions (Millan et al., 2018;Tinto et al., 2011), we use the dive with the maximum depth within that window. If the 2-sigma value is lower than 100 m (i.e., better than the accuracy of the gravity inversion), we use the average depth within the 10x10 km window. The MEOP grid, at a resolution of 10x10km, is gridded at 2.5 km using a bilinear interpolation and merged with BMv3.7, where no MBES data is available. We choose to resample BMv3.7, AntGG2022 and other related data to 2.5 km to avoid excessive degradation of MBES and BMv3.7 data.
Gravity inversion. We calculate a forward model of the gravity anomalies associated with an initial bedrock elevation from IBCSOv2 using the approach of Parker's 1967 (Oldenburg, 1974) in the Geosoft GM-SYS 3D software. The forward gravity signal is compared to AntGG2022, the bedrock is iteratively migrated until the misfit between modeled and observed gravity anomaly is less than a threshold of 0.1 mGal (Millan et al., 2018 and 2020). To account for spatial variations in gravity signal due to changes in bedrock geology, we used the approach of An et al., (2019), which uses the misfit between observed and modeled gravity at locations of known seafloor depths and interpolates the signal in between. Since the AntGG2022 is complete, we obtain a bathymetry around the entire continent, including all ice shelf cavities.
Error assessment. To assess the precision of the inversion, we quantify the misfit between observed and modeled gravity and convert it into an uncertainty in bedrock elevation with a conversion factor of 5.8 mGal per 100 m (Millan et al., 2017 and 2020). Additionally, we perform an inversion in a rectangular region of 200 km by 125 km at latitude 170.215 and longitude -76.656 , with dense MBES coverage (Fig S1) using only data the periphery and compare the results with actual MBES data and with seismic data, where available.