Probabilistic projections and past trends of sea level rise in Finland

We explore past trends and future projections of mean sea level (MSL) at the Finnish coast, in the northeastern Baltic Sea, during the period 1901–2100. We decompose the relative MSL change into three components: regional sea level rise (SLR), postglacial land uplift, and the effect of changes in wind climate. Past trends of regional SLR can be calculated after subtracting the other two components from the MSL trends observed by tide gauges, as the land uplift rates obtained from the semi-empirical model NKG2016LU are independent of tide gauge observations. According to the results, local absolute SLR trends are close to global mean rates. To construct future projections, we combine an ensemble of global SLR projections in a probabilistic framework. In addition, we use climate model results to estimate future changes in wind climate and their effect on MSL in the semi-enclosed Baltic Sea. This yields probability distributions of MSL change for three scenarios representing different future emission pathways. Spatial variations in the MSL projections result primarily from different local land uplift rates: under the medium-emission scenario RCP4.5/SSP2-4.5, for example, the projected MSL change (5 % to 95 % range) over the 21st century varies from −28 (−54 to 24) cm in the Bothnian Bay to 31 (5 to 83) cm in the eastern Gulf of Finland.

1 Introduction

Sea level rise (SLR) is an existential threat to many coastal communities worldwide, but the magnitude and rate of future SLR are still shrouded in uncertainty. Over the last decade, advances in sea level projections have shifted the focus from global to regional assessments (e.g. Slangen et al., 2014; Kopp et al., 2014; Grinsted et al., 2015; Fox-Kemper et al., 2021) and from best estimates with uncertainty ranges to full probability distributions (Jevrejeva et al., 2014; Grinsted et al., 2015; Goodwin et al., 2017; Le Bars et al., 2017). A probability distribution of expected local mean sea level change allows different risk levels to be determined for different operations and types of infrastructure. Nuclear power plants, for example, need to be protected against much more unlikely risks than conventional buildings.

This paper examines past and future mean sea level (MSL) at the Finnish coast, located in the northeastern Baltic Sea (Fig. 1). The Baltic Sea is a semi-enclosed basin with a narrow ocean connection and a characteristic sea level behaviour that differs from open oceans. Wind conditions create decadal variations in sea level at the northeastern Baltic coasts, in Finland and Estonia (Johansson et al., 2014; Suursaar et al., 2006). Moreover, Fennoscandia is an area of substantial postglacial land uplift, as the glacial isostatic adjustment (GIA) following the last ice age is still an ongoing process (Poutanen and Steffen, 2014). Land uplift, in particular, causes the MSL change in the study area to deviate markedly from the global mean.

Figure 1 Map of the Baltic Sea region, with Finnish tide gauges marked with triangles. The star denotes the point (55^∘ N, 15^∘ E) with the best correlation between the zonal geostrophic wind speed component and mean sea level at the Finnish coast.

Despite major advances in the understanding of global mean sea level rise (GMSLR) and its causes, it has proved very difficult to nail down the upper end of the projection range. Upper limits of published GMSLR estimates have varied widely over the past 40 years when modern sea level projections have been available, and the total range of estimates in the literature has expanded rather than contracted (Garner et al., 2018). There is no consensus about the upper extreme tail of the probability distribution. This is mainly because of poorly understood instability mechanisms of marine ice sheets – continental ice resting on ground below sea level and in direct contact with the ocean, mainly found in West Antarctica (Fox-Kemper et al., 2021). Currently, one of the main disputes concerns the so-called marine ice cliff instability (MICI) hypothesis, the structural failure of ice cliffs after the loss of buttressing by floating ice shelves (DeConto and Pollard, 2016; Edwards et al., 2019; DeConto et al., 2021). Incorporating ice cliff failure in ice sheet models increases the Antarctic SLR contribution considerably, especially after 2100 (Fox-Kemper et al., 2021).

The problem for coastal planners and other end users of SLR projections is that there are numerous projections available, which differ considerably from each other. The Intergovernmental Panel on Climate Change (IPCC) Assessment Reports provide an authoritative overview of the state of climate science, and this information is widely used as a basis for national SLR assessments, e.g. in Norway (Simpson et al., 2015) and Sweden (Hieronymus and Kalén, 2020). A survey covering 32 European countries found IPCC to be the primary source of information for SLR planning (McEvoy et al., 2021). There are possible pitfalls in relying solely on IPCC projections, especially in sectors where risk aversion is critical. The consensus-based approach of IPCC tends to produce rather conservative projections, and the upper limits of the IPCC projections have been consistently lower than the upper limits reported in individual studies (Garner et al., 2018). The Fifth Assessment Report (AR5; IPCC, 2013) was criticized for disregarding low-confidence information and omitting the upper tail of the probability distribution from the numerical SLR projections (e.g. Hinkel et al., 2015; Siegert et al., 2020). In the recently published Sixth Assessment Report (AR6, 2021), IPCC for the first time provides local projections spanning the whole probability range. In addition, the marine ice sheet uncertainty has been incorporated as a separate “low-likelihood, high-impact” storyline.

As Behar et al. (2017) point out, there has been a need to use multiple analyses and probability distributions in adaptation planning. In Finland, the approach chosen in coastal management and high-risk applications, such as the safety analyses of coastal nuclear power plants (Jylhä et al., 2018), has been to draw from a wide body of research in an attempt to create local SLR projections covering the full probability range (Johansson et al., 2014; Pellikka et al., 2018). In this study, we build on previous work to calculate projections of mean sea level change in Finland. It has become timely to update earlier national projections after the publication of AR6 and other recent studies. Even though IPCC now provides localized projections, it is still worthwhile to use regional models and perform local analyses to improve local projections. As a small inland sea with a narrow ocean connection, Baltic Sea and its water exchange with the Atlantic Ocean is not well represented in global climate models. Moreover, the probability distributions presented in Johansson et al. (2014) and Pellikka et al. (2018) do not differentiate between future greenhouse gas emission pathways. Merging all uncertainties into a single probability distribution was a choice based on practical needs, but this approach leads to a loss of information on how climate policies and other societal changes affect future SLR outcomes.

