Sea level variability at Adriatic coast and its relationship to atmospheric forcing

. Sea level (SLH) variability at the Adriatic coast was investigated for the period 1872–2001 using monthly average values of observations at 13 tide gauge stations. Linear trends and seasonal cycles were investigated ﬁrst and removed afterwards from the data. Empirical Orthogonal Functions (EOF) analysis was used further on remaining anomalies (SLA) to extract the regional intermonthly variability of SLH. It was found that the leading EOF and its principal component (PC) explain a major part of SLA variability (92%). The correlation between the reconstructed SLA, based on leading EOF and its PC, and overlapping observed SLA values for selected tide gauge stations is between 0.93 and 0.99. Actual SLH values at tide gauge stations can be reconstructed and some gaps in the data can be ﬁlled in on the basis of estimated SLA values if reasonable estimates of long-term trends and seasonal cycles are also available. A strong, seasonally dependent relationship between SLA at the Adriatic coast and atmospheric forcing, represented by sea level pressure (SLP) ﬁelds, was also found. Comparing the time series of leading PC and gridded SLP data for the period 1948–2001, the highest correlation coefﬁcients ( r ) of − 0.92 in winter, − 0.84 in spring, − 0.66 in summer, and − 0.91 in autumn were estimated for a SLP grid point located in northern Italy. The SLP variability on this grid point contains information about the isostatic response of the sea level at the Adriatic coast, but can also be treated as a sort of teleconnection index representing the large-scale SLP variability across central and southern Europe. To some extent the large-scale SLP variability that affects the SLA at the Adriatic coast can be related to the North Atlantic Oscillation (NAO), because signiﬁcant corre-lations were found between the NAO index and the ﬁrst PC of SLA ( r winter =− 0 . 56, r spring =− 0 . 45, r summer =− 0 . 48, and r autumn =− 0 . 43) for the period 1872–2001. The use of partial least-squares (PLS) regression between large-scale SLP ﬁelds and SLA only slightly improved the description of the SLA dependence on SLP forcing in comparison to the single grid point approach. A strong relationship between atmospheric pressure and the sea level could represent an additional possibility for ﬁlling in the gaps in the tide gauge data.


Introduction
The sea level is influenced by a large number of factors operating at different temporal and spatial scales. Real changes in the sea level are related to varying meteorological and oceanographic forcings (Woodworth et al., 1999), such as changes in ocean and atmospheric circulation, changes in temperature and the salinity of surface, intermediate and deep waters, changes in average air pressure and wind fields, changes in regional hydrology, etc. There are essentially two types of observations that measure the sea level directly and can be used to study the variability of the sea level in different time scales: tide gauge measurements, which in a few cases go back over 100 years, and satellite altimeter measurements, which are available only since early 1990s (Cazenave and Nerem, 2004). The main dissadvantages of tide gauge data, which were used in our study, are their poor spatial distribution and their attachment to the land. As tidal gauges measure local sea level height (SLH) relative to a benchmark on land, sea level changes measured by the coastal tide gauges contain contributions both from real changes in the sea level and from vertical crustal movements of the land upon which the gauges are situated (Woodworth et al., 1999). Zerbini et al. (1996) noted that the crustal movements are small compared to the decadal and multidecadal sea level variability, but can be of the same order as the long-term trends in the sea level. To extract the real changes in the sea level from the coastal tide gauge data, the variability due to the crustal movements needs to be averaged (Zerbini et al., 1996), which Table 1. Selected tide gauge stations together with its PSMSL label and geographical data (λ g -longitude and φ g -latitude), period with available data, and number of months with available data (N). is not a trivial problem. When focusing on intermonthly to interannual variability of sea level and its spatial patterns, as in our case, the tide gauge data can be detrended prior to further analysis, which removes the combined effect of the real long-term trends in the sea level and crustal movements.
The main goals of our study were a) to investigate the intermonthly variability of the sea level at the Adriatic coast on the basis of incomplete tide gauge data, b) to extract the patterns of coherent regional sea level variability, and c) to use these patterns for filling in the gaps in the data time series of single stations. Since some recent studies (Tsimplis and Baker, 2000;Tsimplis and Josey, 2001;Woolf et al., 2003) report a strong relationship between SLH in the Mediterranean and atmospheric circulation, especially North Atlantic Oscillations (NAO), we also investigated the relationship of SLH variability on intermonthly scale with atmospheric forcing represented by the variability of large-scale sea level pressure fields across Europe and by the NAO index.

