Absolute Baltic Sea Level Trends in the Satellite Altimetry Era: A Revisit

1. Introduction

Coastal societies are forced to constantly adapt to changes in sea level (SL). Global SL products such as those produced by the European Space Agency's Sea Level Climate Change Initiative (SLCCI) (Legeais et al., 2018), the Integrated Multi-Mission Ocean Altimeter Data for Climate Research (Beckley et al., 2017), and the Copernicus services (Von Schuckmann et al., 2016), are proving to be instrumental in tracking the global SL rise, one of the most severe impacts of climate change. These synoptic and objective SL products are generated using the fleet of satellite-based altimeter sensors in orbit for over two decades. However, for regional to local coastal adaptation and planning for future scenarios, regionally-tailored SL information is required.

In the Baltic Sea (BS), satellite observations are particularly important given that the network of tide gauge (TG)s, which measures relative SL, is strongly affected by Vertical Land Motion (VLM) and, in particular, due to the Glacial Isostatic Adjustment (GIA) (Ludwigsen et al., 2020). For example, relative SL trends in the northern part of the BS over the last few decades have been shown to be strongly negative, while absolute SL trends display significant positive trends (Olivieri and Spada, 2016; Madsen et al., 2019). While global altimetry products have been used to study SL in the area (Karabil et al., 2018), they are affected by data gaps that are smoothed out by the typical strong interpolation in space and time (Madsen et al., 2019). In particular, the BS includes the two main features that limited the use of satellite altimetry since the start of the “altimetry era”: the presence of sea-ice and the proximity of the coast. For example, the average annual maximum extent of sea-ice in March covers up to 40% of the water surface (Leppäranta and Myrberg, 2009), and there are around 200,000 islands in almost 400,000 km² of water surface. However, an advantage for using altimetry in the BS is that the tidal component is limited, which mitigates the known problems of tidal modelling in areas with complex coast and bathymetry.

The advances in altimetry processing have shown that dedicated signal processing techniques are able to enhance the quality and the quantity of the retrievals (Benveniste et al., 2019). This is particularly significant when a better fitting of the radar return (retracking) is combined with a dedicated selection of corrections to the altimetric range to enhance coastal data (Benveniste et al., 2020). In the sea-ice covered areas, classification techniques are able to identify water apertures (leads), which when combined with retracking allow for the retrieval of SL. Such processing has recently driven an improved SL analysis in the Arctic Ocean (Rose et al., 2019), but has never been applied to the BS yet. This gap presented an opportunity for the European Space Agency's Baltic+ Sea Level (ESA Baltic SEAL) Project, to produce a regional gridded SL product that incorporates observations from altimetry measurements acquired from sea-ice leads, and coastal waters.

The use of dedicated altimetry SL products, combined with external datasets, can contribute to characterise local differences in SL trend and variability (Passaro et al., 2016). The SL variability in the BS has been found to be highly correlated with the variability of westerly winds (Andersson, 2002). Wind patterns are modulated by large scale variability of the atmospheric pressure, which can be described by climate indices, and in particular by the North Atlantic Oscillation (NAO) index (Jevrejeva et al., 2005). The objective of this study, therefore, is to analyse the SL trend in the BS during the altimetry era and to characterise the relationship between wind patterns, NAO and variability in absolute SL at a sub-basin scale. To do so, we use the ESA Baltic SEAL monthly gridded SL, summarising the methodology followed to create this dataset and assessing its performances against in-situ data and global gridded SL products. Subsequently, trends are computed providing statistical uncertainty that takes into account the serial correlation in the time series.

2. Methods

2.1. Altimetry Data Processing

This study is based on the analysis of the ESA Baltic SEAL project, whose documentation is freely available from http://balticseal.eu/outputs/. In the context of this project, along-track data from most of the altimetry missions operating in the last two decades are reprocessed to generate a dedicated monthly gridded product from May 1995 to May 2019. In particular the following conventional Low Resolution Mode (LRM) altimetry missions are used: TOPEX-Poseidon (TP, from May 1995), Jason-1 (J-1), Jason-2 (J-2), ERS-2, Envisat, and SARAL. Considering the latest Delay-Doppler (DD) altimetry technology, data from the following missions were acquired: Cryosat-2 (CS-2), Sentinel-3A, and Sentinel-3B (S3-A/B).

High-frequency data are downloaded, i.e., distributed at 20-Hz rate for most of the missions, except 10-Hz for TP, 40-Hz for SARAL, and 18-Hz for Envisat. A detailed list of the data source and the version of the altimetry data acquired for the reprocessing is provided in Ringgaard et al. (2020).

The overall process from radar pulse to SSH estimate delivered as a gridded product is outlined in Figure 1. Also outlined are the enhancements which the ESA Baltic SEAL Project implements to tailor the data produced to BS regional stakeholders, and further develop best practice for coastal altimetry globally. These enhancements are summarised in the following sections and build on the overall pulse-to-SL process.

FIGURE 1

Figure 1. From measurement to corrected sea level estimates delivered as along-track and gridded datasets. A simple flowchart of the origins of the waveform data, and how it is processed to produce the ESA Baltic SEAL estimate of sea level. Process steps which are enhanced and tailored for use in BS, are highlighted in orange.

2.1.1. Classification

In the winter months, parts of the BS, especially the regions Bay of Bothnia and Gulf of Finland as defined in Figure 2, are covered by a dynamic changing sea-ice cover, which makes continuous, gapless SL estimations difficult. Moreover, observations during the winter season are limited to leads (i.e. narrow cracks within the sea-ice) enabling only brief, spatially limited opportunities for open water measurements.

FIGURE 2

Figure 2. The area of study divided in different boxes representing major sub-basins. Color shades display the bathymetry.

The behaviour of the reflected radar signals is used to classify lead returns and to identify open water areas within the sea-ice. The applied open water detection is based on unsupervised, artificial intelligence machine-learning algorithms. The algorithm used in this study has been described in Müller et al. (2017) for LRM altimetry and its versatility to be applied also to DD missions has been shown in Dettmering et al. (2018).

In order to group the reference datasets automatically into a specific number of clusters representing different waveforms types, the waveforms samples from a training set are introduced to a K-medoids clustering algorithm (Xu and Wunsch, 2008; Celebi, 2014). In principle K-medoids searches for hidden similarities within the data based on a given input feature space by minimizing the distance between the individual features and most centrally positioned features (medoids) from the feature space itself. At first the algorithm defines randomly K-medoids followed by the computation of the distance between all features and the initially selected K centres. In the next steps K-medoids rearranges the location of the medoids as long as there is no motion among the features. K-medoids belongs to the “partitional” clustering algorithms, which require a pre-defined number of clusters K.

After clustering, the clusters are assigned manually to the different surface types, by using background knowledge about the physical backscattering properties of the individual surface conditions and feature statistic per each cluster.