The Finnish tide gauge network has been measuring sea levels at the Finnish coast for over a century. Historically, land uplift has overruled sea level rise, causing declining trends in relative sea level and expansion of land area. It has, however, been difficult to separate the different drivers of MSL change, as estimates of land uplift have partially relied on tide gauge observations. In recent years, the development of geodetic measurement techniques, such as Global Navigation Satellite System (GNSS) observations, has enabled independent estimates of land uplift rates. The current land uplift model for Fennoscandia, NKG2016LU (Vestøl et al., 2019), has already been used to estimate the rates of absolute SLR in the Baltic Sea by subtracting land uplift from observed relative sea level trends (Suursaar and Kall, 2018; Madsen et al., 2019; Passaro et al., 2021).

We use the long time series of Finnish tide gauge data to analyse past trends in MSL and its components. This analysis contributes to the understanding of regional SLR trends in the Baltic Sea, which has already been addressed in a number of studies. Using satellite altimetry, Stramska and Chudziak (2013) reported a trend of 3.3 mm a⁻¹ over 1992–2012, while Madsen et al. (2019) arrived at 3.4 mm a⁻¹ over 1993–2014, in accordance with the GMSLR trend. Over the 20th century, the Baltic Sea average was 1.3 mm a⁻¹ according to Madsen et al. (2019), but Suursaar and Kall (2018) reported 2.5 mm a⁻¹ in 1901–2010 based on Estonian tide gauge data. The difference may be explained by significant spatial variability within the Baltic Sea, with larger rates in northern and eastern parts of the basin and lower rates in the southwestern areas (Gräwe et al., 2019; Madsen et al., 2019; Passaro et al., 2021). This pattern is attributed to wind and air pressure changes (Gräwe et al., 2019).

In this study, we utilize tide gauge data, modelling results, and SLR projections from scientific literature to analyse past and future MSL in Finland. The following research questions summarize the aims of this paper:

1. What was the rate of absolute sea level rise in the northeastern Baltic Sea over the past century and the satellite altimetry era?

2. What are the expected changes in mean sea level at the coast of Finland by 2100?

3. How do the probability distributions of SLR provided by IPCC AR6 compare with other recently published projections, and can they satisfy the needs of high-risk coastal management?

While the results of this study are mostly applicable for national planning purposes in Finland, the third question is not confined to the study area, and our analysis of the probability distribution of GMSLR may be of interest to coastal planners in other countries, as well.

2 Data and methods

2.1 Mean sea level components at the Finnish coast

Throughout this paper, mean sea level (MSL) refers to long-term sea level relative to land, while sea level rise (SLR) refers to the absolute sea level change relative to the centre of the Earth. SLR is a global process driven by thermal expansion of seawater and melting of glaciers and ice sheets, but local and regional SLR rates may deviate from the global mean. In the terminology of this paper, SLR is one component of the local MSL change. The MSL change in the study area may be positive or negative, depending on whether the land uplift is strong enough to overrule SLR.

We construct probability distributions of future MSL in Finland by combining three components: regional SLR, post-glacial land uplift, and the effect of changes in wind climate. Johansson et al. (2014) have shown that observed annual mean sea levels at the Finnish coast can be reproduced with good accuracy by computing the sum of these three components. More formally,

$$\begin{array}{}\text{(1)}& h_{\mathrm{m}}(i,t)=r(i,t,t_{\mathrm{0}})-d\left(i\right)(t-t_{\mathrm{0}})+w(i,t)+R\left(i\right),\end{array}$$

where h_m is the estimate of local MSL at the tide gauge i and time t, $r(i,t,t_{\mathrm{0}})$ is the regional SLR relative to some reference year t₀, d is the rate of land uplift, w is the wind-induced MSL component, and R is a constant used to level h_m in relation to the Finnish N2000 height reference system (Saaranen et al., 2009) and the MSL of the reference year t₀. All sea level data used in this analysis are taken from the database of the Finnish Meteorological Institute (FMI).

Figure 2 Schematic illustration of the method. A set of global mean sea level rise (GMSLR) projections is normalized to a common time span. Probability distributions are fitted to each projection and combined to produce a probability distribution of sea level rise (SLR), then combined with land uplift and wind-induced changes to produce local relative mean sea level (MSL) projections. This calculation is performed separately for different emission scenarios.

Our method is presented schematically in Fig. 2. To obtain a probability distribution of r, we use a set of GMSLR projections published over the last decade (Sect. 2.2), also considering regional anomalies (Sect. 2.3). Land uplift rates d are based on geodetic observations and modelling (Sect. 2.4). The wind-induced component w is estimated using climate model projections of geostrophic wind (Sect. 2.5). Finally, we calculate the sum of the three MSL components to obtain a probability distribution of future MSL at a selected location at the Finnish coast. We treat the components as independent random variables, so the MSL distribution is the convolution of the distributions of the three components.

2.2 Probabilistic projections for sea level rise

To take into account the full range of potential future SLR, including low-probability, high-impact events such as marine ice sheet disintegration and controversial mechanisms such as MICI, we use a wide range of GMSLR projections published over the last decade. The projections are listed in Table 1. They are based on different approaches: (i) process-based models, i.e. physical models that simulate the individual components contributing to SLR, including thermal expansion, ocean dynamics, and the melting and disintegration of land-based ice; (ii) semi-empirical models, which project GMSLR indirectly by constructing a statistical relationship between global mean sea level and some other climate variable, e.g. the global mean temperature; and (iii) expert surveys that assess the uncertainty in SLR projections by polling experts in the field and reviewing their understanding.

Jevrejeva et al. (2012)Jevrejeva et al. (2014)Grinsted et al. (2015)Kopp et al. (2014)Kopp et al. (2016)Mengel et al. (2016)Goodwin et al. (2017)Kopp et al. (2017)Kopp et al., 2014DeConto and Pollard, 2016Le Bars et al. (2017)DeConto and Pollard, 2016

Table 1 Sea level rise projections used in this study: median values and 5 %–95 % ranges (in parentheses). All values are in centimetres in 2100 relative to 1995–2014. The low-emission scenario refers to RCP2.6/SSP1-2.6, medium to RCP4.5/SSP2-4.5, and high to RCP8.5/SSP5-8.5. AR6MC and AR6LC are the medium- and low-confidence projections of IPCC AR6 (Fox-Kemper et al., 2021), respectively.