Sea level heights and sea level pressure data
Monthly mean sea level heights (SLH) that were used in our study represent a part of the tide gauge data for the Adriatic Sea collected by thePermanent Service for Mean Sea Level (PSMSL) 1 (Woodworth and Player, 2003). Since only tide gauges spanning at least 20 years of data with missing values amounting to less than 25% were included in our study, data from 13 tide gauge stations covering different time periods (see Table 1) were used to investigate the sea level variability at the Adriatic coast. Data for other Adriatic tide gauge stations were used in our case only for filling in some short gaps in the data of selected stations using simple linear regression. These interpolated values were used as a first approximation for some missing values in the estimation of linear trends and seasonal cycles (Sects. 2.2 and 2.3), but not in the procedure for finding the coherent regional variability with Empirical Orthogonal Function analysis (Sect. 2.6) and the final interpolation of missing data. Stations with data available for less than 20 years or with data of questionable quality were not used in the study.
For the investigation of the relationship between the sea level variability and atmospheric forcing on an intermonthly scale, the NCEP/NCAR (Kalnay et al., 1996;Kistler et al., 2001) sea level pressure (SLP) reanalysis data 2 were used. The SLP reanalysis data with a horizontal resolution of 2.5 • ×2.5 • were extracted for the area with latitudes between 15 • N and 75 • N, and with longitudes between 80 • W and 50 • E. The period for which the SLP data were used, is 1948-2001. Since the NAO is one of the most important atmospheric variability modes in the Northern Hemisphere (Hurrell, 1995;Wanner et al., 2001;Marshall et al., 2001;Woolf et al., 2003) we also investigated the relationship between the NAO index and the SLH on an intermonthly scale. We used the monthly NAO index, defined by Hurrell (1995) as the normalized pressure difference between Reykjavik (Iceland) and Ponte Delgado (Azores) 3 for the period 1872-2001.
The data used in our study were prefiltered according to the methodology accepted within the ESEAS-RI project (Shaw et al., submitted, 2005)   afterwards before extracting the coherent regional patterns of sea level variability. Extracted patterns were used for interpolation of missing data and for investigation of the relationship between the sea level variability and atmospheric forcing represented by the variability of large-scale SLP fields and the NAO index. German, V., Kjaer, N., Pasarić, M., Plag, H.-P., Martinez Benjamin, J. J., Martinez-Garcia, M., Sušnik, M., Tel, E., Raicich, F., and Yildiz, H

Linear trend
Linear trends in the SLH data for selected stations were estimated on the basis of annual mean values. Since there was a shift in the sea level trend in Mediterranean in the early 1960s (Tsimplis and Baker, 2000), we separated the SLH data into two periods. After the examination of the SLH data we selected the year 1961 as a shift year for our study. Linear trends for selected tide gauge stations were estimated separately for the period with years ≤1961 and period with years >1961. The problem is that some of the stations have no (Koper and Bar) or only few data points (Porto Corsini, Dubrovnik, Rovinj) for one of the periods (see Fig. 1). For those stations we did not estimate trends for the short period and we excluded the data for this period from further analysis. When estimating linear trends, the annual mean values of the SLH data were fitted to a straight line as a function of time by a simple linear regression (Wilks, 1995;von Storch and Zwiers, 1999), separately for the two periods. The leastsquares fit to the annual mean of each year is defined as at least data for 10 out of the 12 months available, and the 2 missing months are not consecutive (Shaw et al., submitted, 2005 4 ). For further analysis the estimated trends were subtracted from the data. When investigating the relationship between the sea level variability and atmospheric forcing, the NAO index and SLP data were also detrended the same way before relating their variability to the SLH variability.