The waveform features considered are the usual ones describing the shape of the echo: waveform maximum, trailing edge decline (an exponential function fitted to the trailing edge of the waveform), waveform noise, leading edge slope, trailing edge slope. The features are applied to all required satellite missions and altimetry datasets.

The second part of the classification is related to the classification of the remaining waveforms. Therefore, a K-Nearest Neighbour (KNN) classifier is applied. KNN searches for the closest distance between the reference model and the remaining waveforms (Hastie et al., 2009). The majority of clusters among the K nearest neighbours defines which class has to be assigned to the waveform.

As a final result of the classification, each high-frequency waveform is classified as open water or sea-ice. Waveforms classified as sea-ice returns are flagged and not considered for the generation of the gridded product.

2.1.2. Retracking

All the waveforms from the altimetry missions used in this study, except for TP, are re-fitted by means of dedicated algorithms (retrackers) that are able to track the signal in open ocean, coastal, and sea-ice covered conditions.

ALES+ is the retracker that has been further developed and extended to all the missions considered in this project. ALES+ is based on the Brown-Hayne functional form that models the radar returns from the ocean to the satellite. The Brown-Hayne theoretical ocean model is the standard model for the open ocean retrackers and describes the average return power of a rough scattering surface (i.e., what we simply call waveform) (Brown, 1977; Hayne, 1980). A full description of ALES+ retracker for LRM missions is provided in Passaro et al. (2018a).