Scenarios of future greenhouse gas emissions play a crucial role in SLR projections. In IPCC AR6, the Representative Concentration Pathway (RCP) scenarios used in AR5 were replaced by a new set of emission scenarios called Shared Socioeconomic Pathways (SSPs). The SSP and RCP scenarios are not directly comparable, but both are labelled by the level of radiative forcing reached in 2100. In the pre-AR6 literature, the most widely used RCP scenarios were RCP2.6, RCP4.5, and RCP8.5, which we also use in this paper since there is a sufficient number of GMSLR projections available for these scenarios. The corresponding SSP scenarios are SSP1-2.6, SSP2-4.5, and SSP5-8.5. The SSP scenarios generally result in slightly stronger global warming than the corresponding RCP scenarios (Lee et al., 2021), and there are other differences as well, but for the purposes of this study, it is reasonable to group together the SSP and RCP scenarios with the same nominal radiative forcing. Hence, we construct probability distributions for three emission scenarios, which we call low (2.6), medium (4.5), and high (8.5).

Our ensemble of GMSLR projections includes 8 members for the low-emission scenario, 9 for the medium scenario, and 10 for the high scenario (Table 1). There are two sets of IPCC AR6 projections, medium confidence and low confidence (hereafter, AR6MC and AR6LC), which are included in our ensemble separately. In contrast to AR6MC, AR6LC incorporates information from structured expert judgement (Bamber et al., 2019) and ice sheet simulations (DeConto et al., 2021) that explore the potential for rapid ice discharge from Antarctica through MICI. The probability of these high-end projections cannot be robustly quantified; therefore, IPCC has decided to include them in a separate low-confidence projection to serve the needs of stakeholders with a low risk tolerance. AR6MC and AR6LC differ from each other particularly under the high-emission scenario.

