Long-Term Changes, Inter-Annual, and Monthly Variability of Sea Level at the Coasts of the Spanish Mediterranean and the Gulf of Cádiz

Abstract

One of the effects of climate change is the rise of sea level, which poses an important threat to coastal areas. Therefore, the protection and management of coastal ecosystems as well as human infrastructures and constructions require an accurate knowledge of those changes occurring at a local scale. In this study, long time series of sea level from tide gauges distributed along the southern (Atlantic) and eastern (Mediterranean) Spanish coasts were analyzed. Linear trends were calculated for two periods, from early 1940s to 2018 and from 1990 to 2018. Values for the former period ranged between 0.68 and 1.22 mm/year. These trends experienced a significant increase for the second period, when they ranged between 1.5 and 4.6 mm/year. Previous research analyzed the effect of atmospheric forcing in the Mediterranean Sea by means of 2D numerical models, and the steric contribution was directly evaluated by the integration of density along the water column. In this study, the effect of atmospheric forcing and the thermosteric and halosteric contributions on coastal sea level were empirically determined by means of statistical linear models that established which factors affected sea level at each location and what the numerical response of the observed sea level was to the contributing factors. Atmospheric pressure and the west–east component of the wind hada significant contribution to the sea level variability at most of the tide gauges. The thermosteric and halosteric components of sea level also contributed to the sea level variability at all the tide gauges, with the only exception of Alicante. Atmospheric forcing and the steric components of sea level experienced long-term trends. The combination of such trends, with the response of sea level to these factors, allowed us to estimate their contribution to the observed sea level trends. The part of these trends not explained by the atmospheric variables and the steric contributions was attributed to mass addition. Trends associated with mass addition ranged between 0.6 and 1.2 mm/year for the period 1948–2018 and between 1.0 and 4.5 mm/year for the period 1990–2018.

1. Introduction

Sea level rise is one of the direct and clear effects of climate change [1]. At a regional scale, variations in sea level can be caused by different factors, such as meteorological forcing (atmospheric pressure, wind, and heat fluxes) or changes in the ocean circulation. Nevertheless, when sea level is considered at a global scale, the mass addition linked to the melting of continental glaciers [2,3] and the thermal expansion of seawater are the main contributions. Several efforts have been devoted to the reconstruction of global mean sea level using time series obtained from tide gauges and spatial patterns of variability inferred from satellite altimetry data (Reduced Space Optimal Interpolation, RSOI; [4,5,6,7,8]). Such reconstructions have shown that global mean sea level has increased during the twentieth century, although with a strong decadal variability [9]. During the twentieth century, global sea level rose at a rate of 1.1 mm/year until 1990. This trend increased to 3.1 mm/year since 1993 [5,7]. However, the changes observed for the global ocean do not necessarily reflect those changes operating at regional scales, as the former are dominated by the signal of the Pacific Ocean [10] (Gomis et al., 2012). Hence, specific studies need to be conducted in order to understand the response of each region of the world ocean to climate change as well as to understand the influence of other forcing factors, such as meteorological or circulation variability.

One clear example of a differentiated regional response is the case of the Mediterranean Sea. Several studies have reconstructed the Mediterranean Sea level during the twentieth century and at the beginning of the twenty-first century using similar techniques as those described for the global ocean (RSOI, [11,12,13]). These studies have evidenced that from 1960 to the beginning of the 1990s decade, the sea level in the Mediterranean increased at a much lower rate than that observed forthe global ocean, or that it even decreased in some parts of this region [12,14]. The explanation for this anomalous behavior seemed to be the increase in atmospheric pressure over the Mediterranean Sea during this period [14,15] and a possible decrease of the halosteric contribution [16]. Nevertheless, the sea level trends in the Mediterranean Sea would have recovered since the beginning of the 1990s decade [14,17,18].

Although some models cast doubts about the future evolution of the salinity of the Mediterranean Sea [19], most of the atmospheric ocean global circulation models (AOGCMs) predict an increase in the sea level in the Mediterranean during the twenty-first century, with an acceleration of current trends caused by the mass addition of freshwater and the contribution of the steric component [20,21,22]. These projections pose an important threat to coastal areas, whose protection requires the accurate assessment of changes occurring at a local scale. However, even those changes observed at the regional scale could not be representative of the variability at coastal locations [4,10]. The main goal of the present study was to analyze the inter-annual and long-term variability of local sea level along the southern and eastern coasts of Spain and to assess the contributions to such variability of atmospheric forcing, both atmospheric pressure and wind velocity, and the thermosteric and halosteric components of sea level. Previous research has analyzed atmospheric forcing by means of barotropic 2D models [14,23,24]. The thermosteric and halosteric components of sea level have been estimated using the available temperature and salinity profiles collected in different databases [23,25,26]. In the present study, a different and complementary approach was followed. Atmospheric pressure and Cartesian components of the wind were obtained from reanalysis data. The thermosteric and halosteric contributions of sea level were estimated from available temperature and salinity compilations. Then, a linear statistical model was fitted to the data for determining which factors really affect sea level on a local scale in a significant way. Then, the numerical response of sea level to such factors was obtained for each of the analyzed locations. These relationships between local sea level and the atmospheric and steric components were used to infer their contribution to the observed long-term trends and the sea level rise caused by mass addition.

