Trends and Uncertainties of Regional Barystatic Sea-level Change in the Satellite Altimetry Era

Ocean mass change is one of the main drivers of present-day sea-level change (SLC). Also known as barystatic SLC, it is driven by the exchange of freshwater between the land and the ocean, such as melting of continental ice from glaciers and ice sheets, and variations in land water storage. While many studies have quantified the present-day barystatic contribution to global mean SLC, fewer works have looked into regional changes. This study provides a comprehensive analysis of regional 5 barystatic SLC trends since 1993 (the satellite altimetry era), with a focus on the uncertainty budget. We consider three types of uncertainties: intrinsic (the uncertainty from the data/model itself); temporal (related to the temporal variability in the time series); and spatial-structural (related to the location/distribution of the mass change sources). We collect a range of estimates for the individual freshwater sources, which are used to compute regional patterns (fingerprints) of barystatic SLC and analyse the different types of uncertainty. When all the contributions are combined, we find that the barystatic sea-level 10 trends regionally ranges from −0.43 to 2.55 mm.year−1 for 2003-2016, and from −0.39 to 2.00 mm.year−1 for 1993-2016, depending on the choice of dataset. When all types of uncertainties from all contributions are combined, the total barystatic uncertainties regionally range from 0.62 to 1.29 mm.year−1 for 2003-2016, and from 0.35 to 0.90 mm.year−1 for 1993-2016, also depending on the dataset choice. We find that the temporal uncertainty dominates the budget, although the spatial-structural also has a significant contribution. On average, the intrinsic uncertainty plays a small part in the uncertainty budget. The main 15 source of uncertainty is the temporal uncertainty from the land water storage contribution, which is responsible for at least 50% of the total uncertainty, depending on the region of interest. The second main contributions come from the spatial-structural uncertainty from Antarctica and land water storage, which show that different locations of mass change can lead to trend deviations larger than 20%. As the barystatic SLC contribution and its uncertainty vary significantly from region to region, better insights into regional SLC are important for local management and adaptation planning. 20 Plain Language Summary The ice melt from Antarctica, Greenland and glaciers, and variations in land water storage cause sea-level changes. Here, we characterise the regional trends within these sea-level changes, taking into account mass variations since 1993. We take a 1 https://doi.org/10.5194/esd-2021-80 Preprint. Discussion started: 1 November 2021 c © Author(s) 2021. CC BY 4.0 License.

height, and bilinearly interpolated to a 1°by 1°grid. The SLE model then computes how the source mass change is redistributed over the oceans, taking into account the gravitational, deformation and rotational response of the Earth to these mass changes (Milne and Mitrovica, 1998;Mitrovica et al., 2001;Tamisiea and Mitrovica, 2011). The SLE model uses a pseudospectral approach (Mitrovica and Peltier, 1991) up to spherical harmonic degree and order 180 (equivalent to a spatial resolution of one degree). We assume a purely elastic solid-Earth response to the mass redistribution, based on the Preliminary Reference Earth 130 Model (Dziewonski and Anderson, 1981). The model computes relative SLC, which is the difference in height between the geoid and the solid Earth surface.

Trend and Uncertainty Assessment
As described above, we use a linear trend to compute the barystatic SLC, and we define three independent uncertainties (temporal, intrinsic and spatial-structural) to compute the uncertainty budget. Our trends and associated temporal uncertainty   (Table 1), purple the intermediate products, and green the final products. The yellow boxes indicate steps of the methodology, and the blue the main modules.
We use the following acronyms and abbreviations: OLS: ordinary least-squares; SLE; Sea-level equation; IC: Information Criteria; unc: uncertainty; NM: noise model; Hector: software package by Bos et al. (2013). can be selected to describe the autocorrelation between the residuals of the regression. The uncertainty of the regression model, representing one standard deviation, is then used as our temporal uncertainty.
Based on previous studies (Bos et al., 2013;Royston et al., 2018;Camargo et al., 2020), we test eight noise models to find 140 the best descriptor of the uncertainties in our data: white noise (WN), in which no autocorrelation between the residuals is considered; auto-regressive of orders 1, 5, and 9 (AR(1), AR(5), and AR(9), respectively), in which the order represents the number of previous observations influencing the next one; autoregressive fractionally integrated moving average of order 1 (ARF), which combines an AR(1) model with a fractional integration and a moving average of the noise; generalized Gauss-Markov (GGM), a generalized form of the ARF model.