Seasonal cycle
The coastal recordings of sea level are likely to be biased by local effects, which indicates the potential benefits of the removal of seasonal cycles in the SLH time series before analyzing a coherent regional sea level variability. In many coastal areas, a large fraction of the seasonal cycle is steric, due to density changes as a consequence of temperature and salinity variations (Zerbini et al., 1996). Also, the non-steric component of the seasonal cycle, such as the influence of regional hydrology, was found to be important in some cases (Cazenave et al., 2002;Stanev and Peneva, 2002). The seasonal cycle and its amplitude of coastal SLH present the local response to a seasonal pattern of large-scale sea level variability and hydrology. We used only the SLH data for the full years to estimate the seasonal cycle, since the use of all available data, also the data for years with missing values, could bias the estimation of the seasonal cycle (e.g. Shaw et al., submitted, 2005 4 ). In comparison to the SLH data for the Adriatic Sea, the NCEP/NCAR reanalysis SLP and NAO index data have no missing values for the selected area and time period.
Before making the estimation of the coherent regional patterns of sea level variability and their relationship to atmospheric forcing on an intermonthly scale, the seasonal cycle was removed from the data by subtracting the long-term average monthly values for the full years from the original data.

Spectral analysis
Spectral analysis can be used for the investigation of potential periodicities in the time series (Ghil et al., 2002). The discreet Fourier transform (DFT) was used for the calculation of the spectrum (periodogram) of detrended and deseasonalised SLH data -sea level anomalies (SLA). The missing values in the data were replaced by zeros, which could introduce some artificial signal in the spectra, especially in the case of a large number of missing data. If y n represents the n-th observation of SLA among all N observations at a selected tide gauge station, the corresponding DFT of the time series (von Storch and Zwiers, 1999) can be written as The values of periodogram can be calculated as (Wilks, 1995;von Storch and Zwiers, 1999) where a k and b k are the Fourier coefficients for the selected SLH time series, Since exactly N Fourier coefficients are needed, k in Eq. (3) runs from 1 to K. In the case of an even N , K is N−2 2 , and in the case of an odd N , K is N −1 2 , but b K is 0. The periodogram plots on the figures are presented in logarithmic mode as 10 log 10 P k .

Cross-correlation and autocorrelation
The Pearson's correlation coefficient (Wilks, 1995;von Storch and Zwiers, 1999) was used as a measure of autocorrelation and cross-correlation of SLA time series: The termsū andv in Eq. (4) are the averagesū= 1 M M m=1 u m andv= 1 M M m=1 v m . In the case of the estimation of the cross-correlation between the SLA time series for two tide gauge stations, u m and v m in Eq. (4) represent the pairs of available overlapping SLA values for the same month for these two stations, and M represents the total number of available pairs. The autocorrelation function for a single SLA time series with lags from 0 to 12 months (l∈[0, 12]) was performed only on pairs of lagged SLA data without missing values. In this case M=N −l in Eq. (4) is the total number of such pairs, and (u m =y m , v m =y m+l ) is the pair of lagged data.