2. Data and Methods

2.1. Dataset

We used a monthly averaged time series of sea level from the Revised Local Reference (RLR) dataset from the Permanent Service for Mean Sea Level (PSMSL; [27]). Only those data considered as reliable by the PSMSL were used. Initially, data from 31 tide gauge stations covering the Spanish Mediterranean (including the Balearic Islands) and the Gulf of Cádiz were obtained. As the main goal of this study was to analyze inter-annual and long-term changes in sea level, only those time series longer than 30 years were finally considered for their analyses. In addition, some stations with a shorter length were used for completing the gaps present in the long time series. Other time series were not considered for their analysis or for completing the long-time series but were used for comparison and validation of the reconstructed time series. The dots in Figure 1 show all the tide gauges initially considered. For some locations, more than one tide gauge is operated by different institutions. The numbers of available time series are shown in parentheses. The final set of tide gauges used for the analysis, for filling the gaps in the long time series, or for comparison are shown in

Table 1. In bold, tide gauge stations with long time series used for the analysis of inter-annual and long-term variability. In normal font are the names of those time series used for filling the gaps in the long time series. Shaded squares correspond to those years with fewer than three missing monthly values. The interruption in the years is due to the lack of data from 1939 to 1943.

	Cádiz II	Cádiz III	Tarifa	Tarifa 2	Algeciras	Algeciras 2	Algeciras b	Ceuta	Málaga	Málaga II	Alicante in	Alicante out	Valencia	L’Estartit
1880–1924
1927–1938
1939–1943
1944
1960
1975
1990
1995
2000
2010
2015

and are circled in Figure 1. The names used to identify the tide gauge stations are those used in the PSMSL (www.psmsl.org accessed on 2 December 2020).

In the case of sea level at Cádiz, Marcos et al. (2011, [28]) conducted a very exhaustive study of data archeology, collecting very complete information from different sources, including logbooks. These authors found that, because of leveling problems, an offset of 37.5 mm should be added to the data corresponding to the period 1880–1924 in order to make them comparable with modern data (after 1961). This correction was applied in the present study. Concerning sea level time series from the tide gauges in the outer and inner harbors of Alicante, Marcos et al. (2021, [29]) reconstructed those time series considering historical information and leveling problems. The two mean sea level time series reconstructed by Marcos et al. (2021) are more complete and accurate than those available from the PSMSL. These time series are freely available and were downloaded from the British Oceanographic Data Centre (BODC, www.bodc.ac.uk accessed on 9 July 2021).

2.2. Reconstruction of Long Time Series

Following Calafat et al. (2014), gaps shorter than three months were filled by means of natural cubic splines. Then the seasonal cycle was obtained for each time series, calculating the mean value of each calendar month (Figure 2). These seasonal cycles were estimated using the complete time series available for each tide gauge (see Table 1). The climatological seasonal cycles were then subtracted from the original sea level time series (see Supplementary Materials for details). These series of anomalies or residuals contained gaps longer than three months. Such gaps were filled using linear regression on nearby locations (see Supplementary Materials and [30] for the details). The sea level time series at Algeciras had a gap from 2002 to 2018. There are two more tide gauges located at the same location (Algeciras 2 and Algeciras b in Table 1), but they have been operating since 2007 and 2010, and therefore the Algeciras time series could not be regressed on them, as they do not share a common period of time. For this reason, Algeciras data were filled by means of linear regression on Cádiz, Ceuta, Tarifa, and Málaga. Then the reconstructed time series was correlated with those from the redundant tide gauges operating at Algeciras, offering a means to test the procedure used for filling the gaps. The L’Estartit sea level time series had no gaps. After visual inspection, the year 2018 seemed to be anomalous, as data seemed to reproduce exactly the climatological seasonal cycle. Consequently, this year was substituted by its regression on Barcelona sea level. Statistics from the linear regression used for filling gaps can be found in Supplementary Materials. Figure 3 shows the original time series (black lines) and the reconstructions (red lines) for the case of Cádiz and Tarifa. Figure 4 shows similar results for Algeciras and Ceuta.The reconstructed time series for Algeciras was correlated with the sea level at Algeciras 2 and Algeciras b (Figure 4C,D). The correlation coefficients were 0.88 and 0.85, respectively. Figure 5 shows the original and reconstructed time series for Málaga, Alicante (outer harbor), and L’Estartit.

Each sea level time series has a different starting date, but all of them end in December 2018.

2.3. Pressure and Wind Data

Monthly atmospheric pressure (P hereafter) and wind data were obtained from the reanalysis of the NCEP/NCAR [31]. Both the west–east component (U hereafter) and south–north component (V hereafter) of the wind were obtained. U is positive when directed eastwards and V is positive when directed northwards. This dataset has a resolution of 2.5° in longitude and latitude and extends from January 1948. The grid point closest to each tide gauge was selected and monthly time series of P(hPa), U (m/s), and V(m/s) were obtained from January 1948 to December 2018. The distance between the three tide gauges located in the Strait of Gibraltar (Tarifa, Algeciras, and Ceuta; Figure 1) is lower than the spatial resolution of atmospheric time series. For this reason, common time series of P, U, and V were calculated for the Strait of Gibraltar and were considered as representative of meteorological conditions in Tarifa, Algeciras, and Ceuta.