150
The goodness of the fit of the models is assessed with the modified Bayesian Information Criterion (BIC tp ; He et al. (2019)), which is an intermediate criterion in relation to the Akaike (AIC; Akaike (1974)) and Bayesian (BIC; Schwarz (1978)) criteria. The best noise model is chosen by minimizing the criterion.
The second uncertainty we consider is the spatial-structural uncertainty (Figure 1b, right column). Studies that combine a large number of datasets often base the structural uncertainty of an estimate on the standard deviation over the individual 155 datasets in relation to the ensemble mean (Palmer et al., 2021;Cazenave et al., 2018). However, the small number of samples in our study (4 estimates for each contribution) could lead to unrealistic structural uncertainties when simply based on the standard deviation, as individual outliers could bias the ensemble mean. Instead, we compute the spatial-structural uncertainty by estimating the standard deviation based on the normalized fingerprint for each contribution. First, we use the trend of each contribution to compute sea-level fingerprints normalized to 1 mm.year −1 of global mean SLC. By doing so, we reduce the 160 effect of outliers on the standard deviation of the fingerprints, and only preserve the effect that the mass source distribution has on the fingerprint shape. We then take the standard deviation across the four datasets for each mass source contribution, which leads to four normalized spatial-structural uncertainties reflecting the uncertainty associated with the different spatial resolutions of the datasets. For example, the spatial-structural uncertainty for AIS, will reflect the differences in the fingerprints due to the fact that GRACE datasets provide observations in a 0.25 degrees resolution, while R19 provides mass changes can therefore not include intrinsic uncertainty in our uncertainty budget. Since this uncertainty represents systematic errors and instrumental noise, we assume no autocorrelation in the errors time-series, and propagate the errors through the ordinary least-squares (OLS) regression: Where y is our mass source change time series, Q yy is an identity matrix with the intrinsic errors on the diagonal, and A is the design matrix. Q xx is used to compute the standard error of the trend (Heij et al., 2004): We then use the standard error of the trend (σ trend ) as input in the SLE model, to see how the mass associated with the intrinsic uncertainty is distributed over the oceans.

Combining Trends and Uncertainties
To compute total barystatic SLC and its uncertainties, we sum the individual contributions (AIS, GIS, LWS and GLA) as Whereas the trends are added together linearly, we add the uncertainties in quadrature, assuming they are independent and normally distributed. For each contribution, we first combine the different types of uncertainty following Equation (3): where σ CON T R is the total uncertainty for each individual contribution (AIS, GIS, GLA, LWS). We then compute the total barystatic uncertainty for all contributions (σ total ) following Equation (4):