In the case of the DD waveforms, the correspondent version called ALES+ SAR adopts a simplified version of the Brown-Hayne functional form as an empirical retracker to track the leading edge of the waveform. The model simplification is achieved by assigning a fixed decay of the trailing edge, instead of a dependency with respect to antenna parameters (beamwidth, mispointing) as in the LRM case. A full description of ALES+ SAR retracker is provided in Passaro et al. (2020a). Data from CS-2, S3-A/B can be reprocessed with ALES+ SAR using the ESA Grid Processing On Demand (GPOD) service (https://gpod.eo.esa.int/, Passaro et al., 2020b).

By means of a leading edge detection, ALES+ and ALES+ SAR retrack only a subwaveform whose width is dependent on the wave height in LRM and fixed in DD. In this way, it is possible to avoid considering signal perturbations typical of the coastal zone when fitting the echoes. In the case of peaky waveforms, typical of leads within sea-ice, both algorithms perform a direct estimation of the trailing edge slope.

Together with the retracking, a dedicated sea state bias (SSB) correction is computed. The SSB correction is computed at 20-Hz rate. This guarantees a better precision of the range retrieval, since it decreases the impact of correlated errors in the retracked parameters (Passaro et al., 2018b).

For the LRM missions, the SSB correction is derived by using the same 2D map from Tran et al. (2010), but computed for each high-frequency point using the high-frequency wind speed and significant wave height (SWH) estimations from ALES+.

In the original DD altimetry products, the SSB correction is either missing (Cryosat-2) or computed using the Jason model. Here instead, a first model is computed specifically for the ALES+ SAR retracker. As a reference parameter on which the model is built, we take the rising time of the leading edge, which is taken as a proxy for the significant wave height.

The corrections are derived by observing the SL residuals (with no correction applied) at the points of intersections between satellite tracks (crossover points). A wider region covering the North Sea and the Mediterranean Sea is used in order to have more open ocean crossover points, which are scarce in the BS. The residuals are modelled with respect to the variables influencing the sea state (here the rising time of the leading edge) in a parametric formulation.

In particular, at each crossover m:

$\begin{array}{ll}\Delta S{L}_m=\widehat{\alpha }{\sigma }_{c_o}-\widehat{\alpha }{\sigma }_{c_e}+ϵ& (1)\end{array}$

where o and e stand for odd and even tracks (indicating ascending and descending tracks respectively), ϵ accounts for residual errors, σ_c is the rising time of the leading edge. We, therefore, have a set of m linear equations, which is solved in a least square sense. The chosen α is the one that maximises the variance explained at the crossovers, i.e., the difference between the variance of the crossover difference before and after correcting for the SSB using the computed model.

2.1.3. Choice of Range Corrections

Once the ranges have been obtained by retracking the waveforms, the following altimeter equation is implemented to derive the sea surface height (SSH):

$\begin{array}{l}SSH=H_{orbit}-(R+WT+DT+IONO+SSB+DAC\\ +SET+PT+ROC)& (2)\end{array}$

The atmospheric and geophysical corrections applied are listed in Table 1. H_orbit and R stands for the orbital height above the TOPEX/POSEIDON ellipsoid and the retracked range between the satellite and the sea surface. To generate the gridded product, the SSH is also corrected for ocean tide and load tide using the FES2014 tidal model (Carrere et al., 2015).

TABLE 1

Table 1. Altimeter corrections/parameter applied to the along-track processing.

2.1.4. Multi-Mission Cross-Calibration

In order to ensure a consistent combination of all different altimetry missions available, a cross-calibration is necessary. We follow the global multi-mission crossover analysis (MMXO) approach described by Bosch et al. (2014) in order to produce a harmonized dataset and a consistent vertical reference for all altimetry missions.

For all crossover locations, a radial correction for the observations of both intersecting tracks is estimated by a least squares approach based on crossover differences without the application of any analytic error model. These corrections are later interpolated to all measurement points of all missions included in the analysis.

The approach was developed for global calibration and is adapted for regional applications within ESA Baltic SEAL. This comprises the following points of change in comparison to (Bosch et al., 2014):

1. The maximum acceptable time difference for the crossover computations is increased from 2 to 3 days, in order to ensure enough crossover differences in the BS region.

2. For the same reason, all crossover points are used, including coastal areas.

3. For the computation of crossover differences, high frequency data are used. This is realised by changing the interpolation of along-track heights to crossover locations from point-wise to distance-wise.

4. All missions are equally weighted. The weighting between crossover differences and consecutive differences is adapted in order to account for the smaller region.

2.1.5. Gridding

After the multi-mission cross-calibration, the along-track SL estimates undergo an outlier detection whose consecutive steps are listed in Passaro et al. (2020a). After this flagging, the observations are interpolated on an unstructured triangular grid (i.e., geodesic polyhedron). The grid has a spatial resolution of 6–7 km.

The gridded monthly SL estimates are obtained by fitting an inclined plane to each grid node by means of weighted least square interpolation, considering SSH along-track information within 100 km radius around the grid node centre. Distance-based Gaussian weights are defined in a diagonal matrix W. The median absolute deviation of the along-track SSH observation within an area in the open sea without complex topographic features is computed per mission and used as an estimation of the uncertainty. This is placed as variance on the main diagonal of the uncertainty matrix Q_bb. The uncertainty information and W are combined to form the least-squares weighting matrix P_bb, following the equation:

$\begin{array}{ll}{P}_{bb}=W*(\frac{1}{Q_{bb}})& (3)\end{array}$

In order to eliminate still existing outliers among the along-track data within the cap-size, the weighted least square estimation is performed iteratively. At each iteration, the difference between the monthly average and all along-track SSH values is evaluated. The along-track estimates whose residual exceeds 3 times the standard deviation of all residuals (3-sigma criterion) are flagged as outliers and the weighting matrix P_bb is consequently updated, before a new least-squares adjustment is performed.

Finally, an additional outlier rejection based on a Student distribution is performed: the standardised residuals of each remaining observations within the search radius are tested against quantiles of the Student distribution (t), setting the 99th percentile as boundary condition. If a standardised residual is smaller than the value of the t distribution, the corresponding observation is used in a last iteration of the least-squares adjustment.

2.2. Validation of the Altimetry Product

To validate the gridded SL product from altimetry, we perform a comparison against SL data from TGs. The main source of data is the Copernicus Marine Environment Monitoring Service (CMEMS) service and some data are complemented from the national datasets of the Danish Meteorological Institute (DMI), the Finnish Meteorological Institute (FMI) and the Swedish Meteorological and Hydrological Institute (SMHI). A full list of TG data sources used in this study is available in Ringgaard et al. (2020).

To avoid gaps in the time series, we consider only grid points with at least 250 months of valid data. We also divide the BS in different sub-basins whose naming and geographical extensions are provided in Figure 2.

The TG and altimeter SL measurements are not equivalent and hence both data sets were further processed before they were compared. In particular, to allow the comparison, the DAC was added back to the altimetry data, since TG data is not corrected for it.

The SL reference frame of the altimeter SL height is tied to the TOPEX ellipsoid, while TG SL height data are referred to the Normaal Amsterdams Peil (NAP) reference frame. To allow for comparison of the TG-altimetry pairs, the mean of the gridded SL was removed and set equal to the mean of the corresponding TG.

In order to validate the gridded dataset, the grid points within 20 km from every TG are considered. As the gridded dataset has a frequency of a month, the TG data are monthly averaged. TGs measure relative SL and altimeters measure absolute SL, hence the effect of land uplift is removed from the TG data using the GIA model NKG2016LU (Vestøl et al., 2019). NKG2016LU is the latest GIA model for the BS. The closest absolute land uplift values are located for each TG and used for trend removal.

The root mean square error (RMSE) and the Pearson correlation coefficient (r) are computed for all TG-grid pair time series and the results are displayed in Figure 3. Out of 67 TGs and gridded altimetry pairs, 62 show a correlation higher than 0.6 and 61 have a RMSE lower than 9 cm. The lowest performing area is located north of the Danish straits. A possible reason lies in the performances of the tide corrections, which are much more important north of the straits, than in the BS.

FIGURE 3

Figure 3. Correlation coefficient and root mean square error between every TG considered and the altimetry grid points used for trend computation which are located within 20 km. The circles showing the statistics are co-located with the TGs.

The results of the validation considering the sub-basin division used in this study are summarised in Table 2. For this purpose, TG and altimetry gridded data are averaged in space across each sub-basin and in time every 3 months. The comparison shows that the correlation is never lower than 0.75 and the RMSE is never higher than 0.10 m.

TABLE 2

Table 2. Pearson Correlation Coefficient and RMSE between TGs and gridded altimetry averaged in space over the sub-basins considered in this study and in time every 3 months.

2.3. Methods for SL Analysis

2.3.1. Trend Computation

We estimate the seasonal cycle, the linear trend and the parameter uncertainties by fitting multi-year monthly averages (to approximate the seasonality) and a linear trend to the monthly gridded data. In terms of formulation, this means fitting the time series d(t) with the model y(t), which for every monthly step t is defined as:

$\begin{array}{ll}y(t)=o+at+m_i+ϵ& (4)\end{array}$

Where o is an offset term, a is the linear trend, m_i is the multiyear monthly mean for the month i corresponding to the time step t and ϵ is the residual noise. The trend estimate is found solving the fitting by linear least squares. The standard error σ of the least square solution would be nevertheless unrealistic, since it would not consider the autocorrelation of the time series. Therefore, to account for the autocorrelation, σ is found by an iterative maximum likelihood estimation (MLE), as described in [6]. This requires the definition of an appropriate covariance matrix of the observations, including a formulation of the residual noise. In particular, we investigate the fit of a variety of different stochastic noise model combinations as done in e.g., Royston et al. (2018): These are an autoregressive AR(1) noise model, a power law plus white, a generalized Gauss Markov (GGM) plus white, a Flicker noise plus white and an auto-regressive fractionally-integrated moving-average (ARFIMA) model. For the considered domain we find that on average the AR(1) has the lowest mean (or median) values of the Akaike Information Criterion (AIC, Akaike, 1998) and the Bayesian Information Criterion (BIC, Schwarz, 1978). Finally, the uncertainties provided in this study are scaled as 1.96 * σ to obtain a 95% confidence interval.

2.3.2. Principal Component Analysis

We perform a Principal Component Analysis (PCA) (Preisendorfer, 1988) to investigate the major modes of SL variability in the BS. For this purpose, we consider a set of multiple SL anomalies x_k(t), where k and t describe the dimensionality of the data in space and time, respectively. To identify the maximum modes of joint space and time variations, we determine a set of linear combinations in form of Principal Component (PC) u_m(t) and associated eigenvectors or Empirical Orthogonal Function (EOF) e_km. The linear combinations, or the modes are arranged such that the higher-order modes m = 1, 2, 3, … explain the highest variance fractions of the data. The PCs u_m(t) are equal to the projection of the data vector onto the m^th eigenvector e_km (e.g., Wilks, 2006):

$\begin{array}{ll}{u}_m(t)={\displaystyle \sum _{k=1}^K}{e}_{km}{x}_k(t),m=1,\dots ,M& (5)\end{array}$

In this manner the data is explained by a set of PCs, which represent time series (which are uncorrelated or independent from each other), as well as the EOFs (or eigenvectors) which represent the geographical coherence of the individual modes.

We compute the EOF and their PC from monthly gridded deseasoned SL. The latter is called sea level anomaly (SLA) in this description, since it is the anomaly with respect to a monthly-based average (for example, based on the average for all Januaries in the period of record at a particular grid point). In this way we capture the “full-year” monthly variability and no seasonal variations. Because monthly SL variability is generally most pronounced in winter, the derived full year EOF-pattern are very similar to the ones derived only over the winter season (DJF). EOF patterns are given as point-wise correlations of their PCs with SLAs.

2.3.3. Regression Analysis

We use a simple statistical approach to understand the relation of surface winds and SL trends: We compute point-wise linear regressions of the deseasoned, monthly and basin-averaged surface winds (zonal U component and meridional V component) and SLAs by solving for: SLA(t) = aU(t) + bV(t) + η, where a and b are the first order partial regression coefficients to be estimated and η is the residual (e.g., Storch and Zwiers, 1999; Dangendorf et al., 2013). Based on these point-wise linear regressions we estimate a linear trend (without seasonal component) which is explained by the individual wind components as well as the explained variance of SL variability by the components (as for example in Dangendorf et al., 2013).

3. Results and Discussion

3.1. Absolute Sea Level Trends

Figure 4 shows the map of SL trends estimated using the ESA Baltic SEAL dataset. Superimposed in circles along the coast are the estimations of the TGs, which are corrected for GIA. In accordance to previous studies based on the altimetry era (e.g., Madsen et al., 2019), it is found that the absolute SL has been rising throughout the region. The rate of SL rise increases from the South West of the Baltic Sea (S-W) to the Gotland Basin, and from the Gotland Basin to the Bay of Bothnia and the Gulf of Finland.

FIGURE 4

Figure 4. SL trends (A) and corresponding uncertainties (B) estimated by altimetry (shading) and TGs (circles) from May 1995 to May 2019. TGs are corrected for the GIA using the NKG2016 model. Uncertainties are reported as 95% confidence interval.

By retrieving the SL information from the leads within sea-ice in winter, we are able to extend the analysis to areas characterised by seasonal sea-ice coverage, i.e., Bay of Bothnia and Gulf of Finland. Nevertheless, gaps are still present in such regions along the coast. Other remaining gaps involve locations in the Danish Archipelago, where the predominant presence of land within the search radius of each grid point hinders the possibility to find enough data in particular during years in which few altimeters were in orbit. Finally, our SL analysis does not provide results in some parts of the Turku Archipelago (south-western Finnish coast). The presence of numerous islets in this area means that the vast majority of SL retrieval are located at distances below 1 km from the nearest land. This is well below the possibilities of any LRM altimeter, even using coastal retracking to avoid land contamination.

These data gaps could be artificially mitigated by means of heavier interpolation and different weighting in the gridding process, nevertheless the choice in this study is to avoid generating information that is indeed not available. The comparison of the agreement between the SL trends from different altimetry dataset and TGs presented in Figure 5 is a proof of the validity of our solution. In the histograms, the SL trend estimates from the TGs are compared with the closest estimates from altimetry using data from this study (Figure 5A) and data from CMEMS (Figure 5B, Taburet et al., 2019). The time series are restricted to the interval May 1995-December 2018 to enable the comparison with CMEMS. In Figure 5C, the length of the time series of this study is May 1995–December 2015, to enable the comparison with the gridded product of the SLCCI (Figure 5D, Legeais et al., 2018). In both pairs of comparison, the comparability between trends from altimetry and from TGs improves by 9% using the Baltic+ data in terms of root mean square of the differences. All the altimetry dataset show a median of the trends that is about 0.2 mm/yr lower than in the TG records. We acknowledge that the nonlinear elastic uplift from present day deglaciation, which is not taken in consideration by the GIA model, may affect the bias, although GIA has been shown to be the dominating source of vertical deformation in the region (Ludwigsen et al., 2020).

FIGURE 5

Figure 5. Histograms of the differences in estimated SL trends from gridded altimetry (SAT) and TGs, compared using the closest point. Each panel correspond to different SAT dataset: the altimetry dataset from May 1995 to December 2018 developed in this study (A, Baltic+), the altimetry dataset from May 1995 to December 2018 CMEMS (B, Copernicus), the altimetry dataset from May 1995 to December 2015 developed in this study (C, Baltic+), the altimetry dataset from May 1995 to December 2015 of the SLCCI (D, SLCCI).

Figure 4B also shows the uncertainty of the computed trends. This is a purely statistical function of the number of samples in the time series and their SL variability, taking in to consideration the serial correlation. The same method is used also to estimate the uncertainty of the trend estimated from TGs, which also do not have an uncertainty value associated to each measurement. This is in line with most of the studies estimating trends from altimetry measurements, for example (Benveniste et al., 2020). The possibility to associate an uncertainty to the single altimetry measurement has been explored by (Ablain et al., 2016) and analysed in the BS by Madsen et al. (2019), but requires a large amount of assumptions concerning every single correction added to the altimetric range. Nevertheless, our statistical uncertainties show a similar pattern and range of the ones shown in Madsen et al. (2019).

By grouping the grid points according to their location, Figure 6 displays the averaged SL trends of each sub-basin with their statistical uncertainty. In Figures 6B,C, the monthly time series for the Bay of Bothnia and the S-W are shown as examples, since they present the largest discrepancies in the linear trend estimations. The rise in SL is statistically significant in all sub-basins. The spatial variation of the best estimate of the linear trend is confirmed, although the uncertainties due to the larger variability of the SL time series in most of the sub-basins cannot ensure statistical significance to this assertion yet.

FIGURE 6

Figure 6. (A) SL trends from gridded altimetry averaged across different sub-basins of the BS from May 1995 to May 2019. Corresponding uncertainties are reported as black error bars. (B,C) Monthly SL time series of the S-W Baltic Sea and the Bay of Bothnia. The linear trend is shown as a red line, with the shading representing its uncertainty.

Figure 7 shows the trends in SL considering only the winter months (Figures 7A,C) and only the summer months (Figures 7B,D). Positive trends are found in the winter SL, with a difference of over 4 mm/year comparing the winter trends in their minimum and maximum values. A similar gradient in SL trend estimates is seen for the full time series in Figure 4A. This spatial variation in the trend is less pronounced in summer. Due to the relatively short duration of the time series, the seasonal trend uncertainties are comparatively large. Thus, we investigate whether a similar pattern can also be found in TG records, which spans much longer periods than the altimetry time series.

FIGURE 7

Figure 7. SL trends with uncertainties from gridded altimetry computed using only the winter months (December, January, and February, A,C) and the summer months (June, July, and August, B,D).

Figures 8A–D report the best estimate of the linear trends in SL from the longest TG time series in the region, spanning from 1920 to 2020 at intervals of 25 years. Indeed, the same gradient of about 4 mm/yr in SL trend estimates from altimetry across the basin is observed in the most recent TG record (Figure 8A). It is observed not only that the SL rise in the BS is evident in the last 50 years, but also that the spatial gradient in trends has been increasing in time.

FIGURE 8

Figure 8. SL trends from TGs computed at 25-year intervals from 1920 to 2020. (A) Linear trends 2020–1995, (B) linear trends 1995–1970, (C) linear trends 1970–1945, (D) linear trends 1945–1920.

In the next section, we analyse the possible role of wind patterns and the NAO in shaping this spatial gradient.

3.2. Discussion

3.2.1. Relationships With Wind Pattern

To analyse the spatial and temporal pattern of SL variability and how they differ locally across the basin, we perform an EOF analysis on the deseasoned altimetry time series at each grid point (as described in section 2.3.2). Figures 9A,D show the spatial patterns of the first and second EOFs. We find that 87.4% of the variance in the entire domain is explained by the first EOF, which is associated with a uniform SL pattern across the basin. The second EOF, while representing 3.1% of the variance, is connected to a SL variability with a strong gradient from S-W toward Bay of Bothnia and Gulf of Finland, generating SL anomalies of opposite sign.

FIGURE 9

Figure 9. (A,D) Empirical orthogonal functions of SL variability expressed as correlation of their corresponding PC against the SL time series at each grid point. (B,C) Correlation between the first PC of SL variability and the zonal (U) and meridional (V) component of the surface winds. (E,F) Same for the second PC.

To characterise these two modes, we correlate the accompanying PC to the zonal (U) and meridional (V) components of the surface wind from the ERA5 reanalyses (Hersbach et al., 2020). The results displayed in Figures 9B,C,E,F, show that the PC1 is correlated with U in the South of the basin, with the correlation degrading toward North. The predominance of the zonal component in shaping the SL variance of the region is in accordance with previous studies, for example Johansson et al. (2014). PC2 is instead well described by the variability of V, with correlation values over 0.5 in all our region of study.

To further study how the wind variability may affect the estimates of SL trend during the altimetry era, we perform a multiple regression analysis of the SL time series using the U and V wind components (as described in section 2.3.3). We consider the large scale wind field by spatially averaging the monthly wind speed over the entire domain. The results are shown in Figure 10, in which the explained trend for each component of the wind and for the sum of the two components is presented, with its uncertainty. The average variance explained is 31% by U and 2% by V, but the latter explains over 15% of the variance in the Bay of Bothnia (not shown). Despite the high variance explained, the U regression shows a very small, homogenous trend in the whole BS, while V is responsible for a gradient of over 1 mm/yr from South to North. Although a conclusive statement with statistical relevance cannot be drawn, given the uncertainties, both EOF and regression analysis point out to the same role of the meridional wind component to shape a North-South imbalance in the SL anomalies. Recently, spatial gradients of the SL trend within the BS have been attributed, based on circulation models, to an increase in the days of westerly winds, which increase transport toward the east (Gräwe et al., 2019). Our analysis suggests that the meridional wind component also affects the differences in SL trend among the sub-basins.

FIGURE 10

Figure 10. Trends resulting from the regression of the zonal (U, A) and meridional (V, B) component of the surface winds on the SL time series. (C) Shows the trend obtained summing the two components of the regression. (D–F) Show the corresponding uncertainties of the respective trend estimations in (A–C).

3.2.2. Relationships With North-Atlantic Oscillation

The large scale variability of both SL and wind patterns in our area of study can be well-described using teleconnections. There is a strong agreement that the NAO is the leading mode of atmospheric circulation in the region (Andersson, 2002; Jevrejeva et al., 2005). Interconnections with other climate patterns and the corresponding indices have been shown to play a role in the area, such as the East Atlantic (EAP) pattern, Scandinavian (SCAN) pattern (Chafik et al., 2017) and the BS and North Sea Oscillation index (BANOS) (Karabil et al., 2018).

We focus on the local effects of NAO variability, since the relationship of NAO to the Baltic SL variability has been previously reported to be spatially heterogeneous (Jevrejeva et al., 2005 with TG observations, Stephenson et al., 2006 with global models). Figure 11A shows that in the altimetry era the correlation between the SL variability and the NAO index is dominant in winter, as expected, and uniform in the whole domain except for the S-W.

FIGURE 11

Figure 11. (A) Correlation of the NAO index with SLA from gridded altimetry. (B) Normalized time series of NAO index (green) and SLA difference between Bay of Bothnia and S-W sub-basins (orange). Each point represents the time average of the quantities of the winter months December, January, and February.

The possibility given by our dataset to observe local SL changes in winter in the sea-ice covered areas allow a basin-wide comparison of the Bay of Bothnia against the S-W, which present the largest discrepancies in the linear trend estimations. In Figure 11B, the difference in SL (Bay of Bothnia–S-W) in the winter months is plotted against the NAO index. The correlation between the two curves is 0.71. From the comparison it is seen that positive NAO phases are related to winters in which the SLA are higher in Bay of Bothnia than in S-W. As seen in the next section, this is linked to the action of stronger southerlies and westerlies winds during positive NAO phases, which push the water north and east of the basin through Ekman transport. The intensity of the NAO phase, which is linked to the wind forcing (Dangendorf et al., 2013), is here shown to drive differences of SLA at a sub-basin scale in the BS, with interannual variations that have an effect on the linear trend of the SL estimated on time series spanning two decades of observations. In particular, as seen in Figure 10, the effect of the positive NAO phases in our period of observation results in a wind-related SL trend increasing toward north.

3.2.3. Ekman Currents

Winds affect the surface circulation of water masses through Ekman transport. We show in Figure 12 the average winter wind speed direction (Figures 12E–H) and the resultant Ekman currents (Figures 12A–D) in selected winter seasons to observe the mechanism that may regulate the SL differences among different sub-basins within the BS. The years are chosen based on the highest and lowest differences between the winter SLA of Bay of Bothnia and the S-W as reported in Figure 11B. We analyse the Ekman transport at 15 m depth from the GlobCurrent ocean product (Rio et al., 2014). The Ekman currents are distributed on a 1/4 of a degree grid and they are derived using wind stress from ECMWF, Argo floats and in-situ surface drifter data.

FIGURE 12

Figure 12. Velocity of the Ekman currents (A–D) and surface wind velocity (E–H) during the winters of four selected years. The arrows represent the direction of the velocity vectors and are scaled according to their magnitude.

When considering these results, we are dropping the hypothesis of a fully developed Ekman spiral, in which case the transport would be perpendicular to the wind direction. Nevertheless, given the low depths of the S-W of the BS, the 15 m-depth Ekman transport should be a good approximation at least in this sub-basin. Since the Bay of Bothnia is covered by sea-ice for most of the winters, which hinders formation of an Ekman spiral, and since sea-ice is not taken into account in the GlobCurrent product, we are mostly interested in the effect of the Ekman currents in the southern part of the domain. The results are consistent with the Ekman transport pushing the surface waters to the right of the wind direction.

Winters with SLA higher in Bay of Bothnia than in S-W (e.g., 2000, 2014) are either characterised by strong westerlies in S-W whose intensity decrease toward the North (e.g., 2000), or by a marked southerly component of the wind (e.g., 2014). Years in which the differences are very low, or even flip (e.g., 1996, 2010) are characterised by much lower wind speed.

In conclusion, the years in which the winter SLA in the Bay of Bothnia are higher than in the S-W are characterised by a strong Ekman transport, which affects the sea-ice free part of the domain. This mechanism fails during the negative NAO phases, directly affecting the SL difference between the two sub-basins.

4. Concluding Remarks

This study analysed the SL trend in the BS during most of the altimetry era (1995–2019). A new reprocessing in the framework of the ESA Baltic SEAL project enables the retrieval of more data in the coastal zone and among sea-ice. A trend analysis based on this dataset improves the agreement with trends estimated using GIA-corrected TGs (Figure 5). The information retrieved from the leads among sea-ice covered areas enhances the possibility to study SL variability and its differences across the basin during winter, which is the season with the largest SL rise in our observation period (Figure 7).

The absolute SL rise is statistically significant in the entire domain, since the uncertainties are lower than the trend estimates (Figure 4). Differences in trends among the sub-basins are not statistically significant, but are seen in both the TGs and the altimetry dataset (Figure 6).

The absolute SL rise is a year-round phenomenon, although trends are higher in winter than in summer. The gradient in SL rise across the basin mainly occurs during winter (Figure 7). SL differences between the North and the South West of the BS are shown to be well-correlated with the NAO index in winter (Figure 11). In particular, winter positive NAO phases trigger lower SL anomalies in the S-W, as strong south-westerly winds transport surface water away from the sub-basin (Figure 12).

The NAO drives not only the SL in the entire domain, but is shown to also affect internal sub-basin gradients. A part of it can be explained by wind forcing, which accounts on average for about 40% of the SL variability. Other factors can contribute to the observed spatial gradient of the SL trend, which we plan to consider in a future study. Karabil et al. (2018) for example observes that a possible driver can be the freshwater flux, but this would be particularly pronounced in summertime, therefore would not explain the larger trend differences found in winter. The increasing use of GRACE data to compute mass SL changes at a regional scale (Kusche et al., 2016) and the availability of sea surface temperature and salinity datasets can be combined with measurements from the Argo floats (Guinehut et al., 2004; Boutin et al., 2013) (particularly in a low-depth basin like the BS). This suggests that a regional SL budget based on observational data shall be the subject of a future study and the next step to increase our knowledge of the Baltic SL variability and drivers.

This study highlights the value of developing regionalised SL products, using satellite altimetry measurements. It has improved the efficacy of retrieving meaningful SL observations from areas featuring complex coastlines, and those affected by sea-ice contamination of the altimeter footprint. While current efforts in the exploitation of altimetry in the coastal zone are focused on the analysis of along-track data, this work for the first time employs a coastal-dedicated reprocessing to produce gridded sea level data. Moreover, we have demonstrated that such techniques are able to obtain reliable sea level time series also in areas and seasons interested by sea-ice coverage. The BS has proved to be an excellent region to explore these issues concerning coastal altimetry. Using the best practice advances developed here, comparative analyses can be conducted in more tide-prone regions to test the further applicability of our approach.

Data Availability Statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found at: http://balticseal.eu/data-access/.

Author Contributions

MP designed the study, wrote the manuscript, and was the principal investigator of the research group. He was also the author of the retracking algorithm and of the sea state bias correction developed for this study. FM was responsible for the classification and the gridding. He was the main responsible for the altimetry database organisation, with CS and DD. FM also contributed in the coordination of the activities of the research group. JO was responsible for the trend and sea level variability estimation and interpretation. DD was responsible for the multi-mission calibration. AA and OA were responsible for the mean sea surface used to obtain the sea level anomalies. MH-D contributed in the interpretation of the sea level trends. FS provided the basic resources making the study possible and coordinates the activities of the research group at TUM. RS helped in the coordination of the research group and reviewed the manuscript. JH, KM, and IR contributed to the interpretation of the results. LR, JS, and LT were responsible for the validation of the along-track data. MP, MR, and JB conceptualized the research project. MR and JB reviewed the manuscript and all deliverables of the research project and supported the validation activities of the ALES+ SAR retracking algorithm. All authors read, commented, and reviewed the final manuscript.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

This study is a contribution to the ESA Baltic+ Sea Level project (ESA AO/1-9172/17/I-BG - BALTIC +, contract number 4000126590/19/I/BG). We use the python distribution eofs (Dawson, 2016) for computation of the EOF. We used the Hector software (Bos et al., 2013) to estimate the trend uncertainties from MLE and study the impact of different noise models, as described in section 2.3.1. We thank Samantha Royston for the useful discussion concerning noise models in time series analysis.

References

Ablain, M., Legeais, J., Prandi, P., Marcos, M., Fenoglio-Marc, L., Dieng, H., et al. (2016). Satellite altimetry-based sea level at global and regional scales. Surveys Geophys. 38, 7–31. doi: 10.1007/s10712-016-9389-8

CrossRef Full Text | Google Scholar

Akaike, H. (1998). Information Theory and an Extension of the Maximum Likelihood Principle. New York, NY: Springer New York, 199–213. doi: 10.1007/978-1-4612-1694-0_15

CrossRef Full Text | Google Scholar

Andersson, H. C. (2002). Influence of long-term regional and large-scale atmospheric circulation on the Baltic sea level. Tellus A 54, 76–88. doi: 10.3402/tellusa.v54i1.12125

CrossRef Full Text | Google Scholar

Beckley, B. D., Callahan, P. S., Hancock, D. III, Mitchum, G., and Ray, R. (2017). On the Cal-Mode correction to TOPEX satellite altimetry and its effect on the global mean sea level time series. J. Geophys. Res. 122, 8371–8384. doi: 10.1002/2017JC013090

CrossRef Full Text | Google Scholar

Benveniste, J., Birol, F., Calafat, F., Cazenave, A., Dieng, H., Gouzenes, Y., et al. (2020). (The Climate Change Initiative Coastal Sea Level Team) Coastal sea level anomalies and associated trends from Jason satellite altimetry over 2002-2018. Sci. Data 7:357. doi: 10.1038/s41597-020-00694-w

PubMed Abstract | CrossRef Full Text | Google Scholar

Benveniste, J., Cazenave, A., Vignudelli, S., Fenoglio-Marc, L., Shah, R., Almar, R., et al. (2019). Requirements for a coastal hazards observing system. Front. Mar. Sci. 6:348. doi: 10.3389/fmars.2019.00348

CrossRef Full Text | Google Scholar

Bos, M. S., Fernandes, R. M. S., Williams, S. D. P., and Bastos, L. (2013). Fast error analysis of continuous GNSS observations with missing data. J. Geod. 87, 351–360. doi: 10.1007/s00190-012-0605-0

CrossRef Full Text | Google Scholar

Bosch, W., Dettmering, D., and Schwatke, C. (2014). Multi-mission cross-calibration of satellite altimeters: constructing a long-term data record for global and regional sea level change studies. Remote Sens. 6, 2255–2281. doi: 10.3390/rs6032255

CrossRef Full Text | Google Scholar

Boutin, J., Martin, N., Reverdin, G., Yin, X., and Gaillard, F. (2013). Sea surface freshening inferred from SMOS and ARGO salinity: impact of rain. Ocean Sci. 9, 183–192. doi: 10.5194/os-9-183-2013

CrossRef Full Text | Google Scholar

Brown, G. (1977). The average impulse response of a rough surface and its applications. IEEE Trans. Anten. Propag. 25, 67–74. doi: 10.1109/TAP.1977.1141536

CrossRef Full Text | Google Scholar

Carrere, L., Lyard, F., Cancet, M., and Guillot, A. (2015). “FES 2014, a new tidal model on the global ocean with enhanced accuracy in shallow seas and in the arctic region,” in EGU General Assembly Conference Abstracts, Vol. 17 (Wien), 5481.

Google Scholar

Celebi, M. (2014). Partitional Clustering Algorithms. EBL-Schweitzer. Springer International Publishing. doi: 10.1007/978-3-319-09259-1

CrossRef Full Text | Google Scholar

Chafik, L., Nilsen, J. E. Ø., and Dangendorf, S. (2017). Impact of north atlantic teleconnection patterns on northern European sea level. J. Mar. Sci. Eng. 5:43. doi: 10.3390/jmse5030043

CrossRef Full Text | Google Scholar

Collecte Localisation Satellites (CLS). Dynamic Atmospheric Corrections Are Produced by CLS Space Oceanography Division Using the MOG2D Model From Legos and Distributed by AVISO+ With Support from CNES. AVISO+.

Dangendorf, S., Mudersbach, C., Wahl, T., and Jensen, J. (2013). Characteristics of intra-, inter-annual and decadal sea-level variability and the role of meteorological forcing: the long record of Cuxhaven. Ocean Dyn. 63, 209–224. doi: 10.1007/s10236-013-0598-0

CrossRef Full Text | Google Scholar

Dawson, A. (2016). EOFS: a library for EOF analysis of meteorological, oceanographic, and climate data. J. Open Res. Softw. 4:e14. doi: 10.5334/jors.122

CrossRef Full Text | Google Scholar

Dettmering, D., Wynne, A., Müller, F. L., Passaro, M., and Seitz, F. (2018). Lead detection in polar oceansÑA comparison of different classification methods for Cryosat-2 SAR data. Remote Sens. 10:1190. doi: 10.3390/rs10081190

CrossRef Full Text | Google Scholar

Fernandes, M. J., and Lazaro, C. (2016). GPD+ wet tropospheric corrections for cryosat-2 and GFO altimetry missions. Remote Sens. 8:851. doi: 10.3390/rs8100851

CrossRef Full Text | Google Scholar

Fernandes, M. J., Lázaro, C., Ablain, M., and Pires, N. (2015). Improved wet path delays for all esa and reference altimetric missions. Remote Sens. Environ. 169, 50–74. doi: 10.1016/j.rse.2015.07.023

CrossRef Full Text | Google Scholar

Gräwe, U., Klingbeil, K., Kelln, J., and Dangendorf, S. (2019). Decomposing mean sea level rise in a semi-enclosed basin, the Baltic sea. J. Clim. 32, 3089–3108. doi: 10.1175/JCLI-D-18-0174.1

CrossRef Full Text | Google Scholar

Guinehut, S., Le Traon, P., Larnicol, G., and Philipps, S. (2004). Combining Argo and remote-sensing data to estimate the ocean three-dimensional temperature fieldsÑa first approach based on simulated observations. J. Mar. Syst. 46, 85–98. doi: 10.1016/j.jmarsys.2003.11.022

CrossRef Full Text | Google Scholar

Hastie, T, Tibshirani, R, Friedman, J. (2009). The Elements of Statistical Learning, Vol. 2. New York, NY: Springer-Verlag. doi: 10.1007/978-0-387-84858-7

CrossRef Full Text | Google Scholar

Hayne, G. S. (1980). Radar altimeter mean return waveforms from near-normal-incidence ocean surface scattering. IEEE Trans. Antenn. Propag. 28, 687–692. doi: 10.1109/TAP.1980.1142398

CrossRef Full Text | Google Scholar

Hersbach, H., Bell, B., Berrisford, P., Hirahara, S., Horányi, A., Mu noz-Sabater, J., et al. (2020). The ERA5 global reanalysis. Q. J. R. Meteorol. Soc. 146, 1999–2049. doi: 10.1002/qj.3803

CrossRef Full Text | Google Scholar

Jevrejeva, S., Moore, J., Woodworth, P., and Grinsted, A. (2005). Influence of large-scale atmospheric circulation on European sea level: results based on the wavelet transform method. Tellus A 57, 183–193. doi: 10.3402/tellusa.v57i2.14609

CrossRef Full Text | Google Scholar

Johansson, M. M., Pellikka, H., Kahma, K. K., and Ruosteenoja, K. (2014). Global sea level rise scenarios adapted to the finnish coast. J. Mar. Syst. 129, 35–46. doi: 10.1016/j.jmarsys.2012.08.007

CrossRef Full Text | Google Scholar

Karabil, S., Zorita, E., and Hünicke, B. (2018). Contribution of atmospheric circulation to recent off-shore sea-level variations in the Baltic sea and the north sea. Earth Syst. Dyn. 9, 69–90. doi: 10.5194/esd-9-69-2018

CrossRef Full Text | Google Scholar

Kusche, J., Uebbing, B., Rietbroek, R., Shum, C., and Khan, Z. (2016). Sea level budget in the Bay of Bengal (2002-2014) from GRACE and altimetry. J. Geophys. Res. 121, 1194–1217. doi: 10.1002/2015JC011471

CrossRef Full Text | Google Scholar

Landskron, D., and Böhm, J. (2018). VMF3/GPT3: refined discrete and empirical troposphere mapping functions. J. Geodesy 92, 349–360. doi: 10.1007/s00190-017-1066-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Legeais, J.-F., Ablain, M., Zawadzki, L., Zuo, H., Johannessen, J. A., Scharffenberg, M. G., et al. (2018). An improved and homogeneous altimeter sea level record from the ESA climate change initiative. Earth Syst. Sci. Data 10:281. doi: 10.5194/essd-10-281-2018

CrossRef Full Text | Google Scholar

Leppäranta, M., and Myrberg, K. (2009). Physical Oceanography of the Baltic Sea. Springer Science & Business Media. doi: 10.1007/978-3-540-79703-6

CrossRef Full Text | Google Scholar

Ludwigsen, C. A., Khan, S. A., Andersen, O. B., and Marzeion, B. (2020). Vertical Land Motion from present-day deglaciation in the wider Arctic. Geophys. Res. Lett. 47:e2020GL088144. doi: 10.1029/2020GL088144

CrossRef Full Text | Google Scholar

Madsen, K. S., Höyer, J. L., Suursaar, Ü., She, J., and Knudsen, P. (2019). Sea level trends and variability of the Baltic sea from 2d statistical reconstruction and altimetry. Front. Earth Sci. 7:243. doi: 10.3389/feart.2019.00243

CrossRef Full Text | Google Scholar

Müller, F. L., Dettmering, D., Bosch, W., and Seitz, F. (2017). Monitoring the arctic seas: How satellite altimetry can be used to detect open water in sea-ice regions. Remote Sens. 9:551. doi: 10.3390/rs9060551

CrossRef Full Text | Google Scholar

Olivieri, M., and Spada, G. (2016). Spatial sea-level reconstruction in the Baltic Sea and in the Pacific Ocean from tide gauges observations. Ann. Geophys. 59:0323. doi: 10.4401/ag-6966

CrossRef Full Text | Google Scholar

Passaro, M., Dinardo, S., Quartly, G. D., Snaith, H. M., Benveniste, J., Cipollini, P., et al. (2016). Cross-calibrating ALES Envisat and Cryosat-2 Delay-Doppler: a coastal altimetry study in the Indonesian Seas. Adv. Space Res. 58, 289–303. doi: 10.1016/j.asr.2016.04.011

CrossRef Full Text

Passaro, M., Mueller, F., and Dettmering, D. (2020a). Baltic+ SEAL: Algorithm Theoretical Baseline Document (ATBD). Technical report delivered under the Baltic+ SEAL project, European Space Agency. doi: 10.5270/esa.BalticSEAL.ATBDV1.1

CrossRef Full Text

Passaro, M., Nadzir, Z., and Quartly, G. D. (2018b). Improving the precision of sea level data from satellite altimetry with high-frequency and regional sea state bias corrections. Remote Sens. Environ. 218, 245–254. doi: 10.1016/j.rse.2018.09.007

CrossRef Full Text | Google Scholar

Passaro, M., Restano, M., Sabatino, G., and C., O., and J, B. (2020b). “The ALES+ SAR service for Cryosat-2 and Sentinel-3 at ESA GPOD,” in Presented at the Ocean Surface Topography Science Team (OSTST) Meeting. Available online at: https://meetings.aviso.altimetry.fr/index.html

Google Scholar

Passaro, M., Rose, S., Andersen, O., Boergens, E., Calafat, F., and J, B. (2018a). ALES+: Adapting a homogenous ocean retracker for satellite altimetry to sea ice leads, coastal and inland waters. Remote Sens. Environ. 211, 456–471. doi: 10.1016/j.rse.2018.02.074

CrossRef Full Text | Google Scholar

Petit, G., and Luzum, B. (2010). Iers Conventions 2010, IERS Technical Note; 36. Frankfurt Am Main: Verlag des bundesamts für kartographie und geodäsie. Available online at: http://www.iers.org/TN36

Google Scholar

Preisendorfer, R. W. (1988). Principal Component Analysis in Meteorology and Oceanography. Amsterdam: Elsevier.

Google Scholar

Ringgaard, I., Madsen, K., Müller, F., Tuomi, L., Rautiainen, L., and Passaro, M. (2020). Dataset Description, Version 1.1. Technical report delivered under the Baltic+ SEAL project, ESA Contract 4000126590/19/I/BG - BALTIC+ SEAL (Sea Level), European Space Agency. doi: 10.5270/esa.BalticSEAL.DDV1.1

CrossRef Full Text

Rio, M., Piollé, J.-F., J.A., J., and Donlon, C. (2014). GlobCurrent Product Format and Content. European Space Agency. Available online at: http://globcurrent.ifremer.fr/products-data/products-overview (accessed November 10, 2020).

Rose, S. K., Andersen, O. B., Passaro, M., Ludwigsen, C. A., and Schwatke, C. (2019). Arctic ocean sea level record from the complete radar altimetry era: 1991-2018. Remote Sens. 11:1672. doi: 10.3390/rs11141672

CrossRef Full Text | Google Scholar

Royston, S., Watson, C. S., Legresy, B., King, M. A., Church, J. A., and Bos, M. S. (2018). Sea-level trend uncertainty with Pacific climatic variability and temporally-correlated noise. J. Geophys. Res. 123, 1978–1993. doi: 10.1002/2017JC013655

CrossRef Full Text | Google Scholar

Scharroo, R., and Smith, W. H. F. (2010). A global positioning system based climatology for the total electron content in the ionosphere. J. Geophys. Res. 115:A10318. doi: 10.1029/2009JA014719

CrossRef Full Text | Google Scholar

Schwarz, G. (1978). Estimating the dimension of a model. Ann. Stat. 6, 461–464. doi: 10.1214/aos/1176344136

CrossRef Full Text | Google Scholar

Stephenson, D., Pavan, V., Collins, M., Junge, M., Quadrelli, R., and Participating CMIP2 Modelling Groups (2006). North Atlantic Oscillation response to transient greenhouse gas forcing and the impact on European winter climate: a CMIP2 multi-model assessment. Clim. Dyn. 27, 401–420. doi: 10.1007/s00382-006-0140-x

CrossRef Full Text | Google Scholar

Storch, H. V., and Zwiers, F. W. (1999). Statistical Analysis in Climate Research. Cambridge: Cambridge University Press.

Google Scholar

Taburet, G., Sanchez-Roman, A., Ballarotta, M., Pujol, M.-I., Legeais, J.-F., Fournier, F., et al. (2019). DUACS DT2018: 25 years of reprocessed sea level altimetry products. Ocean Sci. 15, 1207–1224. doi: 10.5194/os-15-1207-2019

CrossRef Full Text | Google Scholar

Tran, N., Labroue, S., Philipps, S., Bronner, E., and Picot, N. (2010). Overview and update of the sea state bias corrections for the Jason-2, Jason-1 and TOPEX missions. Mar. Geodesy 33, 348–362. doi: 10.1080/01490419.2010.487788

CrossRef Full Text | Google Scholar

Vestøl, O., Ågren, J., Steffen, H., Kierulf, H., and Tarasov, L. (2019). NKG2016LU: a new land uplift model for Fennoscandia and the Baltic Region. J. Geodesy 93, 1759–1779. doi: 10.1007/s00190-019-01280-8

CrossRef Full Text | Google Scholar

Von Schuckmann, K., Le Traon, P.-Y., Alvarez-Fanjul, E., Axell, L., Balmaseda, M., Breivik, L.-A., et al. (2016). The Copernicus marine environment monitoring service ocean state report. J. Oper. Oceanogr. 9(Suppl. 2), s235–s320. doi: 10.1080/1755876X.2016.1273446

CrossRef Full Text | Google Scholar

Wilks, D. S. (2006). Statistical Methods in the Atmospheric Sciences, 2nd Edn. International Geophysics Series, Vol. 91. Burlington, MA: Academic Press.

Xu, R., and Wunsch, D. C. (2008). Clustering. Hoboken, NJ: John Wiley & Sons, Inc. doi: 10.1002/9780470382776

CrossRef Full Text | Google Scholar