A seasonal cycle was calculated for each atmospheric time series, estimating the mean value for each calendar month. These seasonal cycles were subtracted from the original time series, producing time series of anomalies or residuals (see Supplementary Materials).

2.4. Steric Component of Sea Level

Temperature and salinity gridded fields were downloaded from the Met Office Hadley Centre observations datasets (version EN.4.2.1; [32]). This product offers monthly temperature and salinity profiles on a 1° × 1° grid at 42 vertical levels. An area adjacent to each tide gauge was selected and all the profiles within that area were averaged in order to obtain monthly time series of temperature and salinity profiles representative of the conditions close to Cadiz, Málaga, Alicante, and L’Estartit (see red boxes in Figure 1). In the case of Tarifa, Algeciras, and Ceuta, the spatial resolution of the gridded fields did not allow us to obtain temperature and salinity profiles representative of each location. An area corresponding to the Strait of Gibraltar conditions was selected (see Figure 1). These time series extend from January 1940 to December 2019. As sea level time series end on December 2018, this was also the final date of the temperature and salinity dataset.

Following Jordá and Gomis (2013a, 2013b, [25,26]), the observed sea level (corrected for land movements) can be decomposed into a mass and a steric component. First, mean profiles of temperature, salinity, and density were calculated for the complete period 1940–2018. These mean profiles were subtracted from each monthly profile, and the variations of sea level caused by the steric component were estimated as:

$$\Delta {\eta }_{St}=-\frac{1}{{\rho }_s}{\displaystyle \underset{-H}{\overset{\eta }{\int }}\Delta \rho (z)dz}$$

(1)

where ρ is the density of sea water, and Δρ is the deviation of the monthly density profile from the long-term mean value estimated for the complete period of time (1940–2018). ρ_S is the density for the sea surface (vertical coordinate z = η). The steric contribution was decomposed into its thermosteric (Δη_T) and halosteric components (Δη_H):

$$\Delta {\eta }_T=-\frac{1}{{\rho }_S}{\displaystyle \underset{-H}{\overset{\eta }{\int }}{\rho }_0\alpha \Delta T(z)dz}$$

(2)

$$\Delta {\eta }_H=-\frac{1}{{\rho }_S}{\displaystyle \underset{-H}{\overset{\eta }{\int }}{\rho }_0\beta \Delta S(z)dz}$$

(3)

where ρ₀ is the reference density (mean value of the complete series). The time series of steric, thermosteric, and halosteric contributions contain a seasonal cycle. Once again these cycles were estimated calculating the mean value for each calendar month (see Figure 2). The seasonal cycles were subtracted from the original time series to obtain time series of anomalies or residuals.

In addition to the EN4 product, gridded temperature and salinity data were obtained from the NCAR/UCAR Research Data Archive [33,34]. This dataset extends from January 1954 to December 2012, not covering the complete period analyzed in the present work. For this reason, it was only used to check the sensitivity of our results to the dataset used.

2.5. Estimation of Linear Trends

Time series of anomalies of sea level; P, U, and V components of the wind; and the steric component of sea level can be affected by long-term changes with time scales comparable to the length of the analyzed period. Such changes are usually modeled as a straight line fitted to the data by means of least squares. The slope of the fitted line can be interpreted as a mean rate of change. The importance of trend estimation is twofold. First, as already stated, it expresses the rate of change of sea level and of other factors that could be involved in sea level variability. Second, the existence of linear trends in the sea level and in any of the factors contributing to its variability does not allow us to establish any relationship between both trends. These relationships (if existing) should be extracted from de-seasoned and de-trended time series.

The confidence intervals for the trends (slope of the straight lines) were calculated using Student’s t distribution. Nevertheless, the monthly sea level values were auto-correlated. For this reason, the equations used were corrected for the loss of degrees of freedom caused by this autocorrelation [35,36].

Several periods of time were considered for the estimation of linear trends. As most of the analyzed time series started during the 1940s decade, first, linear trends were estimated for the complete period covered by such time series. These time series and the exact period of time are included in

Table 2. Linear trends for the time series of de-seasoned sea level, atmospheric pressure, and the U, V components of the wind. Those trends statistically significant at the 0.05 significance level are typed in bold.

	Relative Sea Level		P Atmospheric Pressure		U-Wind	V-Wind
	Period	b ± CI	Period	b ± CI	b ± CI	b ± CI
		(mm/year)		(hPa/year)	(m s−1/year)	(m s−1/year)
Cádiz GIA −0.23	1880–2018	1.28 ± 0.08
	1943–2018	1.08 ± 0.17	1948–2018	0.02 ± 0.01	−0.01 ± 0.01	0
	1990–2018	3.6 ± 0.6	1990–2018	−0.02 ± 0.03	0.03 ± 0.02	−0.03 ± 0.02
Tarifa GIA −0.22	1943–2018	1.22 ± 0.19
	1990–2018	4.2 ± 0.5
Algeciras GIA −0.21	1943–2018	0.89 ± 0.13	1948–2018	0.02 ± 0.01	0	0
	1990–2018	2.6 ± 0.5	1990–2018	−0.01 ± 0.03	0.03 ± 0.02	0.01 ± 0.01