Empirical Orthogonal Functions (EOF) analysis for incomplete data
To remove the noise from climatological and oceanographic data sets, the Empirical Orthogonal Functions (EOF) analysis (Preisendorfer, 1988;von Storch and Zwiers, 1999) is commonly used. The classical EOF analysis assumes that the matrix Y, containing the time series of SLA (y) for different tide gauge stations, is perfect and completely known. In such a case, the EOFs can be calculated as the eigenvectors of the The columns of matrix E represent the eigenvectors of covariance matrix R -the EOFs. The EOFs are orthogonal and normalized. The elements of the diagonal matrix represent the eigenvalues of R containing the information about the variability of Y data explained by the single EOFs. The variability explained by the p-th EOF (the p-th column of E) can be estimated as V AR p =λ p / P p=1 λ p , where P is the number of tide-gauge stations, as the number of all EOFs, and λ p is the p-th eigenvalue. The time amplitudes of the EOFs, often called the principal components (PCs-columns of the matrix C), can be estimated as The data in matrix Y can now be written as Usually the first few EOFs and corresponding PCs contain most of the variability in Y, and the remaining are considered as noise in Y data. Due to the missing data in the available SLA time series, the classical approach to the EOF analysis was not adequate. For the calculation of EOFs and interpolation of missing data, the approach proposed by Beckers and Rixen (2003) and further applied by Alvera-Azcárate et al. (2005) was used in our case. Zeros were chosen for the first approximation of the missing values in the SLA time series. The EOF analysis was performed afterwards, and the new approximation of the missing values was the value reconstructed by the calculated EOFs and their time amplitudes (PCs). This procedure was repeated until the desired precision of data was reached (1 mm in our case, which is also the accuracy of the data). The question remains as to how many EOFs are needed for the reconstruction of measured Y data and the interpolation of missing values. In our case, the cross-validation approach (CV) (Wilks, 1995;von Storch and Zwiers, 1999) was performed for the estimation of the optimal number of EOFs. One-tenth of non-missing data was randomly excluded from Y, and treated as missing -replaced by zeros. The described EOF analysis for the missing data was performed on the adjusted Y, first by retaining only the first EOF. The number of retained EOFs was increased gradually afterwards, and the mean square error (MSE cv ) between the excluded and estimated values was calculated before every addition. The addition of new EOFs is meaningful while the MSE cv decreases. The optimal number of retained EOFs is where MSE cv reaches its minimum. To prevent the convergence in local minimum, this procedure was performed several times (10). The Y data was then reconstructed and the missing data interpolated with the optimal number of EOFs and corresponding PCs.

Partial Least-Squares (PLS) regression
Since the time series of SLP data on different grid points are commonly strongly correlated, a Partial Least-Squares regression (PLS) (Bro, 1998;Helland, 2001;Martens, 2001;Abdi, 2004) can be used to extract large-scale SLP features that are important for intermonthly SLH variability at the Adriatic coast. The general idea behind the PLS method is to extract a few important, non-correlated features (S=[s 1 . . . s Q ]) of a predictor (X), in such a manner that the covariance between this features (S) and the predicted (y) will be as high as possible. If every column of matrix X contains a time series of SLP for each grid point and y contains the time series of the coherent regional sea level variability at the Adriatic coast, both without linear trends and seasonal cycle, the idea of PLS can be written in a mathematical form as where X and y are centred. The main task is to estimate the weight matrix W responsible for the extraction of the important features, while the regression coefficient vector b can be estimated afterwards in a simple least-squares manner. For the estimation of W the algorithm proposed by Bro (1998) with the modification by de Jong (1998) was used. It can be summarized as follows: 1. Set the initial values of y q to the original values of y (q=0 and y 0 =y).
2. Calculate the covariance vector z q =X T y q .
3. Determine the weight vector 5 w q = z q z q that maximizes the covariance between new feature s q and y q .
4. Store w q as an additional column of weight matrix W=[w 1 . . . w q ].
5. Use the weight matrix W for feature extraction as S=XW.
6. Calculate the regression coefficients' vector as b=(S T S) −1 S T y.
7. Deflate y as y q =y−Sb to subtract the part explained by previously determined features.
8. Increase q to q+1 and continue from step 2 until a proper description of y is reached.
The procedure of feature extraction and estimation of regression coefficients was stopped, when adding an additional feature increased the mean square error of the cross-validation (CV) -MSE cv (Wilks, 1995;von Storch and Zwiers, 1999). The data were divided into 10 parts, and each of them was excluded once. The rest of the data were used for the development of the PLS model. The excluded y data were estimated with the developed model by using excluded X data. 5 In this step we have to maximize the covariance between new features s q and y q , with a condition of w q =1. Due to the properties of scalar product, the covariance s T q y q =w T q X T y q =w T q z q will be maximal, if the vector w q will have the same direction as the vector z q , meaning w q =ζ z q . Considering the condition of |w q =1 this leads to the solution of w q = The difference between the CV estimated and the original y was used to estimate the MSE cv . As in the case of CV for selecting the optimal number of retained EOFs of SLA, as well as in the case of PLS regression, the procedure was repeated several times (10) to avoid the convergence to the local minimum of MSE.
After selecting the optimal number of features for PLS model, the final model can be written aŝ 3 Results and discussion