Results
In this Section we first present the noise model selection (Section 3.1) used to compute the barystatic SL trend and temporal 200 uncertainty in Section 3.2. We then present the spatial-structural (Section 3.3) and intrinsic uncertainties (Section 3.4). Lastly, we show the regional patterns of the total barystatic trends (i.e., the sum of the different contributions) and uncertainties (i.e., the sum of the different contributions and types of uncertainties) and zoom in on a few coastal examples (Section 3.5). There is a general preference for AR(1) in areas with smaller LWS changes (i.e., not the large drainage basins). On the other hand, over the large drainage basins, the same model preference mentioned above is maintained (Figure 2, right column). This suggests that GRACE observations and the hydrological models might not always be capturing the same processes.
Different noise models are selected as optimal for the two GRACE datasets: CSR datasets (  On the other hand, internal regions of the ice sheets, where there is little ablation, are better described by the PL model.

Trend and temporal uncertainty
The mass source trend and uncertainties obtained with the selected noise models (Section 3.1) are used to compute the sea-level fingerprints with the SLE model ( Figure 3). To illustrate the difference between the fingerprints based on GRACE and those based on GRACE-independent datasets, we show the trends and uncertainties for the JPL estimates (Figure 3a rotation-deformation patterns are visible in all fingerprints: regions closer to a freshwater source present a lower SLC, due to the mass loss that causes land uplift and reduced gravitational attraction, while in the far-field the sea level rises more than average. While all trends strongly depend on the dataset (Figure 3, first and third column), the uncertainty patterns are rather con-240 sistent. This suggests that, even though different noise models were used to compute the trend for each dataset, the temporal uncertainty is characteristic of each contribution. For glaciers and the ice sheets, the GRACE-independent datasets estimate a higher trend than the GRACE observations. The temporal uncertainties for ice sheets and glaciers are relatively small, especially for the UCI datasets. This indicates that these contributions do not exhibit strong autocorrelations, and as a consequence the uncertainty of the trend will be small. On the other hand, the temporal uncertainty for the LWS is larger than the trend itself, 245 and therefore the LWS trend is not statistically significant. This is probably related to the large internal and decadal variability of the time series, in combination with the relatively short period under study.
The largest inter-dataset differences are displayed in the regional patterns of the LWS contribution. Despite the similar global mean LWS trend value for both JPL and WGP, the regional trend patterns and uncertainty values are very different. This may partially be related to the coarse resolution of GRACE (300 km) in comparison to the hydrological models (0.5 • by 0.5 • grid 250 (55 km by 55 km at the Equator)). This difference can also be related to the difficulty in modelling the complex processes Another significant inter-dataset difference is in the regional trend pattern as a consequence of AIS mass change ( Figure   3a,e). This is mainly related to the location of ice mass changes in each dataset. GRACE observes mass accumulation in East 255 Antarctica, resulting in a positive sea-level trend in the region. This accumulation is not captured by the UCI data set. GRACE has a higher spatial resolution, and thus provides more detail of where the mass change is taking place. The UCI dataset provides estimates on a basin scale, so more detailed changes may be averaged out. The effect of the location of mass change at the source of the contribution is further investigated with the spatial-structural uncertainty (next section).  The regional SLC fingerprints directly reflect the differences in the spatial distribution of the mass change sources of the datasets . Over the ice sheets, for instance, IMBIE provides one time series for the entire Greenland Ice Sheet, which is subdivided into dynamic and surface mass balance changes, and the Antarctic Ice Sheet is divided into three To account for the uncertainties arising from the differences in location of the mass change between datasets, we first normalize 265 the fingerprints and then combine them into estimates of the spatial-structural uncertainty (Figure 4).
For all contributions, the largest spatial uncertainties are concentrated closer to the mass change sources, while the uncertainties are reduced in the far field. The effect of differences resulting from Earth rotational effects (typically leading to four large quadrants) is visible in the far field of the AIS (in the Northern Pacific) and of hotspots of LWS (around the Southern Ocean).
As was the case for the trends (Figure 3a), the AIS shows the strongest spatial differences, as the underlying datasets strongly 270 differ in their spatial detail. The spatial uncertainties represent the error introduced by using datasets that have insufficient resolution to solve the processes being analysed. In addition, it also shows that different physical processes are captured by the different datasets, as is the case for the LWS estimate. The LWS models have higher resolution than the GRACE observations, nonetheless the spatial-structural uncertainty of LWS component (Figure 3d) is the second largest. Black dashed contour line and number indicates the spatial average of the regional uncertainty. The final type of uncertainty considered here is the intrinsic uncertainty, which represents noise related to the dataset itself.
This type of uncertainty arises from the processing of the data, and needs to be provided with the model/observation. This information is only available for the JPL and IMBIE datasets ( Figure 5). All intrinsic uncertainties are fairly small (note that the colorbar ranges only up to 0.10 mm.year −1 ). The maximum uncertainty is seen in the LWS contribution (Figure 5a Overall, the intrinsic uncertainty, which is dependent of the observational technique, is relatively small when compared to the spatial-structural and temporal uncertainties, which are related to the physical processes represented.

Total Barystatic Trend and Uncertainty
Combining the different contributions, as explained in Section 2.2.3, leads to the total barystatic trends and uncertainties shown in Figure 6. Although we analysed six barystatic combinations, here we show only two (JPL and IMB+WGP) to discuss the patterns and the total uncertainty fields. We show these specific combinations because they present the most complete uncertainty budget (as only JPL and IMB had intrinsic uncertainties). Additional combinations are presented in Supplementary 290 Figure A2, with the global mean values listed in Supplementary Table A1. We recall that the aim of this study is not to provide one final ensemble of the barystatic contribution, but rather to focus on the uncertainty budget. Figure 6 shows For most regions of the world, we find that the regional barystatic sea-level trend is higher than the 1-sigma total uncertainty, with exception of the regions near the polar areas (indicated by stipples in Figure 6). Comparing the JPL trend to the IMB+WGP 300 trend, the shape of the pattern is similar, but the global mean (and thereby the regional SLC) is larger for JPL. This is also reflected by the histograms of the regional trend (depicted below the maps), which indicate larger regional SLC values for the JPL dataset (locally ranging from −0.43 to 2.19 mm.year −1 ) than for the IMB+WGP combination (locally ranging from −0.10 to 1.69 mm.year −1 ). The regional histograms also show a clearly skewed distribution of the trend, with mainly positive values. When we compare the two periods of IMB+WGP (Figure 6c, e), the regional histogram is slightly narrower for the 305 longer period (i.e., less divergence for the regional values), ranging from −0.39 to 1.50 mm.year −1 . This is probably because the local effect of internal variability plays a smaller role in the longer period. Nonetheless, the regional pattern is similar for both periods.
The uncertainty patterns ( Figure 6, right panels) are similar for the different dataset combinations (JPL vs. IMB+WGP) and periods (2003-2016 vs. 1993-2016). JPl regional uncertainties range from 0.62 to 0.98 mm.year −1 , while the IMB+WGP 310 combination ranges from 0.62 to 1.02 mm.year −1 , both for the 2003-2016 period. Just like for the trend, the longer period IMB+WGP uncertainties have a similar pattern but with lower values than for the shorter period, with regional values ranging from 0.37 to 0.75 mm.year −1 . Although the total uncertainty is dominated by the temporal uncertainty (see Figure 7), the similarity of the uncertainty pattern for both periods is influenced by the fact that the spatial-structural errors are based on the 2003-2016 period and extended to 1993-2016. On average, the spatial-structural uncertainty represents 16% (25%) of the total 315 uncertainty, while the temporal represents 80% (70%), for the 2003-2016 (1993-2016) period.
To further illustrate how the different contributions and uncertainties contribute to the total uncertainty budget, we selected ten coastal cities around the world in which we break down the total uncertainty of barystatic SLC from 1993-2016 into the four  contributions (Figure 7a), and into the three types of uncertainties (Figure 7b). We also show the different types of uncertainties for each of the contributions (Figure 7c). As in in Figure 6, we show the IMB+WGP combination.

320
The large contribution of the LWS and temporal uncertainty to the uncertainty budget is highlighted on Figure 7. Figure   7a shows that the LWS uncertainty plays an important role at all locations, being responsible for at least 50% of the total uncertainty. While the temporal uncertainty is the main contribution of the LWS uncertainty (Figure 7c), in some locations, such as Washington (US, location 3) and Tokyo (Japan, location 9) the spatial uncertainty is also important. Even without the contribution of LWS to the total uncertainty (Supplementary Figure A5), the temporal uncertainty is still the main contributor.

325
The second main contribution to the uncertainty budget comes from the AIS and glaciers (GLA). In these examples, the relative importance of AIS and GLA is fairly similar, with exception of Vancouver (Canada, location 1), for which the glaciers contribute about 5 times more than AIS. For the GLA, the dominant type of uncertainty strongly depends on the location ( Figure 7c). For example, in Vancouver (Canada, location 1), Lima (Peru, location 2) and Rotterdam (The Netherlands, location 5), the spatial-structural uncertainty dominates the contribution from GLA. In all other locations, the intrinsic and temporal 330 uncertainties play a more important role in the GLA contribution to the uncertainty budget. The AIS uncertainty is mainly dominated by the temporal uncertainty, with exception of Cape Town (South Africa, location 6), which is located within the large uncertainty bulge of the spatial-structural uncertainty from AIS (see Figure 4a).
The GIS and intrinsinc uncertainty play a small role in the uncertainty budget. On average, the GIS is only responsible for about 10% of the total uncertainty (panel a), and its uncertainty is generally dominated by the intrinsic and spatial uncertainties 335 rather than temporal uncertainties (panel c). The intrinsic uncertainty (panel b) is fairly small in all locations. Even for the JPL combination (Supplementary Figure A4), which has intrinsic uncertainty estimation for all contributions, the intrinsic uncertainty is responsible, on average, for only 5% of the total uncertainty. Table A1. Table with global mean barystatic sea-level changes and uncertainties from the original global mean timeseries. Note that these numbers may be different compared to the histograms of Figure 6, which represent the spatial average of the regional trend and uncertainty.
The difference between the trends is due to the use of noise-models for the regional trend, against an ordinary least-squares fit for the global mean trend.