Ceuta GIA −0.21	1944–2018	0.68 ± 0.15
	1990–2018	1.5 ± 0.5
Málaga GIA −0.17	1943–2018	1.08 ± 0.18	1948–2018	0.02 ± 0.01	0	0
	1990–2018	4.6 ± 0.6	1990–2018	−0.03 ± 0.03	0.04 ± 0.02	0
Alicante GIA −0.31	1914–2018	0.92 ± 0.09	1948–2018	0.02 ± 0.01	0	−0.01 ± 0.00
	1948–2018	0.80 ± 0.17	1990–2018	−0.03 ± 0.03	0.01 ± 0.02	0.03 ± 0.02
	1990–2018	2.2 ± 0.6
L’Estartit			1948–2018	0.02 ± 0.01	0	−0.01 ± 0.01
GIA −0.38	1990–2018	2.8 ± 0.7	1990–2018	−0.04 ± 0.04	0.01 ± 0.01	−0.02 ± 0.02

. As the Cadiz sea level time series started in 1880, the linear trends were estimated from 1880 to 2018 and also from 1943 to 2018, for comparison with the rest of the time series. Finally, the time series of sea level at L’Estartit started in 1990, and consequently, the linear trend could only be calculated for the period 1990–2018. This date is also coincident with an increase in sea level trends in the Mediterranean Sea [12,14,24] and in the global ocean [1,7,37] (Oppenheimer et al., 2019; Dangerdorf et al., 2017; Church and White, 2011). For this reason, and for comparison with the other locations analyzed in this work, all the trends were also computed for the period 1990–2018 (see Table 2).

2.6. Statistical Model

To explore the relationship between the observed sea level at each location and the atmospheric variables (P, U, V), and the thermosteric and halosteric contributions, a linear model was used. Instead of imposing the variables that explain the variability of sea level (de-seasoned and de-trended), an explanatory model was selected by means of a forward stepwise linear regression (see for instance [30] or Supplementary Materials for the details). The variance explained by the selected model and the coefficients corresponding to each explanatory variable or predictor were calculated. This procedure was applied for two different periods: 1948–2018 and 1990–2018.

3. Results

3.1. Seasonal Cycles

The observed sea level seasonal cycle showed maximum values in October and November in all the tide gauges analyzed. The steric contribution of sea level also had a maximum value in October, with very similar values in November at Cádiz, at the Strait of Gibraltar, and at L’Estartit station. In the cases of Málaga and Alicante, the maximum value of the steric contribution was observed in September, although the October mean value was very close to the September one. No differences were observed between the steric components calculated using the EN4 or NCAR/UCAR temperature and salinity datasets.

The seasonal cycle of the atmospheric pressure is characterized by maximum values in January, then a decrease with a minimum value in April, and a secondary minimum in August. The atmospheric pressure increases during the rest of the year until December.

3.2. Linear Trends

Table 2 shows the sea level trends for the reconstructed time series. Rates of vertical land movement caused by Global Isostatic Adjustment are included in Table 2 for the position of each tide gauge. These values were obtained from those listed in the PSMSL according to the ICE-5G v1.3 ice model [38]. As already explained in the Data and Methods section, sea level trends were calculated for the complete length of the time series and for the period 1990–2018. In all the cases, sea level showed positive trends. Nevertheless, the values obtained are very sensitive to the chosen period of time. Trends estimated from the beginning of the 1940s decade were much lower than those obtained for the 1990–2018 period. In the case of Cadiz sea level, linear trends estimated from 1880 to 2018 were lower than those corresponding to 1990–2018, but higher than those estimated from 1943. These results simply indicate that, although long-term changes reflect a clear increase in sea level, this behavior is modulated by a strong decadal variability (see the Discussion section for more details).

Despite the observed changes in the sea level trends, all the results obtained show positive values. On the contrary, trends obtained for atmospheric pressure, and the U, V components of the wind show positive or negative values depending on the period analyzed.

Table 3. Linear trends for the steric, thermosteric, and halosteric contributions of sea level estimated from EN4 gridded temperature and salinity profiles. Those trends statistically significant at the 0.05 significance level are typed in bold.

		Steric Component	Thermosteric Component	Halosteric Component
	Period	b ± CI (mm/year)	b ± CI (mm/year)	b ± CI (mm/year)
Cádiz	1940–2018 1990–2018	0.48 ± 0.17 3.2 ± 0.8	1.12 ± 0.10 0.9 ± 0.5	−0.64 ± 0.20 2.2 ± 0.9
Strait of Gibraltar	1940–2018 1990–2018	0.47 ± 0.17 3.4 ± 0.8	1.15 ± 0.11 2.2± 0.6	−0.68 ± 0.21 1.2 ± 0.9
Málaga	1940–2018 1990–2018	0.16 ± 0.16 1.9 ± 0.8	1.15 ± 0.11 2.0 ± 0.6	−0.99 ± 0.20 −0.1 ± 0.9
Alicante	1940–2018 1990–2018	−0.54 ± 0.12 −1.4 ± 0.6	1.12 ± 0.09 2.7 ± 0.4	−1.63 ± 0.12 −4.0 ± 0.5
L’Estartit	1940–2018 1990–2018	−0.82 ± 0.10 −2.7 ± 0.4	1.03 ± 0.08 2.8 ± 0.3	−1.84 ± 0.11 −5.4 ± 0.4