IPCC AR6 sea level projections (Fox-Kemper et al., 2021) were downloaded from Physical Oceanography Distributed Active Archive Center (PO.DAAC; https://podaac.jpl.nasa.gov/announcements/2021-08-09-Sea-level-projections-from-the-IPCC-6th-Assessment-Report, last access: 20 April 2023) on 29 September 2021. Other data related to SLR projections are extracted from the tables or supplementary material of the referenced publications.

For each emission scenario, the ensemble of projections is transformed into a probability distribution of GMSLR in a three-step process (Fig. 2). First, the projections are normalized to a common baseline, 1995–2014. Second, a probability distribution is fitted to each projection. Third, the ensemble of distributions is combined to form a single probability distribution.

1. Normalization. All projections are normalized to a common baseline, 1995–2014, by assuming an initial GMSLR rate of 3.25 ± 0.37 mm a⁻¹ (Fox-Kemper et al., 2021 for 1995–2018) and constant acceleration from the initial rate. In practice, to normalize a projection having a baseline year y₀, we remove the rise realized between y₀ and 2005 after making a second-order fit to the projected SLR from y₀ to 2100 and the initial SLR rate.

2. Fitting distributions. To create an ensemble of continuous probability distributions, we fit distributions to the percentiles of the published GMSLR projections. The number of percentiles available for the fit vary from 3 (median and the 5 %–95 % range) to 31, depending on projection. There is no theoretical basis for choosing what type of distribution might best represent future SLR, but it is well known that the probability distributions of GMSLR are positively skewed with a fat upper tail. The tail represents the low-probability, high-impact scenarios that are hard to quantify. We experiment with three types of distributions that allow a positive skew: Weibull, Fréchet, and skew normal. It turns out that no single distribution gives a good fit to all projections, as the shape of the distribution varies between projections. We choose the best fit for each projection, defined as the smallest residual sum of squares (RSS) between the original data points and the fitted distribution in semi-log space (x, log (1−F), where x is sea level rise and F is the corresponding cumulative probability).
   The fits to individual projections are shown in Fig. A1 for the high-emission scenario. For IPCC AR6, which provides a large number of percentiles, we do not make a fit but instead interpolate and extrapolate the distribution linearly in semi-log space to make it continuous. Also, Jevrejeva et al. (2014) provide a continuous probability distribution, which we use as such (digitized from their Fig. 3).

3. Combining distributions. Finally, we calculate a mixture distribution to obtain the probability distribution of GMSLR for each emission scenario:

   $$\begin{array}{}\text{(2)}& F\left(x\right)=\sum _{i=\mathrm{1}}^n{\mathit{\lambda }}_i{P}_i\left(x\right),\end{array}$$

   where F is the cumulative distribution function (CDF) of GMSLR by 2100, P_i are the CDFs of individual projections (the n ensemble members), and λ_i are weights attributed to each projection. Setting the weights is a matter of expert judgement. It is reasonable to assign different weights to the ensemble members, if there are plausible arguments for some projections and projection methods to be more credible than others. We give a lower weight (0.5) to projections that incorporate MICI (Goodwin et al., 2017; Le Bars et al., 2017; and AR6LC) to factor in the low confidence in these projections. AR6MC is given a large weight (4) as it can be thought to represent the latest scientific consensus on SLR projections, and other projections are given equal weight (1). We explore different choices of weights, and the effects of the choices on the resulting projections, in Sect. 3.2.

2.3 Considerations of regional SLR anomaly

Various processes create spatial inhomogeneity in SLR, including variations in ocean temperature, salinity, and circulation as well as gravitational, rotational, and solid Earth responses to ice melt. Each melting glacier produces a characteristic geographical pattern (fingerprint) of SLR, resulting from changes in Earth's gravity field and rotation as well as crustal deformation that follow the reduced weight of the ice mass. Sea level rise is smaller than average in the proximity of the melting glacier or ice sheet, even negative up to the distance of ca. 2000 km, because of elastic uplift of the crust and the reduced gravitational pull of the ice mass (Mitrovica et al., 2011).

At the Finnish coast, the fingerprint effect makes the SLR contributions of the Greenland Ice Sheet and other northern glaciers notably smaller compared to the global mean. The contribution of Greenland is essentially zero in the Bothnian Bay and ca. 15 % of the global mean at the southern Finnish coast (Mitrovica et al., 2001; Kopp et al., 2014). On the other hand, according to global climate models, thermal expansion and ocean dynamical effects (hereafter referred to as the ocean component) produce a larger-than-average SLR contribution in the Baltic Sea. Pellikka et al. (2018) analysed the ocean component from Coupled Model Intercomparison Project Phase 5 (CMIP5) models at the Finnish coast: median estimates for different emission scenarios were 20 % to 60 % higher compared to the global mean. This additional SLR was, however, overruled by the reduction caused by glacier and ice sheet fingerprints, and the total SLR (excluding land uplift) was about 80 % of the global mean in end-of-century projections.

Figure 3 Median sea level rise projections of IPCC AR6 (medium confidence, 2100 relative to 1995–2014) divided into different components, for three emission scenarios. The ocean component includes thermal expansion and ocean dynamical effects. The global projection is shown together with three local projections: Vaasa and Helsinki at the Finnish coast, Hanstholm at the North Sea coast of Denmark.

In IPCC AR6, which uses the current CMIP6 climate models, the local ocean component at the Finnish coast is larger than previously estimated: almost twice the global mean in 2100. Figure 3 shows the proportions of different SLR components – the ocean component, different glaciers and ice sheets, land water storage, and vertical land motion – in the AR6MC median projections locally and globally. At the Finnish coast, the large ocean component compensates the reduced glacier components so that the absolute SLR is approximately the same as the global mean. Projections for Hanstholm, at the North Sea coast of Denmark, are plotted for comparison and show no notable difference. Therefore, the large ocean component is not a feature confined to the Baltic Sea basin, whose sea level dynamics the global models are not able to simulate well. The reason for the notably larger ocean component in CMIP6 compared to CMIP5 seems to be the increased dynamic sea level rise in the North Atlantic, associated with a weakening of the Atlantic meridional overturning circulation (Lyu et al., 2020).

It is important to note that the conformity of global and local absolute SLR takes place only because two phenomena happen to cancel each other: the smaller-than-average glacial contribution and the larger-than-average ocean contribution. This may not be true for other, methodologically different projections or projections with different time spans, as the relative proportions of the SLR contributors may be different. Localized projections based on essentially different methods (such as semi-empirical modelling) are not available, however. To study projections with different time spans, we calculated the ratio of local to global SLR from AR6 median projections every 10 years from 2030 to 2150. In Helsinki, for example, the ratio declines from 1.2–1.4 in 2030 to 1–1.1 in 2100 and 0.9–1 in 2150 (Fig. 4), as the contribution of land ice melt increases and gradually tends to dominate the projected SLR.

Figure 4 Ratio of local absolute sea level rise in Helsinki (without the effect of land uplift) to the global mean in IPCC AR6 median projections. Medium-confidence projections (AR6MC) are plotted with dark lines, low-confidence projections (AR6LC) with pale dashed lines, while colours indicate different emission scenarios: low (SSP1-2.6), medium (SSP2-4.5), and high (SSP5-8.5).

Based on this analysis, we apply the GMSLR projections for 2100 to the Finnish coast as such, without regional adjustments. In other words, we approximate the regional SLR $r(i,t,t_{\mathrm{0}})$ in Eq. (1) with the global mean SLR, g(t,t₀). Using a scaling approach similar to Pellikka et al. (2018) to convert GMSLR projections into regional ones would lead to small adjustments in the final projections, which are minor compared to uncertainties involved in such analysis. Therefore, we omit these adjustments for simplicity.

2.4 Land uplift

Vertical land motions (VLMs) in the northern Baltic Sea area are dominated by glacial isostatic adjustment (GIA). Slow crustal recovery from the pressure of the ice sheet that covered the area during the last glacial period causes post-glacial land uplift. At the Finnish coast, the rate of land uplift is of the same order of magnitude (3 to 9 mm a⁻¹) as the current rate of GMSLR (3.7 mm a⁻¹; Fox-Kemper et al., 2021). The uplift rate can be considered constant in time over the next few centuries (Poutanen and Steffen, 2014).

We use land uplift rates calculated by NKG2016LU, a semi-empirical model of land uplift in the Fennoscandian area (Vestøl et al., 2019). The model has been computed in the Working Group of Geoid and Height Systems of the Nordic Geodetic Commission (NKG) and combines two parts: (i) the observation-based part that uses geodetic observations of land uplift, namely levelling and GNSS time series in the area, and (ii) the geophysical GIA model that supplements data to areas where observations are sparse. The GIA model is fitted to the geodetic observations using least squares collocation (LSC). The uncertainty estimate of the model is a combination of the observations' and the GIA model's uncertainty. The levelled version used in the current study gives the uplift relative to the geoid (Fig. 5). Land uplift rates, d, at the Finnish tide gauges are listed in Table 3 in “Results” (Sect. 3.1). We assume a normal (Gaussian) uncertainty distribution for the land uplift, using the 1 SD uncertainty ranges of the model to fit the distributions.

Figure 5 The present-day land uplift (millimetres per year) over Fennoscandia relative to geoid given by the semi-empirical land uplift model NKG2016LU_lev.

The local projections of IPCC AR6 account for VLM by calculating background rates of mean sea level change from tide gauge data, using the statistical model of Kopp et al. (2014). The long-term background trend is scenario-independent and assumed constant over the projection period. This way, other drivers of long-term VLM such as tectonics, volcanism, and anthropogenic subsidence can be considered in addition to GIA. IPCC assigns low to medium confidence in VLM projections and points out that in many regions, more detailed regional analyses would be needed to produce higher-fidelity projections. In Fennoscandia, VLM is dominated by GIA, which is a relatively well known and predictable process. Comparing the VLM rates of AR6 and the NKG2016LU model, the differences are small and result in a difference of a few centimetres at most over 100 years.

2.5 Changes in wind climate

The wind-induced component w (Eq. 1) is related to sea level dynamics within the Baltic Sea basin. The narrow and shallow Danish Straits are the only ocean connection of the Baltic Sea, and water exchange in the Danish Straits plays an important role in Baltic sea level variations. The barotropic water exchange is governed by large-scale wind and air pressure patterns, which also affect the dynamical sea level within the Baltic Sea basin: persistent westerly winds push water through the straits and generate a slope in sea level, elevating sea levels at the eastern Baltic coast. Climatological changes in wind conditions can have a notable effect on long-term MSL at the Finnish coast (Johansson et al., 2014). We treat this component separately, as global climate models are too coarse to fully capture the water exchange in the Danish Straits and the dynamical sea level topography within the Baltic Sea.

Mean sea level in the Baltic Sea is positively correlated with the North Atlantic Oscillation (NAO) index, especially in winter (Johansson et al., 2001; Andersson, 2002). High NAO is related to a large north–south air pressure difference over the North Atlantic, which in turn is associated with strong westerlies. The NAO index is not the best metric to characterize the connection between regional atmospheric conditions and MSL at the Finnish coast, however. Johansson et al. (2014) found even better correlations with the zonal component of the geostrophic wind in the region. The best correlation was found with geostrophic winds calculated at the grid point 55^∘ N, 15^∘ E, close to the Danish Straits (Fig. 1). The zonal geostrophic wind u_g at this grid point explains more than 80 % of the interannual variability of sea level at the Finnish coast over the 20th century. Following Johansson et al. (2014), we calculate u_g from the Northern Hemisphere daily mean sea-level pressure fields for the years 1899–2018 (CISL Research Data Archive, 1979) according to the standard definition

$$\begin{array}{}\text{(3)}& u_{\mathrm{g}}=-{\displaystyle \frac{\mathrm{1}}{f\mathit{\rho }}}{\displaystyle \frac{\partial p}{\partial y}},\end{array}$$

where f is the latitude-dependent Coriolis parameter, ρ the air density (calculated assuming a constant temperature of 283 K), and $\partial p/\partial y$ the pressure gradient in the north–south direction.

To project u_g into the future, we study zonal geostrophic winds derived from the sea level pressure fields produced by the CMIP5 climate model runs (Taylor et al., 2012). We use results from 17 atmosphere–ocean general circulation models (AOGCMs), listed in Table A1. We select all models that (i) are assessed by Luomaranta et al. (2014) as capable of adequately simulating the observed climate (temperature, precipitation) in northern Europe and (ii) provide information of daily sea level pressure for all three emission scenarios.

From the daily mean sea level pressure fields of each model, we calculate u_g using Eq. (3) (but instead of constant temperature, the modelled pressure and 2 m temperature fields are used to calculate air density). These are calculated to the middle of the grid points of each model and interpolated to the actual model grid. For more details, see Ruosteenoja et al. (2019). From the daily geostrophic winds, we calculate 20-year annual means for the historical reference period (1986–2005) and for the end of the 21st century (2081–2100). Finally, the results are interpolated bilinearly to represent the point 55^∘ N, 15^∘ E.

The spread of the model-projected values of u_g for the historical period (1986–2005) is large (Table A1), the averages for this 20-year period ranging from 0.0 to 4.9 m s⁻¹ among the 17-model ensemble, while the observational average (calculated from the CISL dataset) is 3.2 m s⁻¹. To compensate for this large inter-model variability, we use a delta-change approach: instead of using the projected values for the future time periods as such, we calculate the model-specific changes from the historical period to the future period and add these to the observational historical average of 3.2 m s⁻¹. This still results in a large spread of the model projections, ranging from 2.4 to 5.7 m s⁻¹ for 2081–2100. The means, maxima and minima, as well as standard deviations among the 17-model ensemble for each emission scenario are given in Table 2.

Table 2 Projections of zonal geostrophic wind in 2081–2100 at 55^∘ N, 15^∘ E: means, maxima, minima, and standard deviations (metres per second) among the 17-model ensemble under low-emission (RCP2.6), medium-emission (RCP4.5), and high-emission (RCP8.5) scenarios. The values are calculated so that the projected change from the historical period 1986–2005 is added to the observational average of 3.2 m s⁻¹ over that period.

The wind-induced sea level component is then calculated as $w(i,t)=p_i{u}_{\mathrm{g}}\left(t\right)$, where p_i are regression coefficients between u_g and the detrended annual mean sea levels at the tide gauge i, calculated from 20th century observations. For more details, see Johansson et al. (2014). The regression coefficients range between 6.8 and 7.8 cm (m s⁻¹)⁻¹. The probability distribution for future w(i,t) is a product of two normal (Gaussian) distributions, whose standard deviations are the regression fit uncertainty of p_i and the standard deviation of the 17-model ensemble of u_g.

2.6 Past regional SLR trends

The land uplift rates given by the NKG2016LU model are solely based on geodetic observations. Therefore, they are independent of tide gauge observations and other sea-level-related information. This allows us to estimate the rate of SLR at the Finnish coast over the 20th century using the three-component model of Eq. (1), the unknown variable being the regional SLR r. For all Finnish tide gauges, we take the observed annual mean sea levels (h_m) and subtract the wind-induced component (w) and the land uplift (d) from the time series. The trend calculated from the residual corresponds to the local absolute SLR rate. Removing w from the annual mean sea levels reduces the year-to-year variability; therefore, r has narrower uncertainty ranges than the trend of h_m.

3 Results

3.1 Rate of past sea level rise

Past SLR trends for two time spans, t₀–2018 and 1993–2018, are shown in Table 3 for all Finnish tide gauges: t₀ is the beginning of the observational period and, depending on the site, varies between 1901–1933. The historical SLR trends (t₀–2018) at the Finnish coast range from 1.2 to 1.5 mm a⁻¹, in accordance with some recent estimates of the rate of the 20th century GMSLR (Hay et al., 2015; Dangendorf et al., 2017; Frederikse et al., 2020) but slightly lower than the AR6 estimate of 1.73 ± 0.45 mm a⁻¹ (Fox-Kemper et al., 2021). Local trends obtained for the satellite altimetry era, 1993–2018, range from 3.0 to 3.8 mm a⁻¹, in accordance with the global mean rate of 3.25 ± 0.37 mm a⁻¹ (Fox-Kemper et al., 2021).

Table 3 Historical trends (millimetres per year) of observed mean sea level h_m and its components: land uplift d, wind-induced component w, and regional sea level rise r. Error estimates are 1 SD. The trends are calculated from the beginning of the observations (t₀) to 2018; the trend of regional sea level rise r is also shown for 1993–2018.

3.2 Ensemble projection of GMSLR and the effect of weighting

We have combined 10 SLR projections published over the last decade (Table 1) to yield probability distributions of GMSLR for three emission scenarios. This analysis results in the following median estimates (5 % to 95 % ranges) of GMSLR in 2100 relative to 1995–2014: 45 (26 to 78) cm for the low-emission scenario, 60 (36 to 111) cm for the medium scenario, and 85 (54 to 176) cm for the high scenario. Figure 6 shows our probability distributions compared with the two IPCC AR6 projections, AR6MC and AR6LC, plotted as complementary cumulative distribution functions (1 − CDF) in semi-log scale so that differences in the upper tail are clearly visible. Interestingly, the combined distribution is nearly identical to AR6MC under the low-emission scenario but close to AR6LC under the medium scenario. Under the high scenario, the combined projection lies between AR6MC and AR6LC. The results are not drastically different even if the two AR6 projections are dropped out of the ensemble and we only use the eight other projections to calculate the combined projection.

Figure 6 Complementary cumulative distribution functions (1−F) of global mean sea level in 2100 relative to 1995–2014 under different emission scenarios. The red line denotes a weighted combination of 10 sea level rise projections, including the two IPCC AR6 projections, medium and low confidence (AR6MC and AR6LC, respectively). AR6MC and AR6LC are plotted separately for comparison.

As long as unresolved questions around ice sheet instability remain, sea level projections aiming to cover the whole probability range unavoidably rely to some extent on educated guesses, or expert judgement. Expert judgement is involved in the individual GMSLR projections (Table 1). In our method, subjectivity is involved in the choice of the type of distribution fitted in case there are multiple good fits, as well as in the weights given to the individual projections in the mixture (Eq. 2). Next, we examine how sensitive our results are to these choices.

As seen from Fig. A1, the three distribution types (Weibull, Fréchet, and skew normal) diverge substantially in the extreme tail of the probability distribution. In some cases, it is clear which distribution provides the best fit, but often any one of them could be chosen, while for Kopp et al. (2017) no good fit is found. We have tested different fitting strategies, such as always choosing the most pessimistic or the most optimistic option when multiple good fits are available, and the differences in the resulting mixture distribution are small (max. 2 cm in the 99.9th percentile). We conclude that our method is robust to the choice of the distribution type at least up to the cumulative probability of 99.9 %. Extrapolations beyond that are in any case very speculative, and we avoid presenting such results.

What matters more is the weight given to the various members of the projection ensemble when calculating the mixture, especially the weight given to projections that incorporate MICI. We experiment with different weighting strategies to see how much the results are affected by the weights:

1. All members are included with an equal weight.

2. AR6MC is given a 4-fold weight, as it can be considered to represent the latest scientific consensus on SLR projections.

3. Projections that incorporate MICI (Goodwin et al., 2017; Le Bars et al., 2017; and AR6LC) are given a lower weight (0.5), which reflects the low confidence in these projections. AR6MC is given a weight of 4 and the other members 1.

4. Projections that incorporate MICI are excluded from the ensemble, while AR6MC is given a weight of 4.

Increasing the weight of AR6MC somewhat lowers the final projection, but the difference is ca. 10 cm at most compared to the unweighted case. Unsurprisingly, lowering the weight of the MICI projections has a large effect on the upper tail of the distribution under higher-emission scenarios. The 99.9th percentile is lowered by 0.15, 0.5, and 0.9 m under the low, medium and high scenarios, respectively, if MICI projections are left out compared to the unweighted case. For the final results, we choose the fitting strategy 3, which represents a middle road – we want to take into account the MICI uncertainty but not give it too much weight because of the controversy. By adjusting the weights, the projections could be tailored to specific purposes according to the risk tolerance of the application.

3.3 Future mean sea level at the Finnish coast

After adding the land uplift and the wind-induced component related to the sea level dynamics within the Baltic Sea basin, we obtain probability distributions of MSL change at the Finnish coast by 2100. Figure 7 shows the distributions for three locations under the medium-emission scenario. The distributions of the three components are presented in Fig. 7b–d, while the distribution of the total MSL change in Fig. 7a is the sum (mathematically, convolution) of these components. As can be seen, land uplift (Fig. 7c) is the only component that shows a notable difference between locations. In reality, also the SLR distribution (Fig. 7b) varies slightly from place to place, but we apply the same distribution for all locations as explained in Sect. 2.3. Differences in the distribution of the wind-induced component are negligible between locations over the Finnish coast.

Figure 7 Probability distributions of projected mean sea level change (a) and its components (b–d) in Vaasa, Turku, and Helsinki under the medium-emission scenario (RCP4.5/SSP2-4.5), 2100 relative to 1995–2014. The distributions (a) and (b) are complementary cumulative distribution functions (1−F), while (c) and (d) are probability density functions (f).

Table 4 Projections of mean sea level change (centimetres) from 1995–2014 to 2100 for some of the largest coastal cities in Finland: Oulu, Vaasa, Pori, Turku, and Helsinki. The global sea level rise over the same period, calculated as a weighted combination of 10 SLR projections, is shown for comparison. The 5th, 50th, 95th, and 99th percentiles are given for three emission scenarios: low (RCP2.6/SSP1-2.6), medium (RCP4.5/SSP2-4.5), and high (RCP8.5/SSP5-8.5).

Projections of MSL change are tabulated for selected locations in Table 4, while projections for all Finnish tide gauges are provided in the Appendix (Tables A2, A3, A4). The projected MSL values are also given in the Finnish N2000 height reference system (Saaranen et al., 2009). Mean sea level over the reference period 1995–2014 was determined as the average of observed annual mean sea levels in N2000.

Under the low-emission scenario, substantial MSL rise is not expected at the Finnish coast over the 21st century: the median estimate ranges from a −43 cm decline in the Bothnian Bay to a 16 cm rise in the eastern Gulf of Finland (Table A2). However, according to the high-end (95 %) estimate, a rise of 50 cm is possible at the southern Finnish coast. Under the medium scenario, a sea level decline of −20 to −30 cm is expected in the Bothnian Bay, whereas the Gulf of Finland shows a MSL rise of 20 to 30 cm with an upper limit (95 %) of 70 to 80 cm. The high-emission scenario is the only one in which sea level rise is expected over the whole Finnish coastline, according to the median projection. The 95th percentiles are close to 1 m in the Bothnian Bay and 1.5 m in the Gulf of Finland.

Figure 8 shows the observed and projected MSL in Vaasa and Helsinki. The future time series are fitted to the 2100 projections assuming constant acceleration in the MSL trend. Vaasa is located in the region of strongest land uplift (8.8 mm a⁻¹), while in Helsinki the rate of land uplift is considerably weaker (3.7 mm a⁻¹). This is reflected in the past and future MSL trends. In Helsinki, the historical decline of mean sea level is about to reverse, and MSL rise is expected over the 21st century. In Vaasa, the decline will likely still continue for several decades before the turning point is reached by the end of the century. Uncertainties are large, however, as seen in the wide range of different percentiles plotted.

Figure 8 Observed annual mean sea levels and future projections under the medium-emission scenario (RCP4.5/SSP2-4.5) in Vaasa, region of strong land uplift, and Helsinki, where the land uplift is weaker. The reference level is the Finnish N2000 height reference system.

4 Discussion

4.1 Past mean sea level trends

The apparent sea level trend observed at a tide gauge site in Finland is a combination of three trends: land uplift, regional sea level rise, and the trend in the wind-induced component related to changes in the water balance of the semi-enclosed Baltic Sea. To illustrate the relative roles of all these components, we take Helsinki as an example.

The observed trend in annual mean sea levels at Helsinki was −1.88 (−0.37) mm a⁻¹ in 1904–2018 (1993–2018). This consists of land uplift of −3.69 mm a⁻¹, regional sea level rise of 1.32 (3.15) mm a⁻¹, and the trend in the wind-induced component of 0.49 (0.17) mm a⁻¹. Thus, in the recent past, the trends of SLR and land uplift have approximately balanced each other. In particular, we note that the trend in the wind-induced component is large enough that excluding it from the analysis would lead to a significant bias in the estimates of regional SLR rates.

Long-term absolute SLR trends (from the early 20th century to 2018) at the various locations on the Finnish coast vary from 1.22 to 1.48 mm a⁻¹ (Table 3). These are lower than the absolute SLR trends calculated for the Estonian coast by Suursaar and Kall (2018): 2.0 to 2.8 mm a⁻¹ in 1900–2010. The authors attribute the difference between their SLR rates and GMSLR to poorly resolved local land subsidence and the meteorologically driven internal sea level dynamics of the Baltic Sea. Over the satellite altimetry era, trends obtained in the present study (3.04 to 3.76 mm a⁻¹ in 1993–2018) are lower than the trends reported by Passaro et al. (2021) for the northern Baltic Sea based on satellite altimetry (up to 5 to 6 mm a⁻¹ in 1995–2019). Our SLR trend excludes the wind component (0.15 to 0.18 mm a⁻¹ in 1993–2018), which is included in the results of Passaro et al. (2021), but this does not explain the discrepancy comprehensively.

4.2 Mean sea level projections

We have presented probability distributions of mean sea level change over the 21st century along the Finnish coast for three scenarios representing different greenhouse gas emission pathways. Land uplift will continue at a constant rate, while acceleration in global sea level rise is expected together with a small increase in westerly winds in the northeast Atlantic, causing additional wind-induced sea level rise in the Baltic Sea. Due to accelerating SLR, a reversal of the historical declining MSL trend is expected at the southern Finnish coast. In the north, in the Bothnian Bay, postglacial land uplift is probably strong enough to compensate for SLR so that significant relative sea level rise is unlikely over this century, at least if the emissions of greenhouse gases are mitigated.

The uncertainty in the future MSL is high, both because of the unknown trajectory of emissions (reflected in the large differences in SLR outcomes between the emission scenarios) and the uncertainty associated with the SLR projections (reflected in the width of the probability distributions). Under the low-emission scenario, the MSL change would likely remain negative over most of the Finnish coastline during this century, the postglacial land uplift being strong enough to compensate for the SLR. According to the median projection, only the areas of weakest land uplift at the coast of the Gulf of Finland might experience a slightly rising MSL. Under the high scenario, substantial MSL rise is expected on the southern Finnish coast. Even in the areas of strongest land uplift, the historical decline in MSL may well be reversed and the MSL return close to present-day levels by the end of the century. The upper tails of the distributions point to a small but non-negligible risk of a considerably larger SLR, up to 1–2 m.

The previously published MSL projections for the Finnish coast (Pellikka et al., 2018), used in national adaptation planning, are close to the medium-emission projections presented in this study. The probability distribution calculated for the high-emission scenario shows potential for much stronger sea level rise, however. In the earlier estimates, the expected absolute SLR in Finland (excluding land uplift) was lower than the global mean due to the reduction caused by glacier fingerprints. According to current knowledge, this reduction may be cancelled over this century by larger-than-average contributions induced by thermal expansion and ocean dynamics.

Our mixture distribution of GMSLR is in line with the IPCC AR6 projection (Fox-Kemper et al., 2021, medium confidence) under the low-emission scenario, but higher than that in the upper tail under medium and high emissions. The medium-confidence projections of AR6 do not include processes related to marine ice sheets that are characterized by deep uncertainty. These are included in a separate low-confidence projection, which is similar to our mixture projection under medium emissions and even higher under high emissions. It can be concluded that AR6 manages to cover potential SLR futures better than the earlier IPCC reports, and including other projections in the mixture does not change the picture radically. In high-risk applications, it is important to consider the risk of SLR above the medium-confidence projection of AR6.

As can be seen in Fig. 7, the rate of SLR is the dominant source of uncertainty in the local MSL projections. The uncertainty intervals of the wind-induced and land uplift components are much narrower. The rate of land uplift simulated by the GIA model is effectively constrained by observations. Varying changes in the zonal geostrophic wind among the AOGCM simulations stem both from structural differences between the various models and internal variability in the climate system. For simplicity, we have assumed that inter-model differences in the u_g response follow the Gaussian distribution (Sect. 2.5). Owing to the small number of AOGCMs studied (N=17), we could neither confirm nor disprove the normality assumption; the multi-model frequency histograms of the u_g response for the various time spans and greenhouse gas scenarios proved to be significantly affected by stochastic variations. Nevertheless, the contribution of uncertainty in wind-induced forcing and, in particular, land uplift to the total local MSL change is far smaller than the uncertainty originating from the global sea level rise. Evidently, the exact shapes of the probability distributions of these factors have a negligible impact on the total uncertainty.

4.3 Feasibility of the emission scenarios

The low-emission scenario (RCP2.6/SSP1-2.6) represents ambitious climate mitigation in line with the Paris Agreement: global warming stays below 2 ^∘C as the emissions decline to net zero in the second half of the 21st century (Chen et al., 2021). In the medium scenario (RCP4.5/SSP2-4.5), emissions peak around mid-century, resulting in an estimated global warming of 2.7 ^∘C by 2100. The high scenario (RCP8.5/SSP5-8.5) represents a fossil-fuel intensive pathway, where emissions grow strongly during this century and global warming by 2100 can exceed 4 ^∘C.

IPCC does not assess the likelihood or feasibility of the different emission scenarios, but such estimates are often needed by policymakers, and implicit assessments of scenario probabilities may well be made by end users in the absence of relevant expert assessments (Ho et al., 2019). Hausfather and Peters (2020) argue that considering current trends in the energy sector, SSP5-8.5 should be regarded as a highly unlikely worst-case scenario and that current policies point to the modest mitigation pathway SSP2-4.5. On the other hand, the feasibility of the low-emission scenario can also be questioned, as the ambition in climate change mitigation should strongly increase to keep the goals of the Paris Agreement – and hence SSP1-2.6 – within reach (den Elzen et al., 2022). This again points to SSP2-4.5 as the most likely greenhouse gas scenario.

4.4 Notes for the end users of the projections

The projections presented in this study are intended for use in many different applications where information about future MSL is needed. The projections are presented as probability distributions so that they can be widely applied for different purposes. It is up to the end users to choose the projections appropriately so that the probability level of the projection is in line with the acceptable risk level, which depends on the application.

Tables A2–A4 present median projections with three different uncertainty ranges. The narrowest range, 17 %–83 %, corresponds to the likely range used by IPCC in their SLR projections. The widest range, 1 %–99 %, is intended to be used in high-risk applications. Even higher projections could be presented with an extremely small probability; it is impossible to determine a definite upper limit for sea level rise.

The end users of all probabilistic SLR projections need to recognize that the probabilities are always to some extent subjective and dependent on different assumptions. As long as there is deep uncertainty related to ice sheet behaviour in the warming climate, it is recommended not to be overconfident in any seemingly precise SLR projection (Behar et al., 2017). Uncertainties are large, especially in the upper end; median projections are better constrained. The upper-end projections should be viewed as illustrative, reflecting the best scientific understanding at the moment.

In many applications, the short-term sea level maxima (flooding) are of high relevance, not just the mean sea level. At the Finnish coast, short-term maxima vary between 1–2 m above the mean. In Helsinki, for example, the worst coastal flood on record, 1.51 m above the mean sea level, took place in January 2005 during the storm Gudrun and caused significant damage, which was alleviated with temporary flood protection structures. Even under the medium-emission scenario, coastal floods of that magnitude may become common by the end of the century (Pellikka et al., 2018), which would require substantial investments in adaptation.

In the future, the height of coastal floods will increase simply because the mean sea level is rising. In addition, it is possible that climate change affects the weather phenomena causing short-term sea level variations, most importantly through changes in wind and air pressure conditions. The change in water depth as the mean sea level rises also impacts the dynamics of short-term sea level phenomena, such as seiche oscillations, tides, and wind waves. The effect of this interaction between the mean sea level and short-term variability can be significant especially in shallow water areas (Idier et al., 2019), such as the Baltic Sea, but is not very well understood. A more detailed study of short-term sea level variability is beyond the scope of this paper.

5 Conclusions

The main findings of this study can be summarized as follows:

1. Historical trends of absolute sea level rise at the Finnish coast, excluding the effect of land uplift and wind-induced changes in Baltic sea level, are in accordance with global mean rates. This evolution is expected to continue over the 21st century, as well, as the smaller-than-average glacier contribution is compensated by the larger-than-average ocean contribution in the Baltic Sea.

2. Median projections of mean sea level change in Finland (2100 relative to 1995–2014) range from −43 to +16 cm for the low-emission scenario, −28 to +31 cm for the medium scenario, and +1 to +61 cm for the high scenario depending on location. The lowest values occur in the area of strong land uplift in the north and the highest ones in the Gulf of Finland in the south. The upper tail of the probability distribution characterizes the risk of higher sea level rise. The 95th percentiles range from −9 to +50 cm (low emissions), +24 to +83 cm (medium), and +93 to +152 cm (high) depending on location.

3. Our weighted combination of 10 individual GMSLR projections is close to or in between of the two probability distributions provided by IPCC AR6 (Fox-Kemper et al., 2021, medium and low confidence). It seems that AR6 manages to cover the full spectrum of potential sea level rise by 2100 better than earlier IPCC reports. The low-confidence storyline should be considered in high-risk applications, such as in assessing nuclear power plant safety.

Appendix A: Additional results

Here, we provide additional tables and figures that might be beneficial for the end users of this study. Table A1 shows the projections of zonal geostrophic wind from the individual CMIP5 models used in this study. Tables A2, A3, and A4 give the MSL projections for all Finnish tide gauges and different emission scenarios. Figure A1 shows the individual probability distributions that are used to compose the mixture distribution of GMSLR.

Table A1 Twenty-year means of zonal geostrophic wind (metres per second) at 55^∘ N, 15^∘ E, according to CMIP5 model results: historical (1986–2005) and future (2081–2100) under low-emission (RCP2.6), medium-emission (RCP4.5), and high-emission (RCP8.5) scenarios.

Table A2 Projected mean sea level (centimetres) at the Finnish tide gauges in 2100 under the low-emission scenario (RCP2.6/SSP1-2.6).

Table A3 Projected mean sea level (centimetres) at the Finnish tide gauges in 2100 under the medium-emission scenario (RCP4.5/SSP2-4.5).

Table A4 Projected mean sea level (centimetres) at the Finnish tide gauges in 2100 under the high-emission scenario (RCP8.5/SSP5-8.5).

Figure A1 Members of the projection ensemble under the high-emission scenario. Fits for all three distribution types (Weibull, Fréchet, skew normal) are shown as complementary cumulative distribution functions (1−F). No fits were applied to Jevrejeva et al. (2014) and IPCC AR6 (Fox-Kemper et al., 2021). AR6MC and AR6LC are the medium- and low-confidence projections of IPCC AR6, respectively.