Linear trends in the sea level at the Adriatic coast
The findings of Tsimplis and Baker (2000) concerning the shift in the sea level in the Adriatic sea around the year 1960 can also be confirmed by our results (Fig. 1). The shift in the early 1960s is especially evident in the annual averages of sea level at Trieste and Venezia (Punta della Salute). An increasing trend can be seen in the period before the year 1961 and becomes much weaker after this year. Some indications of the shift can also be found in data for locations Bakar and Split (rt Manjana and Harbour), despite the fact that data are available for only a few years before 1961. It can also be seen that after 1961 there are practically no significant trends (at 5% significance level: α=0.05) in the sea level at selected tide gauge stations with the only exception of station Bar in Montenegro (  (Fig. 2), so we included the station Porto Corsini in further analysis.
Estimated trends of the crustal movement (Table 2) are comparable among the selected tide-gauge stations and are slightly less than 1 mm year . If we subtract these values from the trends in SLH, we can see that there was an increasing trend in the sea level before the shift in the early 1960s and a slightly decreasing trend after the shift practically on the entire Adriatic coast.

Seasonal cycle in the sea level at the Adriatic coast
After the linear trends ( Fig. 1 and Table 2) were removed, seasonal cycle on selected tide gauge stations was examined. Long-term averages of mean monthly deviations (CY C) from linear trends, used as estimates for the seasonal cycle, are presented in Fig. 3. The amplitude of the annual cycle is in most cases about 120 mm, with maximum in late autumn (November). With regards to the minimum, stations can be separated into two groups. The first group has the tide-gauge stations from the North Adriatic coast and Istria (Venezia, Trieste, Koper, Rovinj), having one minimum at sea level in March. The second group has the stations from Dalmacia (Bakar, Split, Dubrovnik, Bar), having two almost equivalent minima, one in March and one in July. Seasonal cycle for station Porto Corsini is more similar to the seasonal cycle of the second group. We would expect similarity to the first group, but the deviation of the annual cycle for Porto Corsini from the annual cycle of stations in this group could be related to the problematic quality of data, especially in the second period with measurements. In general, the annual cycles of both groups are very similar, despite the mentioned differences in the second minimum in summer. Table 2. Estimates of the linear trend at the sea level (b est ) and their significance (sig) at the 5% level (α=0.05) for selected tide gauge stations at the Adriatic coast before and after the shift (1961), together with the number of years used in the trend estimation (N) and the trends for the land movement (b lnd ).

Coherent variability of sea level anomalies at the Adriatic coast
For the investigation of coherent variability of the sea level at the Adriatic coast, we used only the sea level anomalies (SLA) -detrended and deseasonalised SLH data. A strong cross correlation (r) can be found between SLA time series for selected tide gauge stations for overlapping periods (Table 3). With few exceptions, the values of the Pearson correlation coefficient are higher than 0.8, which indicates a good possibility for interpolation of missing data on the basis of other simultaneously observed SLH data across the Adriatic coast.   Using the Beckers and Rixen (2003) approach to EOF analysis for estimation of the missing values we found that 92% of the SLA variability on selected tide gauge stations can be explained by the first leading EOF (Fig. 4). Results of cross-validation have shown that additional EOFs do not contribute to the quality of reconstructed SLA data at selected stations but only contribute to the noise. On the basis of the first leading EOF (left upper plot in Fig. 4) and its PC (lower plot in Fig. 4) we can fill in some gaps in the SLA time series of selected tide gauge stations. Consequently, we can reconstruct the sea-level (SLH) data at selected stations if reasonable estimates of trends and seasonal cycle are available (see Sect. 3.4). The correlation between original SLA data for every single tide gauge station and PC of the first EOF are shown in the last column of Table 3. It can be seen that the correlation coefficient is higher than 0.93 in all cases, which means that the first EOF describes well the variability of SLA practically fort the entire Adriatic coast and that the presented approach can be useful for the interpolation of missing SLA values.
The temporal behavior of the PC for the leading EOF corresponds to a weak red noise, which can be seen from the values of autocorrelation coefficients (upper right plot in Fig. 4 -red, solid lines) and the periodogram (upper middle plot on Fig. 4 -blue, solid line). The spectrum of the corresponding red noise with autocoorelation coefficient r(1)=0.33 is presented in Fig. 4 as a solid black line. The periodogram has no obvious peaks, although in some cases, especially for long-term variability (lower frequencies), the values exceed the 5% significance limit for red noise (broken black line). Despite the fact that the theoretical spectrum of the weak AR(1) process (black, solid line) describes well the general behavior of the estimated spectrum of PC (blue line), it can be seen from the autocorrelation function that the order of AR is higher than 1. The autocorrelation function r(l) for l>1 months decreases slower than r l (1), which is the theoretical value of the autocorrelation function of the AR(1) process (upper right plot on Fig. 4 -black, dashed lines). Theoretical white noise spectrum (solid red line) and its 5% significance limit (broken red line) were also added to the upper middle plot of Fig. 4 for easier evaluation of the PC spectrum. Examples

Reconstruction of sea level data
Knowing reasonable estimates for trends (see Fig. 1 and Table 2) and seasonal cycles (see Fig. 3), we can reconstruct the actual sea level (SLH) at the tide gauge stations by adding the trend and seasonal cycle to the reconstructed SLA data (e.g. Figs. 5 and 6). In the case of the tide gauge stations in Istria (Koper and Rovinj) and Dalmacia (Bakar, Split, Dubrovnik, Bar) the reconstruction and filling in of the missing data is possible only for the period after 1961, as the trends were not estimated for the period before this year or were estimated on the basis of data from a short time period and are not useful for the extrapolation to earlier years. Due to the unusual trends in the data for station Porto Corsini, the reconstruction of SLH data for this station is also questionable. On the other hand, the reconstruction for the entire period 1872-2001 can be made for the stations Trieste and Venezia-Diga sud di Lido. An example of reconstructed SLH data, together with the observed SLH data for stations Trieste and Koper, is shown in Fig. 7.

Relationship between sea level anomalies and sea level pressure variability
The results of the EOF analysis show that the regional temporal variability of SLA at the Adriatic coast can be described well by a single PC time series corresponding to the leading EOF. A unique variability of SLA on the entire northern and eastern Adriatic coast suggests that some large-scale forcing, like variability in the sea level pressure and related atmospheric circulation, is probably responsible for intermonthly variability of SLH, in addition to the variability related to the seasonal cycle and long-term trends. It is well known that a substantial fraction of fluctuations in coastal SLH on time scales of months to decades is non-steric and forced mechanically by the atmosphere (Zerbini et al., 1996). Some recent studies (Tsimplis and Baker, 2000;Tsimplis and Josey, 2001;Woolf et al., 2003) report a relationship between SLH in the Mediterranean Sea, including the Adriatic Sea, and atmospheric circulation, especially NAO. Woolf et al. (2003) found that sensitivity to NAO is especially strong in the northern Adriatic. A statistically significant relationship between NAO and SLA in the Adriatic Sea was also confirmed in our study, as the coorelation coefficient (r) between the NAO index and the first PC of SLA is −0.48 when including all the months in the year and the entire period 1872-2001 (N =1548) (e.g. r critical =−0.32 for d.f.=N −2=100, α=0.001). As the NAO signal is the strongest in winter (Hurrell, 1995) the r values are even higher when using only data for the period December -February (r winter =−0.56, N=387). The r values for other seasons are slightly lower (r spring =−0.45, r summer =−0.48, and r autumn =−0.43, N =387). Tsimplis and Josey (2001) stated that the NAO affects the sea level in the Mediterranean via anomalies in freshwater flux (difference between precipitation and evaporation, river runoff) and the surface pressure, the latter occurring by the inverse barometer effect. For example, high NAO index does not indicate only high pressure over the Mediterranean, but also reduced winter precipitation and strong northerly winds (Hurrell, 1995).
The NAO index is only a crude measure of atmospheric variability and regional indicies may be more useful in particular regions (Woolf et al., 2003). Since the SLP is one of the most prominent forcing factors of SLH variability in monthly time scales (Zerbini et al., 1996), we also focused on the large-scale SLP patterns to locate the area of SLP variability that has the most important influence on the Adriatic SLA. Our goal was to estimate to what extent the SLA variability can be explained by the variability of large-scale SLP patterns. The SLP fields influence the SLA variability on a monthly time scale through the isostatic response, but also through the non-local wind effects related to the SLP fields (Ponte, 1994;Zerbini et al., 1996) and impacts on regional hydrology (Tsimplis and Josey, 2001;Stanev and Peneva, 2002). Due to the seasonal dependence of typical circulation patterns, the relationship between SLA and SLP fields in the period 1948-2001 was investigated separately for the months of all four seasons (winter -DJF; spring -MAM; summer -JJA; autumn -SON).
In a first step, the time series of leading PC, containing most of the information about the variability of SLA in the Adriatic region, was correlated to the gridded SLP values across the North Atlantic and European region. The results of the correlation analysis (Fig. 8) for the single grid point approach (SGP) show a similar spatial pattern in the correlation maps for all seasons, with the highest values across northern Italy. The difference is in the amplitude of the correlation patterns. The correlation is the strongest in the winter season, with r values up to −0.92, the weakest in the summer season, with r values up to −0.66, and with r values up to −0.91 for autumn and up to −0.84 for spring. This means that a time series of a single grid point located in northern Italy, representing the variability of SLP in a single cell of 2.5 • ×2.5 • , can explain a great part of the variability of the Adriatic SLA data -up to 84% in winter, up to 83% in autumn, up to 71% in spring, and up to 43% in summer. Such a time series contains information about the isostatic response of the sea level at the Adriatic coast, but can also be treated as a sort of teleconnection index representing large-scale SLP variability across central and southern Europe, because the SLP time series of neighboring grid points are correlated among each other. Variability of tropospheric mass distribution over the North Atlantic -Mediterranean sector assumes the form of a seasonally stationary wave pattern with opposite-sign nodes over Greenland-Iceland and northern Mediterranean (Eshel and Farrell, 2000). Thus, when Greenland pressure is anomalously high/low, anomalous low/high covers the northern Mediterranean, yielding southwesterlies/northeasterlies across the Adriatic that strongly influence the SLH. The correlation maps (Fig. 8) have a form that reflect some similarities with the NAO pattern (Hurrell, 1995), but the r values describing the relationship between the NAO index and the leading PC of the SLA for the same period (−0.55 for winter, −0.48 for spring, −0.54 for summer, and −0.49 for autumn) are much smaller than for a grid point across northern Italy. In order to describe better the relationship between largescale atmospheric forcing and the SLA at the Adriatic coast, we tried to find important SLP features that influence the SLA variability instead of only relating the SLA to a single grid point data series. Using PLS regression the variability of gridded SLP values across the region between 20 • N to 70 • N latitude and 40 • W to 50 • E was related to the leading PC of the SLA across the Adriatic region for the period 1948-2001. The predictor region was chosen as it led to the highest determination coefficients (R 2 ) of developed PLS models. The maps of regression coefficients (b PLS ) for all four seasons are shown in Fig. 9. The developed models for relating the first PC to the large-scale SLP field explain 89% of the variability of the reconstructed SLA in winter, 81% in spring, 55% in summer, and 83% in autumn, in the period 1948-2001. In all cases, a more complicated approach with PLS regression and large-scale SLP fields improves the SLP description of the SLH variability at the Adriatic coast only for a few percentage points in terms of R 2 (for winter 5%, for spring 10%, for summer 12%, for autumn 0%), in comparison to the use of SLP information for a single grid point. The resulting maps of regression coefficients are very complex. The highest importance (the most negative values of b PLS ) is given to the area of northern Italy and the northern Adriatic, similar as in the case of the correlation to a single grid-point SLP variability, with the exception of the summer season. The negative values of b PLS across the southern and central Europe are surrounded with weaker centres of positive b PLS values. According to the authors' opinion, the price of complexity of the PLS method with high calculation demands, and complex b PLS maps that can not be explained easily, is to high to recommend the PLS approach in further analysis in comparison to a simple SGP approach. However, from both, the SGP and PLS approaches we can see that in all seasons, with the exception of summer, a great part of the variability of the monthly average SLH can be related to the variability of the large-scale SLP field. Lower values of explained variability in the summer months could be related to more intensive regional weather processes, resulting in winds that can not be described by large-scale SLP patterns, but influence the regional SLH variability.

Conclusions
Our study confirmed the findings of Tsimplis and Baker (2000)   The average amplitude of the seasonal cycle at the selected tide gauge stations at the Adriatic coast cycle is about 120 mm, with maximum in November. The occurrence of minimum in the seasonal cycle separates selected tide gauge stations into two groups, one including stations located at the coast of the Northern Adriatic and Istria, having minimum at the sea level in March, and the other including stations at the coast of Dalmacia, having two almost equivalent minima, one in March and one in July.
Adjusted EOF analysis proposed by Beckers and Rixen (2003) was found to be successful in the reconstruction of monthly sea level anomalies, treated as deviations of the sea level from the annual cycle and the long-term trend. The method was also found to be useful as a part of the procedure for filling in the gaps in the sea level data from the tide gauge stations at the Adriatic coast if reasonable estimates of the seasonal cycle and long-term trend for the selected station are available. Applying adjusted EOF analysis, a coherent regional sea level variability can be separated from a local noise even in the case of incomplete sea level data covering different time spans. The regional variability of sea level anomalies at the Adriatic coast was found to be strongly related to the atmospheric forcing, represented by sea level pressure, in all seasons with the exception of summer. The forcing that affects the sea level variability at the Adriatic coast can be described fairly well with the sea level pressure variability across northern Italy, represented by the variability of the sea level pressure at a single grid point of the large-scale field with a resolution of 2.5 • ×2.5 • . A more complicated approach with the PLS method did not capture the relationship between the air pressure forcing and the sea level variability much better than a single grid point regression approach. Empirical models based on the strong relationship between the atmospheric pressure and the sea level could represent an additional tool for the adjusted EOF analysis in filling in the gaps in the tide gauge data.
A strong relationship between the large-scale air pressure variability and the sea level indicates that potential changes in atmospheric circulation in the future could introduce some long-term changes in the sea level related to the climate change. But on longer time scales, the expected increase in the global temperature (Houghton et al., 2001) will still play a major role in long-term sea level changes. On the other hand, the atmospheric circulation will remain one of the major forcings of the short-term variability of sea level, responsible for the extremes.