shows the linear trends for the steric, thermosteric, and halosteric components of sea level estimated from the EN4 gridded temperature and salinity fields. Calculations were repeated for the periods 1940–2018 and 1990–2018. These trends were also estimated using the NCAR/UCAR dataset, but no significant changes were observed. The steric contribution to sea level was positive at Cádiz, the Strait of Gibraltar and Málaga, and negative for the Alicante and L’Estartit regions. The negative values of the steric contributions in Alicante and L’Estartit were caused by a strong and negative halosteric contribution for the two periods analyzed. The halosteric component was negative for Cádiz, Gibraltar, and Málaga for the period 1940–2018, but in these cases it was compensated by the thermosteric contribution. Both the halosteric and thermosteric components were positive for the period 1990–2018 at these three locations (see Table 3).

3.3. Inter-Annual Variability—Linear Model

In order to explain the factors that produce the monthly and inter-annual variability of sea level, their linear dependence on the atmospheric variables and the thermosteric and halosteric contributions was assumed. The variables (predictors) that should be included in the model were not considered a priori, and a forward stepwise procedure selected those variables that explained most of the observed variance.

Table 4. Linear model selected by the forward stepwise multiple linear regression for the de-seasoned and de-trended sea level at Cádiz, Tarifa, Algeciras, Ceuta, Málaga, Alicante, and L’Estartit. The multiple correlation coefficient (root square of the explained variance) and the 95% confidence intervals for the coefficients in the model are also shown.

	1948–2018	1990–2018
Cádiz	R = 0.66 η = b0 + b1P + b2U + b5ηH b1 = −13.76 ± 1.14 mm/hPa b2 = −8.4 ± 1.8 mm/(m s−1) b5 = 0.14 ± 0.04	R = 0.79 η = b0 + b1P + b2U + b5ηH b1 = −14.5 ± 1.3 mm/hPa b2 = −8.6 ± 2.2 mm/(m s−1) b5 = 0.10 ± 0.05
Tarifa	R = 0.57 η = b0 + b1P + b2U + b3V + b4ηT + b5ηH b1 = −13.0 ± 1.9 mm/hPa b2 = −8.3 ± 2.1 mm/(m s−1) b3 = 15 ± 7 mm/(m s−1) b4 = 0.23 ± 0.10 b5 = 0.22 ± 0.05	R = 0.79 η = b0 + b1P + b2U + b3V + b4ηT + b5ηH b1 = −11.9 ± 2.1 mm/hPa b2 = −7.0 ± 2.0 mm/(m s−1) b3 = 28 ± 8 mm/(m s−1) b4 = 0.22 ± 0.08 b5 = 0.15 ± 0.05
Algeciras	R = 0.69 η = b0 + b1P + b2U + b3V + b4ηT + b5ηH b1 = −12.3± 1.2 mm/hPa b2 = −7.5 ± 1.3 mm/(m s−1) b3 = 15.2 ± 4.2 mm/(m s−1) b4 = 0.09 ± 0.07 b5 = 0.11 ± 0.03	R = 0.80 η = b0 + b1P + b2U + b3V + b4ηT + b5ηH b1 = −11.0 ± 1.6 mm/hPa b2 = −5.9 ± 1.6 mm/(m s−1) b3 = 24 ± 6 mm/(m s−1) b4 = 0.11 ± 0.07 b5 = 0.13 ± 0.04
Ceuta	R = 0.67 η = b0 + b1P + b3V + b4ηT b1 = −12.8 ± 1.3 mm/hPa b3 = 13.3 ± 4.6 mm/(m s−1) b4 = 0.14 ± 0.06	R = 0.78 η = b0 + b1P + b3V + b4ηT b1 = −12.4 ± 1.6 mm/hPa b3 = 16 ± 6 mm/(m s−1) b4 = 0.09 ± 0.06
Málaga	R = 0.65 η = b0 + b1P + b2U + b5ηH b1 = −14.5 ± 1.3 mm/hPa b2 = −12.4 ± 1.8 mm/(m s−1) b5 = 0.18 ± 0.05	R = 0.70 η = b0 + b1P + b2U + b5ηH b1 = −15.1 ± 2.1 mm/hPa b2 = −11 ± 3 mm/(m s−1) b5 = 0.22 ± 0.08
Alicante	R = 0.72 η = b0 + b1P + b2U + b3V b1 = −13.7 ± 1.1 mm/hPa b2 = −6.4 ± 1.9 mm/(m s−1) b3 = 6 ± 3 mm/(m s−1)	R = 0.61 η = b0 + b1P + b2U + b3V b1 = −11.2 ± 2.0 mm/hPa b2 = −6 ± 4 mm/(m s−1) b3 = 8 ± 6 mm/(m s−1)
L’Estartit		R = 0.87 η = b0 + b1P + b2U + b3V + b4ηT + b5ηH b1 = −13.0 ± 0.9 mm/hPa b2 = −7 ± 3 mm/(m s−1) b3 = 4.8 ± 2.0 mm/(m s−1) b4 = 0.17 ± 0.12 b5 = 0.20 ± 0.07

shows the selected model for each location, the multiple correlation coefficients, and the values of the coefficients that indicate the response of sea level to each predictor. The linear model was fitted for the period 1948–2018 (variables P, U, and V limit the extension of the model), and for the period 1990–2018. Figure 6 shows the time series of sea level anomalies (black lines) and the regression using the models and the predictors in Table 4 for the locations of Cádiz, Tarifa, and Algeciras (red lines). These results are presented for the two periods analyzed. In the case of the long period 1948–2018, annual mean values are presented for the clarity of the plot, whereas monthly values are used for the period 1990–2018. Figure 7 shows the same results for the locations of Ceuta, Málaga, and Alicante. Finally, Figure 8 shows the results corresponding to the L’Estartit time series. In this case, the only available period was 1990–2018. Figure 8A shows the time series of annual mean values (observed and regressed values), Figure 8B shows the same results for the monthly time series, and Figure 8C is a zoom for the period 2010–2018 to highlight the good agreement between the observations and the model.

4. Discussion and Summary

Several studie shave focused on sea level rise at a global scale [5,9,37]. Nevertheless, adequate assessment for coastal protection requires knowledge of sea level variability at a local scale. It should be noticed that important differences between global and local sea level could be caused by several factors, such as atmospheric forcing, mesoscale circulation, river runoff, etc. [1,4,10]. Previous research has addressed the problem of reconstructing global sea level backward in time from coastal tide gauges, using spatial patterns of variability inferred from satellite altimetry data. This methodology has been applied both at a global scale [4,5,39] and at a regional one [11,12,13] for the case of the Mediterranean Sea). In this study, we followed a different approach and, instead of global or regional sea level, local time series were considered. The reconstruction of a time series that was as long as possible was addressed using multiple linear regression on nearby tide gauges and those reconstructions carried out by Marcos et al. (2021, 2011, [28,29]). The reconstructed time series at Algeciras was significantly correlated with sea level time series from two tide gauges at the same location (correlation coefficients 0.88 and 0.85), giving us confidence in the reliability of the final long time series. Time series should be de-seasoned and de-trended prior to regression analysis on any potential explanatory variable. This was the reason for estimating the seasonal cycles for all the sea level time series as well as for the atmospheric variables and the steric contributions. Although the analysis of these cycles was not the objective of this work, it should be noticed that our results are coincident with those obtained by García-Lafuente et al. (2004, [40]). The maximum values of sea level are observed in October and November, coinciding with the maximum values of the steric contribution, which reflects the seasonal cycle of the heat exchange between the atmosphere and the sea. García-Lafuente et al.(2004, [40]) already found that most of the variance associated with the sea level seasonal cycle was explained by the steric component, whereas the atmospheric variables accounted for the semi-annual cycle of sea level.

Our motivation was different in the case of linear trends. This long-term variability should be removed before addressing the linear regression that explains the monthly and inter-annual variability of local sea level. However, the estimation of such linear trends was one of the main goals of this study. The local sea level increased at all the tide gauges analyzed. The linear trends were significantly lower when they were estimated from the beginning of the 1940s than when they were estimated for the period 1990–2018. In the case of Cádiz, where data from 1880 were available, the trend for the period 1880–2018 showed an intermediate value between those corresponding to the other two periods. However, it should be noticed that the trend obtained in the present work for the period 1880–2018 (1.28 ± 0.08 mm/year) was much higher than that presented in Marcos et al. (2011); 0.7 ± 0.1 mm/year). In the case of the Alicante time series, we downloaded the same time series used by Marcos et al. (2021, [29]). On the other hand, we did not use the time series from Marcos et al. (2011, [28]) in the case of Cádiz sea level. Although we included the correction of 37.5 mm for the period 1880–1924, proposed by these authors, we used the time series provided by the PSMSL. The difference between both time series and the possible influence of the method used to fill the gaps would account for the differences found in the trend estimations. For these reasons, trends for the Cádiz sea level should be taken with caution. Nevertheless, we are confident in the qualitative results that show that sea level trends increased during the 1990–2018 period and were lower from the early 1940s to 2018. The trend obtained for the complete period of time would correspond to an intermediate value between the former periods.

The differences between the linear trends estimated for different sub-periods simply reflect the decadal variability super-imposed on the monthly and inter-annual variability. According to Figure 3, Figure 4 and Figure 5, local sea level did not increase significantly or even decreased from the beginning of the 1940s or 1960s decades until the beginning of the 1990s. With sea level trends estimated from the beginning of the 1940s to 1990 (period not included in Table 2), the values obtained showed low positive trends of 0.47 mm/year in Alegeciras and Ceuta, 0.3 mm/year in Alicante, and negative values of −0.88, −0.36, and −1.22 mm/year at Cádiz, Tarifa, and Málaga, respectively. The sea level drop in the Mediterranean Sea from 1960 to 1990 has been extensively studied. Tsimplis and Baker (2000) analyzed time series of temperature and salinity in the deep waters of the Western and Eastern Mediterranean and time series of winter NAO index; they concluded that the sea level reduction could be linked to the temperature and salinity increase in these waters. Such changes could also be related to a strong positive phase of the NAO index. Tsimplis and Josey (2001, [15]) considered that most of the influence of the positive NAO index on sea level during the period 1960–1990 was through the inverse barometer effect, suggesting that atmospheric forcing should be considered in the analysis of sea level variability. Gomis et al. (2008, [24]) and Marcos and Tsimplis (2007, [14]) analyzed the effect of atmospheric forcing using a 2D barotropic model and also found that the inverse barometer effect was the responsible for the low positive or even negative trends observed in the Mediterranean sea level until the beginning of the 1990s decade.

In agreement with these findings, our results show that the atmospheric pressure increased in all the locations from 1948 to 1990 and that the halosteric component also decreased in all the regions considered. As already explained, previous studies have analyzed the effect of pressure and winds on sea level using barotropic 2D circulation models. In this study, we used a different and complementary approach. We analyzed the regression of local sea level on local atmospheric pressure and the U, V components of wind and on the thermosteric and halosteric contributions (Equations (2) and (3)). The forward stepwise regression selected which of the variables hada significant influence on the monthly and inter-annual variability of the observed sea level. The response of sea level to these variables (both atmospheric and steric variability) was expressed by the linear coefficients of the model. These coefficients enabled us to check the impact of long-term changes in atmospheric variables and the steric components on sea level. Jordà and Gomis (2013a; 2013b, [25,26]) estimated the contribution to sea level of mass addition using satellite gravimetry data and subtracting the steric contribution to the observed sea level, and they found a good agreement between both methodologies. However, the use of gravimetry data is not appropriate for the study of local sea level due to the low spatial resolution of such data. On the other hand, determining the exact contribution of the steric component of sea level is a difficult task, as it could be dependent on the reference level used for the integration of Equations (2) and (3) [23], and it is not clear how the steric component calculated from data collected at the open sea would contribute to the coastal sea level [4]. In the case of the halosteric component, a salinity increase would increase the density of sea water, and the volume of water per unit of mass would decrease. However, it is not clear if this salinity increase could be linked to changes in the water and salt exchange through the Strait of Gibraltar and could finally produce an increase in the sea level in the Mediterranean Sea [12]. For this reason, the steric components were calculated using the maximum depth at each region, but it was not assumed that the contributions obtained by means of Equations (2A) and (2B) were the real contributions to the observed coastal sea level. As in the case of the atmospheric variables, the thermosteric and halosteric components estimated from temperature and salinity fields were considered as potential predictors for the observed sea level, and the response of sea level was estimated from the linear regression.

The linear model was fitted for the two periods 1948–2018 and 1990–2018. Sea level reconstructed using the selected linear model explained a large percentage of the observed sea level variance. The agreement between observations and the statistical model is especially good for the period 1990–2018, when the multiple correlation coefficients reached values as high as 0.87 at L’Estartit. The lowest value was 0.61 for Alicante. The increase inthe variance explained by the linear model for the period 1990–2018 could be the result of the higher number of available TS data, which makes more accurate the estimation of the thermosteric and halosteric contributions. This would also explain the better performance of the model for the L’Estartit sea level. This location is close to the Gulf of Lion, which is an area intensively sampled because of its importance in the formation of deep waters.

In all the locations, the predictors selected to explain the observed variance of sea level were the same for both periods (see Table 4). The complete model (P, U, V, Δη_T, Δη_H) was selected for the sea level at Tarifa, Algeciras, and L’Estartit. In all the cases, the atmospheric pressure was one of the selected predictors, with coefficients ranging between −11.0 and −15.1 mm/hPa, close to the theoretical value of −10 mm/hPa. The U component of the wind was another factor contributing to the coastal sea level variability in all the locations, with Ceuta the only exception. The linear coefficients for this variable were negative and ranged between −5.9 mm/(ms⁻¹) and −12.4 mm/(ms⁻¹). The negative sign indicates the well-known effect of westerly winds that induce upwelling in the northern coasts of the Gulf of Cádiz, Strait of Gibraltar, and the Alboran Sea. This direction of the wind also produces the sea level decrease in the Alicante and L’Estartit stations, as it pushes water offshore. On the contrary, southerly winds (positive V) produce the accumulation of water against the coast at Tarifa, Algeciras, Alicante, and L’Estartit. As expected, the coefficients for the thermosteric and halosteric contributions were always positive when these variables were included in the model. Nevertheless, the response of the observed sea level to such contributions was smaller than 1 (Table 4), evidencing the already explained limitations of these calculations.

Relative sea level can be influenced by active tectonics [41,42] and by the redistribution of mass derived from the last deglaciation (GIA). In the present study, the local sea level trends were corrected for the GIA using the values provided in the PSMSL (Table 2; [38]). Then, the observed trends in atmospheric variables (Table 2) and the thermosteric and halosteric contributions (Table 3) were combined with the response of local sea level to these factors (Table 4). Multiplying such trends for the coefficients in the linear model, the influence of each factor on the observed long-term trends was estimated. If the trend estimated for a variable was not significant at the 95% confidence level, we could not be sure that it was different from zero, and its influence on the sea level trend was discarded. In case the trend for a given variable was statistically significant, whether or not it had been included in the linear model was considered. Finally, the influence of those variables that showed significant trends, and that could be considered as predictors for local sea level, were estimated by multiplying their linear trends by the coefficients in the linear model. As both the trends and the coefficients had errors, the errors of the contributions to the sea level trends were calculated using the formula for propagation of errors. These calculations were repeated for the periods 1948–2018 and 1990–2018.

Column 2 in

Table 5. Linear trends in mm/year and the contribution of each variable to the observed local sea level, corrected for the GIA. The second column shows the observed sea level trend corrected for GIA. Columns 3 to 7 show the contributions of P, U, and V and the thermosteric and halosteric components to the observed sea level. Column 8 is the sea level trend caused by mass addition. This column is equal to the second one minus columns 3 to 7. In all the cases, the first line corresponds to the period 1948–2018, and the second line corresponds to the period 1990–2018. When a variable does not contribute to the sea level trend because such variable experienced no trend or because it was not selected in the model, the table is filled with a line.

mm/year	Corrected (GIA)	P	U	V	ηT	ηH	Mass Addition
Cádiz	0.91 ± 0.17 3.4 ± 0.6	−0.27 ± 0.16 -------	0.08 ± 0.10 −0.25 ± 0.24	------- -------	------- -------	−0.10 ± 0.06 0.31 ± 0.24	1.2 ± 0.5 3.3 ± 1.1
Tarifa	1.00 ± 0.19 4.0 ± 0.5	−0.26 ± 0.14 -------	------- −0.21 ± 0.15	0.3 ± 0.3	0.27 ± 0.12 0.48 ± 0.22	−0.15 ± 0.06 0.18 ± 0.15	1.2 ± 0.5 3.3 ± 1.7
Algeciras	0.68 ± 0.13 2.3 ± 0.5	−0.25 ± 0.13 -------	------- −0.18 ± 0.13	0.24 ± 0.24	0.10 ± 0.08 0.24 ± 0.17	−0.08 ± 0.03 0.15 ± 0.13	0.9 ± 0.4 1.9 ± 1.5
Ceuta	0.47 ± 0.15 1.3 ± 0.5	−0.26 ± 0.13 -------	------- -------	0.16 ± 0.16	0.16 ± 0.07 0.20 ± 0.14	------- -------	0.6 ± 0.4 1.0 ± 1.2
Málaga	0.91 ± 0.18 4.5 ± 0.6	−0.29 ± 0.17 0.5 ± 0.5	------- −0.4 ± 0.3	------- -------	---------- ----------	−0.18 ± 0.09 -------------	1.4 ± 0.3 4.5 ± 0.9
Alicante	0.49 ± 0.17 2.2 ± 0.6	−0.27 ± 0.16 0.3 ± 0.4	------- -------	−0.06 ± 0.03 0.2 ± 0.3	------- -------	------- -------	0.8 ± 0.4 1.3 ± 1.4
L’Estartit	------- 2.4 ± 0.7	------- 0.5 ± 0.5	------- −0.07 ± 0.08	------- −0.1 ± 0.1	------- 0.5 ± 0.3	------- −1.1 ± 0.4	------- 2.7 ± 2.1

shows the trends for local sea level corrected for GIA. Columns 3 to 7 show the P, U, and V and the thermosteric and halosteric contributions to the observed sea level trends. Notice that a variable could have been selected as a predictor or explanatory variable for the monthly sea level variability, but if such a variable experienced no change in the long term, it would not contribute to the sea level trend. Conversely, a variable could be affected by a significant trend, but if the linear regression indicated that such a variable had no significant influence on sea level variability, it would not contribute to the sea level trend. When all the different contributions were subtracted from the observed sea level trend (corrected for GIA), the addition of mass could be estimated. For the first period (1948–2018), the mass component ranged between 1.2 ± 0.5 mm/year at Tarifa and 0.6 ± 0.4 mm/year at Ceuta. Considering the large uncertainties of these trends, these trends cannot be considered significantly different. Those locations with higher trends also had larger uncertainties. During the second period, 1990–2018, the contribution of the mass component to sea level trends increased in all the locations. This contribution was not significant at Ceuta (1.0 ± 1.2 mm/year) or in Alicante (1.3 ± 1.4 mm/year) and ranged between 1.9 ± 1.5 mm/year in Algeciras and 3.3 ± 1.7 mm/year in Tarifa. Once again, the values for all the tide gauges (with the exception of Ceuta) were not statistically different.

In summary, no a priori assumptions were made about the influence of atmospheric variables and steric components on observed sea level. The statistical analysis of de-seasoned and de-trended time series established the response of observed sea level to such potential contributions. This analysis shows the ability of a linear model to reproduce most of the variance in observed sea level time series. The performance of this model improves considerably for the period 1990–2018. There is a clear increase in sea level trends (corrected for GIA) in all the locations, being accompanied by an increase in the mass component. The contribution of the halosteric component was negative in the Mediterranean Sea during the two periods analyzed. In the Gulf of Cádiz and Strait of Gibraltar, it was negative during the first period (1948–2018) and positive during the second one (1990–2018). These negative values of the halosteric component simply indicate that there is a reduction in the volume of water per unit mass, but it is not clear yet if the salinity increase is linked to a mass addition associated with changes in the water and salt exchanges through the Strait of Gibraltar. This is a fundamental question that should be addressed in order to understand future changes in Mediterranean sea level.