Coastal Sea Level Trends from a Joint Use of Satellite Radar Altimetry, GPS and Tide Gauges: Case Study of the Northern Adriatic Sea

Abstract

For the last century, tide gauges have been used to measure sea level change along the world’s coastline. However, tide gauges are heterogeneously distributed and sparse in coverage. The measured sea level changes are also affected by solid-Earth geophysics. Since 1992, satellite radar altimetry technique made possible to measure heights at sea independent of land changes. Recently various efforts started to improve the sea level record reprocessing past altimetry missions to create an almost 30 year-long combined record for sea level research studies. Moreover, coastal altimetry, i.e. the extension of altimetry into the oceanic coastal zone and its exploitation for looking at climate-scale variations of sea level, has had a steady progress in recent years and has become a recognized mission target for present and future satellite altimeters. Global sea level rise is today well acknowledged. On the opposite, the regional and local patterns are much more complicated to observe and explain. Sea level falls in some places and rises in others, as a consequence of natural cycles and anthropogenic causes. As relative sea level height continues to increase, many coastal cities can have the local elevation closer to the flooding line. It is evident that at land-sea interface a single technique is not enough to de-couple land and sea level changes. Satellite radar altimetry and tide gauges would coincide at coast if land had no vertical motion. By noting this fact, the difference of the two independent measurements is a proxy of land motion. In this chapter, we review recent advances in open ocean and coastal altimetry to measure sea level changes close to the coasts over the satellite radar altimetry era. The various methods to measure sea level trends are discussed, with focus on a more robust inverse method that has been tested in the Northern Adriatic Sea, where Global Positioning System (GPS) data are available to conduct a realistic assessment of uncertainties. The results show that the classical approach of estimating Vertical Land Motion (VLM) provides values that are almost half of those provided by the new Linear Inverse Problem With Constraints (LIPWC) method, in a new formulation which makes use of a change of variable (LIPWCCOV). Moreover, the accuracy of the new VLM estimates is lower when compared to the VLM estimated from GPS measurements. The experimental Sea Level Climate Change Initiative (SLCCI) data set (high resolution along track) coastal sea level product (developed within Climate Change Initiative (CCI project) that has been also assessed in the Gulf of Trieste show that the trends calculated with the gridded and along track datasets exhibit some differences, probably due to the different methodologies used in the generation of the products.

1. Introduction

Sea level rise is primarily an issue at the boundary line with land. It represents a potential threat to infrastructures and population living in low-elevation coastal areas [1]. The land disappears not only because the rising sea changes the coastline, but also because at a place there could be the land moving up or down, therefore contrasting or accelerating sea level rise [2]. Sea level can change significantly from one coastal location to another, as a result of a number of ocean, atmospheric and land processes that occur at various spatial and temporal scales [3].

In a global change scenario, as speculated in Li et al. [4], a slow rise of the sea level of few cm associated to climate change would make a difference to the coastline. It would not retreat from land, making it permanent. The flooding line of transient events (e.g., storm surges, tsunamis, etc.) would also uplift, increasing the risk of more frequent land inundation and more inland propagation [5]. An example is the City of Venice that has long been vulnerable to short duration flooding during winter [6]. The problem was so important that a system of 78 storm gates, known as MOSE [7] has been constructed to protect the city when high water is expected in Venice [8]. Long-term rising sea levels will represent additional challenges in the future [8].

Understanding the climate-related contribution to the sea level change and how much it will likely affect coastal regions is a major challenge, as it also requires local-scale measurements of the land effects. In this chapter, we review the sea level trend measuring system involving the integration of recent satellite-based observations from radar altimeters and Global Positioning System (GPS) receivers with historical data from tide gauge stations. The latest advances in open ocean and coastal altimetry to measure sea level changes close to the coasts over the satellite radar altimetry era are also summarized. A more robust inverse method to estimate sea level trends is also presented. It has been tested in the Northern Adriatic Sea, where GPS data are available to conduct a realistic assessment of uncertainties.

2. Techniques for measuring sea levels

Since Roman period, sea level has been measured nearby land just sticking a graduated pole within protected piscinae [9]. Since the 19th Century, tide gauges have been used in some coastal places around the world to measure the local change in sea level relative to the adjacent land [10, 11]. The baseline for measuring sea level over time is typically a mean sea level computed by averaging all the measurements over a period of years at each location. This relative sea level that will rise if ocean levels rise and/or land levels fall is the net change in the sea level and is the quantity of interest to the local coastal community in the real-time monitoring.

However, understanding the future coastal sea level changes and their relative significance requires to remove the effect of waves, tides, and other short-term fluctuations. But tide gauges alone cannot determine whether the sea level is rising, the land is sinking, or both. Sea level can rise or retreat in the long-term in response to the natural processes that alter the volume of water, including the climate-related contribution. The land level changing over time (the so-called vertical land motion, VLM or subsidence/uplift) can rise or fall due to natural processes (e.g., tectonic shifts; sediment loading, glacial isostatic adjustment, etc.) but also as a consequence of man-made factors (e.g., ground water extraction; oil and gas pumping, etc.).

There are various techniques for measuring VLM, e.g., geotechnical investigations using spirit levels and borehole extensometers [12]; geodetic surveying with Global Positioning System (GPS) satellite technology [13]; satellite remote sensing observations that use a technique called interferometric synthetic aperture radar (InSAR) [14]. The advantage of satellites is that they ensure almost global coverage in a repeatable manner and consistency of the measuring system over long time periods that is an important requirement for the detection of slow changes over time.

The quantification of VLM before the modern satellite era is difficult due to the poor coverage of geotechnical investigations. The advent of GPS receivers and their co-location with tide gauges made possible to continuously measure land elevation changes with simultaneous sea level reading at the same location [15]. GPS sensors return vertical and horizontal positions. The vertical position is a measure of the elevation of the land surface relative to the center of Earth, also referred as absolute. It is generally two to three times less precise than the horizontal components. The present picture is that, while there are many tide gauges around the world, not all have permanent GPS receivers co-located or near them [16].

The InSAR tool uses repeating multiple satellite radar imagery to create a time profile of land elevation change. The advantage respect to the GPS technology is the much higher density of VLM measuring points in the imaged area. The technique to provide statistically significant results requires a sufficient number of images and reduced scattering over time for the area of interest. The availability of images from the first satellites (e.g., ERS and Envisat) can be very irregular both in time and space [17]. Only the recent Sentinel-1 constellation provides global coverage and more frequent revisiting. Other satellite missions (e.g., COSMO-SkyMed, RADARSAT, etc.) only provide imagery on demand.

Sea level can be also measured with satellites using radar altimeters. These instruments send microwave pulses down along the satellite’s ground track and measure their echoes, revisiting the same place every 10 days or more depending on the mission. The time their echoes take to bounce back allows the system to measure the satellite’s altitude above the sea surface (the so-called range). It can be then corrected for instrumental and environmental effects. Knowing the satellite orbit with respect to Earth’s center of mass, the absolute, not relative, sea level can be thus calculated, and its change tracked over time.

Routine sea level observations began in 1992 with the TOPEX/Poseidon spacecraft on a 10-day repeat cycle, and this has subsequently been followed up by the Jason 1/2/3 series and the recent Sentinel-6 mission, providing a near-global fully consistent along track data set of sea level to understand how sea levels have changed over the past nearly three decades. Over the years, various satellite missions with different orbital configurations and other scientific objectives were launched, e.g., ERS-1/2, Envisat, Sentinel-3, SARAL, CryoSat-2 and HY-2A/B.

But single satellites have limitations. The sea level is tracked along paths whose distance is relatively large. A satellite alone could not fly in the region of interest, as it is for example the case of Venice for the T/P-Jason-Sentinel-6 family. Moreover, it has been difficult to retrieve data near coast where both the presence of land and more complex ocean surface scattering make the standard processing problematic [18].

3. Coastal sea level from satellite radar altimetry – a review

The satellite radar altimetry system has been initially conceived for usage in open ocean. The processing of radar echoes and the development of corrections is now at mature stage in this domain, with the various datasets routinely used for global sea level studies. However, data were normally flagged as bad and therefore rejected in the coastal zone. But the situation rapidly changed in the last ten years for two reasons: (1) the prospect of recovering a valuable long-term sea level data around the global coastline; (2) the improved suitability of the new and future altimeters (like those on CryoSat-2, AltiKa, Sentinel-3, Sentinel-6, Crystal). Therefore, a new domain “coastal altimetry”, i.e. the extension of altimetry into the oceanic coastal zone has been emerging, with a community around it developing a set of coastal altimetry techniques in order to get more and better sea level data closer to the coast.

The analyses of radar echoes revealed that pulse-limited missions, if reprocessed with dedicated models, could provide reliable range measurements to few km from the coastline. An example is the Adaptive Leading Edge Subwaveform (ALES) retracking algorithm, that has been validated and applied successfully to sea level research, demonstrating the ability to increase the quality and the quantity of sea level retrievals in coastal areas [19].

In addition, it was noted that geophysical corrections that must be applied to altimeter range data have a significant impact in coastal altimetry and therefore their constant improvement is crucial. There have been noticeable developments to improve the tropospheric delay [20], the tidal sea level where global models have still large errors [21] and the mean sea surface models, suitable for the observation of the coastal sea level [22]. There have been also improvements in procedures to avoid aliasing of major tidal signals and short-period ocean response to meteorological forcing aliases onto low frequency signals [23].

The wet tropospheric correction is the major source of uncertainty in altimetry budget error, due to its large spatial and temporal variability: this is reason why a multi-channel passive microwave radiometer is on the same platform as the altimeter. Unfortunately, this estimate gets quickly corrupted as soon as land enters the radiometer footprint, i.e. 20–50 Km from the coast. Alternative corrections have been devised and appear to be successful at least in some particular conditions [24]. A very promising approach was the one attempting to estimate the wet tropospheric path delay from GPS measurements known as GPD (GNSS-derived Path Delay), and its latest version called GPD+ (Plus) [25].

The classical data editing used in open ocean was also considered excessively restrictive and revisited with novel editing/re-interpolation approaches (e.g., [26]). The new data from the various reprocessing efforts are now bringing altimetry around the global coastline, with a higher spatial resolution and precision that was previously not available in coastal and shelf sea areas, while constant improvement [18] and validation [27] are still ongoing. The new coastal altimetry datasets open a new opportunity to study sea level change from open ocean to coast and differences in trend and variability at various distances from the coast, also nearby tide gauges [28].

4. A new sea level record from satellite radar altimetry for climate studies

Several radar altimetry missions have been in operation since the first launch in 1973 (see Figure 1). The TOPEX/Poseidon and Jason series (with the addition of the just launched Sentinel-6) is the reference mission for long-term sea level studies, as it is ensured the continuity in the same orbit [29]. However, a single altimeter only provides measurement along a track from open ocean closer to coast. There is always a trade-off between temporal sampling and ground-track spatial coverage. A single altimeter always leaves gaps along the coast between neighboring tracks: tenths to hundred km are not covered, so that the vast majority of the worldwide coast is not sampled. The coverage can be augmented with additional existing altimeters.

Figure 1.

Main characteristics of satellite altimetry missions operating until now and planned for the future.

Data from the various altimeter missions were used to create several datasets. Examples include RADS [30], X-TRACK [28], etc. that also provide sea level estimates. Since 1992, at least two altimeter satellites have been operating simultaneously, and during some periods, even more than two. Such data can be combined in a single product to provide a consistent long-term sea level data set, globally with sufficient spatial coverage over almost three decades. However, altimeter missions need to be accurately homogenized and cross-calibrated to reduce biases and uncertainties [31].

A satellite-based sea level data set to analyze long-term trends that uses the available historic observations from the various radar altimeters is key requirement for the climate community [32]. A recent reprocessing within the European Space Agency (ESA) Sea Level Climate Change Initiative (SLCCI) has produced a gridded altimetry product with a spatial resolution of 0.25° (which is around 25 km resolution) from 1993 to 2015 [33, 34], thus permitting a more detailed view of sea-level change around the world coastlines.

The sea level Environment Climate Variable (ECV) (at global and regional scales) is now operationally produced by the Copernicus Climate Change Service (C3S) [35] by applying the altimetry processing standards developed in the SLCCI initiative. The C3S product ensures a stable number of two altimeters since the beginning and the reference field used to compute sea level anomalies (SLA) is a homogeneous mean sea surface for all missions. The C3S record is a regional product, gridded at 0.125° in the Mediterranean Sea, starting in 1993 and offering ongoing coverage [36]. Both the SLCCI and the C3S datasets are state-of-the-art products designed to be a reference for climate-related sea level studies.

In the case-study illustrated in the chapter, the SLCCI and C3S datasets are used to assess their maturity as state-of-the-art altimetry datasets in climatological studies. The multi-mission gridded products have not still tuned for last 10 km from the coast, where the amount of valid data might decrease. The ESA CCI + Sea Level project, started in 2017, is extending the processing to the coastal zone, and an experimental coastal sea level product is going to be released to the public, in six selected regions: Northern Europe, Mediterranean Sea, Western Africa, North Indian Ocean, Southeast Asia and Australia [37]. This product is along-track and combines the enhanced spatial resolution provided by high-rate data (20-Hz), the post-processing strategy of X-TRACK and the advantage of the ALES retracker [38]. The product relies on the GPD+ wet tropospheric correction [39] and the FES2014 tidal corrections [40]. The X-TRACK/ALES SLCCI 20 Hz along-track dataset will be indicated with SLCCI-AT hereinafter.

5. Methods of estimating sea level trends

The trend is an indicator describing how sea level has changed over long time. It provides a simple predictive scenario if what observed in the past might be representative in the near future. The classical approach is to calculate a straight line through sea level data using a linear regression. The most used method for fitting data is least squares. However, other methods based on more complex models exist to estimate trends from sea level time series [41]. The trend estimation is sensitive to the length of the record and start/end periods. There might be variability at different inter-annual to decadal timescales occurring within the data. Moreover, in addition to the linear trend, there might be autocorrelation of the noise in the data [42].

A single tide gauge cannot explain to what extent the observed trend is related to ocean and/or land changes, without any nearby GPS. With the advent of satellite radar altimetry and the possibility to use altimeter passages nearby tide gauges a new method was proposed by Cazenave et al. [43]. It assumes that both the tide gauge and altimetry system measure the same ocean signal and the difference is a measure of VLM at the gauge: hereinafter we refer to this method as the “direct” or “classical” method. Another assumption is that there are no instrumental errors introducing significant drifts. This direct method provides VLM at the selected tide gauge station only.

Different implementations of the basic idea were successively proposed involving more tide gauges, more rigorous error analysis with mitigation of the uncertainties introduced by the assumptions and taking advantage of longer and improved altimeter-derived time series available at that time (e.g., [44, 45, 46, 47] and others).

An advanced method to estimate VLM that includes supplementary constraints from adjacent tide gauges has been proposed by Kuo et al. [48]. Its solution is based on the inversion of a linear system, formed mixing differences of altimetry- and tide gauge-derived trends, and differences of trends from neighboring tide gauges only, introduced in the linear system through Lagrange multipliers. As the solution of such a system requires its inversion, the method is referred to as Linear Inverse Problem with Constraints (LIPWC), or shortly “inverse” method. The new method optimally combines short-term altimetry records with long-term tide gauge observations. It assumes that absolute sea level change at tide gauges over a long time span is the same. The advantage of the method is that long (>40 years) tide gauge records contribute to reduce the error in the final VLM solution, and random and systematic errors in one or more time series trend are shared among all the other, cutting down the impact on the originating one. The disadvantage is that the method cannot be applied if the absolute sea level change is different from place to place. Nevertheless, this method can be useful in closed and semi-enclosed basins and could be adapted to work also in case a GPS at the coast is used instead of a tide gauge.

Kuo et al. [48] applied the inverse method within a semi-enclosed sea (Baltic Sea region of Fennoscandia). The results showed a significant reduction of uncertainties compared with those from conventional approaches, which are limited to the overlapping periods between altimetry and tide gauges. An extension of the method has been applied to Great Lakes and in open ocean regions, such as Alaskan coast [49]. It has been also extended along the coasts of southern Europe [50] with constraints between pairs of tide gauges based on correlation and overlapping periods. The same method has been extended to open ocean in New Zealand straddles, the Tasman Sea and Pacific Ocean [51]. All studies confirmed the superiority of the inverse method to the classical direct approach.

A new variant of the inverse method considers to difference sea level trends between pairs of tide gauge records and pairs of altimetry records [52]. Another study proposed different mathematical and statistical models, which enable simultaneous estimation of absolute and relative sea level trends and VLM at a tide gauge station merging altimetry and tide gauge records without the aid of geological information or GPS measurements [53].

6. A revisited linear inverse model to estimate sea level trends

The linear inverse model proposed by Kuo et al. [48, 49], and then by Wöppelmann and Marcos [50], assumes that the absolute sea level change rates are similar at all the tide gauge (TG) sites. This assumption is particularly important for the successful inversion. The explanation will be provided in this section.

The difference between the absolute sea level rise (ASLR) and the relative sea level rise (RSLR) rates, i.e. the velocities at which the sea level vertical motion is observed by satellite altimeters and TGs, denoted respectively with $\dot{g}$ and $\dot{s}$, is an estimate of the vertical velocity at which the land beneath TGs is moving. Such vertical crustal velocity, as previously stated, is named vertical land motion (VLM) and indicated by $\dot{u}$. A subscript $i$ is added to denote that the quantities $\dot{g}$, $\dot{s}$ and $\dot{u}$ refer to the $i$-th TG of a group of N:

$${\dot{u}}_i={\dot{g}}_i-{\dot{s}}_{i}i=1,\dots ,N$$E1

Eq. (1) is sufficient to obtain good estimates of the VLM rates at each TG, provided that all the variables in the equation refer to the same period and to coherent geophysical processes and have negligible inherent drifts and errors. Eq. (1) can be expressed in vector–matrix notation:

$$G\dot{u}=d;G=I_N$$E2

where $\dot{u}$ is the column vector of the unknown VLMs $\dot{u}=\left({\dot{u}}_1,\dots ,{\dot{u}}_N\right)$, and $d$ is the column vector whose elements are formed by the right-hand side of Eq. (1): $d=\left({\dot{g}}_1-{\dot{s}}_1,\dots ,{\dot{g}}_N-{\dot{s}}_N\right)$. In this picture all the unknown VLMs are mutually independent, and the linear system is easily inverted, offering the solution component by component. However, the solution is affected by large errors, as the period during which Eq. (1) is valid corresponds to the overlap period of TG and satellite altimetry observations, and thus no more back in time than 1992. In fact, the current time span of satellite altimetry data is less than 30 years. Such a short time span hinders the derivation of accurate trends from altimeter-gauge time series, as they are affected by inter-annual and decadal sea level signals, in particular by the 18.6-year lunar nodal tide, leading to uncertainties of the order of 1–2 mm yr⁻¹ [47, 54, 55]. For this reason, Kuo et al. [48] proposed a more elaborate linear system, in which constraints formed by the differenced time series of TGs over longer time periods (>40 yr) pose strong limits to the magnitude of the final errors thanks to the length of the time series. Such constraints are formed imposing that the rate of relative vertical motion between two TGs must equal the difference of their VLMs:

$${\dot{\mathit{ru}}}_{\mathit{ij}}=\left({\dot{g}}_i-{\dot{s}}_i\right)-\left({\dot{g}}_j-{\dot{s}}_j\right).$$E3

At this stage, the constraints still contain explicitly the ASLR at sites $i,j$, which are not known back in time beyond the beginning of the altimetry era. But if each couple of TGs in Eq. (3) are observing the same ASLR for some reason, for example they can be inside a lake or a semi-enclosed basin, Eq. (3) simplifies to:

$${\dot{\mathit{ru}}}_{\mathit{ij}}=\left({\dot{s}}_j-{\dot{s}}_i\right),$$E4

leaving out any reference to the ASLRs. Containing only the differences of the RSLR, Eq. (4) can be extended to the whole period of overlapping observations of the two TGs, which usually are longer and affected by lower errors. To distinguish the RSLR observed at the TG in the altimetry era from that observed in the common, longer period of observations of TGs $i,j$, we rewrite Eq. (4) as:

$${\dot{\mathit{ru}}}_{\mathit{ij}}={\dot{\mathrm{\zeta }}}_j-{\dot{\mathrm{\zeta }}}_i$$E5

where we used the Greek letter $\mathrm{\zeta }$ to indicate that the TG RSLR difference is calculated over the complete overlapping time span of the two tide gauges, even before the altimetry era.

The two linear systems for the Eqs. (2) and (5) are written in vector–matrix form as:

$$G\dot{u}=d;\dot{u}=\left(\begin{array}{c}\dot{u_1}\\ ⋮\\ \dot{u_N}\end{array}\right);d=\left(\begin{array}{c}\dot{g_1}-\dot{s_1}\\ ⋮\\ \dot{g_N}-\dot{s_N}\end{array}\right);G=I_N,$$E6

$$F\dot{u}=h;h=-F\left(\begin{array}{c}\dot{{\mathrm{\zeta }}_1}\\ ⋮\\ \dot{{\mathrm{\zeta }}_N}\end{array}\right);F=\left(\begin{array}{cccccc}1& -1& 0& 0& \dots & 0\\ 0& 1& -1& 0& \dots & 0\\ & \dots & & & \dots & \\ 0& 0& 0& \dots & 1& -1\end{array}\right).$$E7

The matrix F is the design matrix by which the constraints are formed and introduced in the linear system. The constraints can be chosen arbitrarily, but they must be linearly independent so that the rank of matrix $F$ is L ≤ N – 1. For the system to admit an OLS solution, L, the rank of the matrix $F$, must be < N, so that one degree of freedom is left in the linear system for the OLS procedure to perform the unknowns estimate. Without such degree of freedom, the system would become even-determined, and the N constraints ${\dot{\mathit{ru}}}_{\mathit{ij}}={\dot{\mathrm{\zeta }}}_j-{\dot{\mathrm{\zeta }}}_i$ would automatically determine the solutions for the unknown ${\dot{u}}_i$, and there would be no need for an OLS estimation.

Not always assumption ${\dot{g}}_i-{\dot{g}}_j=0$ is valid for every couple of TGs paired in a constraint. In such case independent linear systems are to be considered for each group of homogeneous TGs. Wöppelmann and Marcos [50] for example, applying this method to the seas of the southern Europe, treated separately the sites inside the Mediterranean Sea from those on the Atlantic coast of the Iberian Peninsula, as the oceanographic behaviors of the two sets of TG sites were observed to be markedly different. In this case the two linear system are totally independent and there cannot be a connection between the two groups.

A unique linear system incorporating the constraints (7) in the system (6) is formed recurring to the Lagrange multipliers technique: the inverse linear problem with constraints [56] (LIPWC). It stems from the minimization of the expression $\mathrm{\Phi }\left(\dot{u}\right)=e^{T}e+2{\lambda }^T\left(F\dot{u}-h\right)$, obtained as the sum of the L₂ norm of the prediction error $e=G\dot{u}-d$ and the inner product of the constraint equations $F\dot{u}-h=0$ by the Lagrange multipliers ${\lambda }^T$. The resulting ordinary least squares (OLS) equation in matrix form is:

$$\left(\begin{array}{cc}I& F^T\\ F& 0\end{array}\right)\left(\begin{array}{c}\dot{u}\\ \mathrm{\lambda }\end{array}\right)=\left(\begin{array}{c}d\\ h\end{array}\right)$$E8

which is a linear system of the form $A·X=Y$, with $A=\left(\begin{array}{cc}I& F^T\\ F& 0\end{array}\right)$. Such system is solved by direct inversion, provided the inverse of matrix A exists:

$$X=\left(\begin{array}{c}\dot{u}\\ \mathrm{\lambda }\end{array}\right)=A^{-1}Y;Y=\left(\begin{array}{c}d\\ h\end{array}\right)$$E9

The standard errors of the $\dot{u}$ are estimated as the diagonal elements of the covariance matrix $A^{-1}\mathrm{\Omega }{\left(A^{-1}\right)}^T$, with Ω given by

$$\mathrm{\Omega }=\left(\begin{array}{cccc}{\mathrm{\sigma }}_1^2& 0& \dots & 0\\ 0& {\mathrm{\sigma }}_2^2& \dots & 0\\ ⋮& ⋮& \ddots & ⋮\\ 0& 0& \dots & {\mathrm{\sigma }}_{N+L}^2\end{array}\right)$$E10

where N and L are respectively the number of parameters ${\dot{u}}_{\mathit{ij}}$ and constraints ${\dot{\mathit{ru}}}_{\mathit{ij}}$. The previous expression for Ω holds assuming no autocorrelation and heteroscedasticity of the regression residuals. The resultant errors of the OLS estimators are generally referred to as heteroscedasticity-consistent standard errors or White–Huber robust standard errors [57]. Finally, the estimated parameters are given by $X_i\pm \delta X_i$:

$$X=A^{-1}Y;X_i=\left\{\begin{array}{c}{\dot{u}}_i\text{if}0\le i\le N\\ {\mathrm{\lambda }}_i\text{if}N+1\le i\le N+L\end{array}\right.$$

$$\delta X_i=\mathit{diag}{\left(A^{-1}\mathrm{\Omega }{\left(A^{-1}\right)}^T\right)}_i$$E11

A possible attenuation of the condition that ${\dot{g}}_i-{\dot{g}}_j=0,$ for each pairs of TGs involved in a constrain, arises from the observation that if ${\dot{\mathit{ru}}}_{\mathit{ij}}=\left({\dot{g}}_i-{\dot{s}}_i\right)-\left({\dot{g}}_j-{\dot{s}}_j\right)$ can be rewritten as ${\dot{\mathit{ru}}}_{\mathit{ij}}=\left(-{\dot{s}}_i^{\prime }\right)-\left(-{\dot{s}}_j^{\prime }\right)$, then the reference to the ASLR has disappeared. Such situation can be achieved by a change of variable, as proposed by De Biasio et al. [58]:

$${\dot{g}}_i\to {\dot{g}}_i^{\prime }={\dot{g}}_i-{\dot{g}}_i=0$$

$${\dot{s}}_i\to {\dot{s}}_i^{\prime }={\dot{s}}_i-{\dot{g}}_i$$

$${\dot{\zeta }}_i\to {\dot{\zeta }}_i^{\prime }={\dot{\zeta }}_i–{\dot{g}}_i$$

$${\dot{u}}_i={\dot{g}}_i-{\dot{s}}_i={\dot{g}}_i-{\dot{g}}_i+{\dot{g}}_i-{\dot{s}}_i={\dot{g}}_i^{\prime }-{\dot{s}}_i^{\prime }.$$E12

In other words, we overcome the limitation of equal ASLR at all TGs by removing from both, the TG time series and the altimetry time series associated with the TG, a linear trend equal to that measured by the altimeter. Such change of variables (COV) does not alter the statistical properties of the TG and altimetry time series but eliminates any difference in relative sea level changes due to different absolute sea level changes. For it to work two assumptions are necessary:

1. All the time series have a linear trend in every period in which they are considered.

2. The absolute sea level rates observed by satellite altimetry in its era can be extended backward in time to cover the timespan of the associated TG relative sea level time series.

While the first assumption can be easily verified by visual inspection or with more precise statistical methods, as the goodness-of-fit R² test [59], the second assumption, needed to permit the third change of variable in Eq. (12), can be more difficult to assess. In general, the linearity of a TG’s RSLR trend can partly corroborate the validity of the second assumption, as the probability that two different, non-linear trends of the local ASLR and VLM perfectly combine by chance, to give an overall linear trend, is obviously low.

7. Case-study of Venice and Adriatic Sea

The method of derivation of VLM described in the previous section will now be applied to a real case. To this end we have chosen the Adriatic Sea for its complexity and for the interest in this area. Indeed, several historical heritage cities and commercial/productive sites lie in the coastal area of the region, not to mention the number of people leaving along the Adriatic Sea coast, which at the end of last century was already higher than 3.5 million [60]. First, we will derive the VLM values in the Adriatic Sea using the classical LIPWC technique. After that, the LIPWC method will be applied to the same data using the change of variable presented in Eq. (12), and the results of the two strategies compared.

The TG in the Adriatic Sea for which long time series of monthly sea level are available are few. Table 1 reports their name, position, and data availability. The records of the TGs have been formed in some case from different sources, principally the Permanent Service for Mean Sea Level (PSMSL) [61], the Venice Tide Forecast and Early Warning Center (Centro Previsioni e Segnalazioni Maree, CPSM) of Venice Municipality, the Istituto Superiore per la Protezione e la Ricerca Ambientale (Italian Institute for Environmental Protection and Research (ISPRA)) and the Institute of Marine Sciences of the National Research Council of Italy (CNR-ISMAR).

Location	TG Name	Lat; (° N)	Long; (° E)	Data; (%)	Time Span; (Years)	Record Length; (Years)
Venice	VENEZIA *	45°25′51.45″	12°20′13.39″	97	1872–2018	147
Venice off-shore	VEPTF	45°18′51.29″	12°30′29.69″	100	1974–2018	45
Trieste	TRIESTE *	45°38′50.00″	13°45′33.90″	89	1875–2018	145
Rovinj	ROVINJ	45°05′01.18″	13°37′44.86″	99	1955–2018	64
Split	SPLIT *	43°30′23.88″	16°26′18.44″	100	1952–2018	67
Dubrovnik	DUBROVNIK	42°39′28.40″	18°03′38.84″	99	1956–2018	63

Table 1.

Principal characteristics of tide gauges considered in Adriatic Sea. Some of the tide gauge records have been formed by collating partial records from different sources. Such situation is marked by an asterisk after the tide gauge name.

Figure 2 shows the position of the TGs on the map of the Adriatic Sea region. Some of the TG records have been formed by collating partial records from different sources. Such TGs are marked by an asterisk. The individual positions of the TGs with respect to the twelve closest nodes of the.

Figure 2.

Positions of six tide gauges in Adriatic Sea. Color bar indicates length of available time series of sea level at tide gauges; shortest time series is about 50 years.

C3S altimetry grid are shown in Figure 3: note that some of the grid nodes are represented over land. This is an artifact of the gridding procedure that partially extrapolates over land the SLA field [36].

Figure 3.

Geographical location of the six tide gauges and position of the twelve nearest grid points of C3S altimetry SLA. Tide gauges are marked by black squares. Altimetry grid nodes are red dots. Blue triangle marks the grid node with best correlation match. Also shown (green circle) the lowest root mean square difference of the TG and altimetry monthly time series.

VEPTF is the shortest record in the set, as it started sea level recordings only in 1974. Nonetheless, its length is almost double that of the altimetry era, and abundantly double the period of the lunar nodal tide. To treat evenly all the TG records, we consider in situ sea level data from 1974 up to 2018 for all the TGs.

Plots of the in situ, as well as of the altimetry sea level anomaly monthly means observed at the six locations in the Adriatic Sea are reported in Figure 4: the seasonal and tidal signals have been removed from both the in situ and the altimetry datasets. The altimetry grid node associated to the TG time series has been chosen as the one whose time series has the higher correlation coefficient with the sea level time series of the TG, among the twelve grid nodes closest to the latter. All the sea level trend errors have been calculated considering serial correlation and are given with a 95% confidence interval.

Figure 4.

Plots of the sea level anomalies registered by tide gauges (1974–2018, in blue) and observed by the satellite altimetry (1993–2018, C3S dataset, in orange).

The altimetry dataset used to represent sea level anomaly in Figure 4 is C3S. The in situ and the remotely sensed sea level records are in good agreement, as the lowest Pearson’s correlation coefficient between altimetry and TG sea level time series is 0.82 at the Rovinj station, while all the others reach values larger than 0.91.

However, in some period a marked difference between in situ and altimetry SLA are seen, as for example in VENEZIA during 2012–2019 (TG sea level higher than altimetry), which is also confirmed by the nearby TG of VEPTF and seems to interest in a lesser extent also TRIESTE and DUBROVNIK, and for ROVINJ in 2014–2015 (TG sea level lower than altimetry), and also SPLIT in 2002–2005 (TG higher than altimetry). On the other hand, common patterns are identified in all the records throughout the observation period.

Global positioning system (GPS) observations are synergistically included in our analysis. Several sources of GPS data, at different elaboration levels, are currently available online for geocentric surface velocity data and trends from continuous GPS (CGPS) stations at selected locations, in particular near TGs: Système d’Observation du Niveau des Eaux Littorales (SONEL)/Université La Rochelle (ULR) [62], Nevada Geodetic Laboratory (NGL, University of Nevada) [63]. Other sources of information are available from local and national public agencies: in this study we used for the VENEZIA TG station also data acquired and processed by ISPRA at the PSAL tide gauge [64], which provides a relevant part of the VENEZIA sea level record. Table 2 reports the vertical velocities registered at five of the six tide gauges considered in this study, with their time span and the values provided by one or more centers for the same TG by one or more GPS stations nearby.

GPS; station	NGL (mm y^−1); Span (yr)	SONEL (mm y^−1); Span (yr)	ISPRA (mm y^−1); Span (yr)	Distance from TG (Km)	Pooled Mean; (mm y^−1)
PSAL; VENEZIA	−1.70 ± 0.80; 2014–2020	—	−1.46 ± 0.09; 2010–2015	0.01	−1.59 ± 0.65; 2014–2020
TRIE; TRIESTE	−0.52 ± 0.45; 2003–2020	0.20 ± 0.26; 2003–2013	—	6.97	−0.25 ± 0.52; 2003–2020
PORE; ROVINJ	−1.51 ± 1.03; 2011–2021	—	—	16.62	−1.51 ± 1.03; 2011–2021
SPLT; SPLIT	0.45 ± 0.68; 2004–2012	−0.25 ± 0.34; 2004–2012	—	< 4.00	0.10 ± 0.64; 2004–2012
DUBR+DUB2; DUBROVNIK	−1.83 ± 0.70; 2000–2020	—	—	4.15	−1.83 ± 0.70; 2000–2020

Table 2.

Geocentric surface vertical velocities at three locations in the Adriatic Sea from GPS stations.

PSAL is almost co-located with the VENEZIA PUNTA DELLA SALUTE TG. For SPLIT, data from the CGPS station of SPLT were acquired. It is worth mentioning that SPLT is, with PSAL in VENEZIA, among the few CGPS co-located with TGs in the Adriatic Sea. The TRIE CGPS station is the nearest to the TRIESTE TG, but 6.9 km far from it, over a hill north-west of Trieste: for this reason, TRIE CGPS station cannot be considered co-located with the TRIESTE TG. Neither can the PORE CGPS station for ROVNIJ be considered as such, and the DUBR and DUB2 CGPS stations in DUBROVNIK: PORE is located 16 km north of Rovinj along the coast, while DUBR and DUB2 are located 4 km away and 400 m in height.

With the data described so far, the VLM can be derived with the classical method, i.e. subtracting the TG RSLR from the ASLR observed by altimetry at the associated grid point. This approach, described in [43], allows to estimate the VLM separately at each location for which RSL and ASL records are available. The error associated to these estimates is drastically reduced when the linear trend of VLM is calculated by differencing the time series of the ASL and RSL, instead of combining the two errors of ASLR and RSLR as they were two independent measurements. From here on, all the errors on the sea level change rate are calculated according to this convention. The results of such approach are shown in Table 3: in column 1 appear the TG locations, in column 2 the ASLR derived by altimetry, in column 3 the RSLR derived by the TG, and in column 4 the VLM (${\dot{u}}_i={\dot{g}}_i-{\dot{s}}_i$) derived by differencing the time series of ASL and RSL monthly time series.

Location	$\dot{g}$; (mm y^−1)	$\dot{s}$; (mm y^−1)	$\left(\dot{g}-\dot{s}\right)$; (mm y^−1)
VENEZIA	3.36 ± 1.45	5.15 ± 1.73	−1.79 ± 0.65
VEPTF	3.38 ± 1.46	5.50 ± 1.73	−2.12 ± 0.67
TRIESTE	3.75 ± 1.58	3.56 ± 1.66	0.18 ± 0.60
ROVINJ	3.33 ± 1.58	1.03 ± 1.85	2.30 ± 1.06
SPLIT	3.60 ± 1.36	2.92 ± 1.65	0.68 ± 0.63
DUBROVNIK	3.34 ± 1.22	3.79 ± 1.48	−0.45 ± 0.55

Table 3.

Results of calculations using C3S altimetry dataset (1993–2018). Column 1 reports the TG location; columns 2 and 3 the absolute and relative sea level rates in the altimetry era; column 4 the VLM calculated with the classical approach (ALT-TG). All data are in mm y⁻¹.

First of all we note that the error of the VLM estimates in the fourth column, obtained as standard error of the trend of the differenced time series (ASL-RSL) are much lower than that provided by the error propagation formula for the difference of the trend estimates of two statistically-independent time series, as in this case the error propagation formula would provide ${\sigma }_{\mathit{gs}}=\sqrt[2]{{\sigma }_g^2+{\sigma }_s^2}$, and in the VENEZIA case, for example, it would determine a standard error of 2.26 mm y⁻¹ instead of the 0.65 mm y⁻¹ resulting by calculating the trend and the standard deviation of the differentiated time series.

A second aspect worth to note is the independence of each VLM determination from all the others. That means if one of the VLM estimates is affected by large errors or relies on data of bad quality (RSL and/or ASL), it does not influence the evaluation of the others.

The third observation about the numbers reported in Table 3 is that while the ASLR is almost constant at all sites of the Adriatic Sea considered in this study, the RSLR observed at the TGs are much more varied, determining VLM estimates going from −2.12 to +2.30 mm y⁻¹. From one side, this means that the vertical velocity field applicable to the Adriatic area is not constant, potentially revealing that different processes could be at the base of the observed crustal motions. From the other side, such numbers reveal also that the VLM is an essential parameter in sea level studies conducted mainly from tide gauge data. Thus, every methodology able to estimate the VLM at the TG is of extreme interest to correct the RSL observed at the TGs themselves, in particular where no geodetic measurements are available to estimate the VLM.

To conclude this section about the classical approach to VLM estimate from sea level data, we report a comparison of the two gridded altimetry datasets: C3S and SLCCI. To provide a fair comparison, both datasets have been limited to the same common period of temporal coverage: 1993–2015. The results of the classical approach to VLM estimate are given in Table 4.

Location	$\dot{g}$ C3S; (mm y^−1)	$\dot{g}$ SLCCI; (mm y^−1)	$\dot{s}$; (mm y^−1)	$\left(\dot{g}-\dot{s}\right)$; C3S; (mm y^−1)	$\left(\dot{g}-\dot{s}\right)$; SLCCI; (mm y^−1)
VENEZIA	4.16 ± 1.71	4.47 ± 2.07	6.08 ± 2.08	−1.93 ± 0.79	−1.61 ± 0.91
VEPTF	4.17 ± 1.73	4.47 ± 2.07	6.44 ± 2.07	−2.27 ± 0.81	−1.97 ± 0.93
TRIESTE	4.51 ± 1.84	3.42 ± 1.78	4.49 ± 1.98	0.02 ± 0.71	−1.07 ± 0.74
ROVINJ	4.09 ± 1.85	4.37 ± 1.86	1.91 ± 2.04	2.18 ± 1.17	2.46 ± 1.04
SPLIT	4.44 ± 1.57	4.15 ± 1.48	4.15 ± 1.87	0.29 ± 0.70	0.00 ± 0.76
DUBROVNIK	4.01 ± 1.42	3.98 ± 1.45	4.67 ± 1.70	−0.66 ± 0.62	−0.69 ± 0.70

Table 4.

VLM estimates from in situ RSL and remotely observed C3S and SLCCI ASL 1993-2015.

Column 4 of Table 4 reports the RSLR, which is common to both ways to calculate VLM, the classical and the LIPWC.. In columns 2 and 3, differences can be seen in the ASLR measured by C3S and SLCCI: the most notable refers to TRIESTE, which appears to observe an ASLR of 4.51 mm y⁻¹ in the C3S dataset, and 3.42 mm y⁻¹ in the SLCCI dataset: these numbers differ by more than 1 mm y⁻¹. The difference in ASLR for TRIESTE is reflected in the final VLM rate. Regardless the marked difference for TRIESTE, the other rates appear in good agreement between the two datasets, even if in general C3S supplies lower errors.

So far, we have shown the results of the classical approach to VLM determination from altimetry and tide gauge. From now on we present and analyze the results of the linear inverse problem with constraints, in the modified version which exploit a change of variable to disentangle the contribution of the ASLR from the system. To do so we examine only the results relative to the C3S dataset, for three reasons: first of all, the final results do not differ much between the two datasets; second, the C3S gridded product has an enhanced resolution in the Mediterranean Sea, and an appropriate regional processing; third, the C3S dataset has a time span longer than SLCCI, and most important, it is intended to be continuously updated in the future. VLM results derived with the LIPWC-COV approach are shown in Table 5, together with the values of the ASLR and RSLR values used for the calculation, and the VLM derived with the classical approach for ease of comparison.

Location	$\dot{g}$; (mm y^−1)	$\dot{s}$; (mm y^−1)	$\dot{\zeta }$; (mm y^−1)	$\left(\dot{g}-\dot{s}\right)$; (mm y^−1)	LIPWC-COV; (mm y^−1)	GPS; (mm y^−1)
VENEZIA	3.36 ± 1.45	5.15 ± 1.73	3.26 ± 0.73	−1.79 ± 0.65	−0.93 ± 0.39	−1.59 ± 0.65
VEPTF	3.38 ± 1.46	5.50 ± 1.73	3.78 ± 0.73	−2.12 ± 0.67	−1.41 ± 0.47	—
TRIESTE	3.75 ± 1.58	3.56 ± 1.66	2.30 ± 0.67	0.18 ± 0.60	0.42 ± 0.33	−0.25 ± 0.52
ROVINJ	3.33 ± 1.58	1.03 ± 1.85	1.36 ± 0.71	2.30 ± 1.06	0.93 ± 0.37	−1.51 ± 1.03
SPLIT	3.60 ± 1.36	2.92 ± 1.65	2.20 ± 0.66	0.68 ± 0.63	0.37 ± 0.33	0.10 ± 0.64
DUBROVNIK	3.34 ± 1.22	3.79 ± 1.48	2.69 ± 0.58	−0.45 ± 0.55	−0.41 ± 0.46	−1.83 ± 0.70

Table 5.

VLM results using C3S altimetry dataset (1993–2018). Location in column 1; ASLR in column 2; RSLR in the altimetry era (1993-2018) in column 3; RSLR 1974-2018 in column 4; VLM calculated with the classical approach in column 5, as in column 4 of Table 3, and with the LIPWC-COV in column 6. Columns 7 reports the VLM values directly detected by the GPS stations associated with three TGs.

The difference between the results obtained in the classical approach and the LIPWC-COV approach is evident; while the classical approach range of VLMs is [−2.12 2.30] mm y⁻¹, that provided by LIPWC-COV is almost half as wide: [−1.41 0.93]. This is to the result of the introduction of the constraints in Eq. (7), which enter the linear system, propagating the structure and values of the relative vertical motion between the TGs in the solution. The effect of the constraints is best seen in Figure 5, where the plot of the LIPWC-COV solution follows the general form of the classical solution, but with a reduced spread.

Figure 5.

Scatterplots of VLM values derived with the classical $\left(\dot{\mathrm{g}}-\dot{\mathrm{s}}\right)$ and the LIPWC-COV approaches using the C3S altimetry dataset (period 1993–2018). GPS estimates (in black) are also reported. Results from the study of Wöppelmann and Marcos (W&M) for the period 1992–2010 are shown in green for comparison. The zero level is drawn in black. (adapted from [58]).

In Figure 5 are reported also the VLM values measured by the CGPS stations of PSAL (VENEZIA), TRIE (TRIESTE) and SPLT (SPLIT), and the values of VLM estimated by Wöppelmann and Marcos [50] with the LIPWC technique without the change of variable.

In the classical approach, as there is no optimization of errors as in the LIPWC technique, we see a wide spread of the VLM values. This is particularly evident for ROVINJ TG, whose $\left(\dot{g}-\dot{s}\right)$ estimates reach more than 2 mm y⁻¹, while the LIPWC-COV approach calculates it as less than 1 mm y⁻¹. The LIPWC solution proposed by Wöppelmann and Marcos [50] presents much lower standard errors than LIPWC-COV solution described in this study. We presume that such low standard errors were attained by a different methodology in calculating the rates of absolute and relative sea level change rates and their formal errors. Moreover, in the years following 2010 the VLM rates at the five common TG locations have remained substantially unmodified with respect to the Wöppelmann and Marcos’ results. As a final step regarding VLM, we have calculated the root mean square difference (RMSD) of the VLM calculated with GPS and those calculated with the classical and the new (LIPWC-COV) approaches and found that the second one is lower: Classic approach: 1.84 mm y⁻¹; LIPWC-COV approach: 1.34 mm y⁻¹.

The discrepancy observed between this study and that of Wöppelmann and Marcos can largely be ascribed to the different periods covered by the altimetry datasets (C3S and SLCCI datasets cover time periods respectively 44% and 23% longer than the study of Wöppelmann and Marcos). Other factors that may contribute to explain the difference between the results of the two studies are the processing of the altimetry data and the inclusion of the VEPTF TG in this study. The rates of absolute sea level change at the TGs, calculated as the sum of relative sea level change and VLMs derived in this study with the LIPWC-COV approach, for the whole period covered by the TG record, are reported in Table 6.

Location	$\dot{\zeta }$-VLMLIPWC-COV
VENEZIA	2.33 ± 0.83
VEPTF	2.37 ± 0.86
TRIESTE	2.71 ± 0.75
ROVINJ	2.29 ± 0.80
SPLIT	2.57 ± 0.74
DUBROVNIK	2.28 ± 0.74
Pooled mean	2.43 ± 0.80
Sample mean	2.43 ± 0.18

Table 6.

ASLR from TG records over whole period 1974–2018, corrected for VLM estimated with the LIPWC-COV approach. All data are in mm y⁻¹.

The uncertainty of the sample mean (last row of Table 6) was obtained as standard error of the sample mean, considering the rates as random and independent variables. The absolute sea level change rates vary in a very narrow interval, 2.33–2.71, with a sample mean of 2.43 mm y⁻¹. The standard deviation of the sample is much lower than the precision of each individual determination of SL change rate at the TGs. As pointed out by Wöppelmann and Marcos [50], such a low dispersion is unlikely to be determined from estimates of independent random variables: it is instead the evidence of the high performance of LIPWC method for determining accurate VLM rates from TG and altimetry differenced time series. The ASLR rates calculated by altimetry in 1993–2018 and through the LIPWC-COV technique (1974–2018) are shown in Figure 6.

Figure 6.

Absolute sea level change rates as calculated by altimetry 1993–2018, and by the LIPWC-COV approach integrating data from TGs in 1974–2018.

Clearly, the ASLR values calculated for the longer period are smaller than those calculated in the shorter one, but the modulation of the rate from TG to TG is apparently reflected in the LIPWC-COV approach. As already noted, the errors associated to the ASLR rates derived in the LIPWC-COV are also smaller, thanks to the introduction of the constraints on the relative vertical land motion between paired TGs. The mean value of the ASLR calculated for the Adriatic Sea with the LIPWC-COV approach, is in general agreement with both regional studies on the Mediterranean Sea (0.7 ± 0.2 mm y⁻¹ (1945–2000) [65]; 1.60 ± 0.35 mm y⁻¹ (1992–2010) [50]; 2.44 ± 0.5 mm y⁻¹ (1993–2012) [66]; 2.87 ± 0.33 mm y⁻¹ (1992–2016) [67]), and at global scale (2.0 ± 0.3 mm y⁻¹ (1971–2010) [68]; 3.0 ± 0.7 mm y⁻¹ (1993–2010) [69]; 2.8 ± 0.5 mm y⁻¹ (1993–2010) [68]; 3.1 mm y⁻¹ (1992–2018) [32]).

Among the ASLR altimetry rates associated with the six TGs in the Adriatic Sea, those for TRIESTE are very different in the C3S and SLCCI dataset. In order to investigate such a large difference (0.02 ± 0.71 mm y⁻¹ SLCCI; −1.07 ± 0.74 mm y⁻¹ C3S; see Table 4) the SLCCI-AT X-TRACK/ALES 20 Hz along-track coastal altimetry dataset has been used.

The analysis focuses on the descending track 196 of the Jason-1 (2001–2013) and Jason-2 (2008–2019) altimetry missions and covers the period 2002–2016 with 532 cycles (from 22 Jan 2002 to 23 Jun 2016 at 350 m resolution along-track) with MSSH computed using cycles from 1 to 517. The position of the track 196 and the geographical setting are shown in Figure 7.

Figure 7.

Gulf of Trieste. The positions of the SLCCI along-track 20 Hz altimetric product version 1.0 samples of the descending track 196 (white circles), and the TRIESTE TG stations (red triangle). Umag in Croatia and Grado in Italy are also shown (green circles).

Altimetry data at the 71 observation points of track 196 are compared to 10′ interval RSL observations of the TRIESTE TG. The TG time series did not undergo any filtering or processing, and the astronomical tide and Dynamic Atmospheric Correction (DAC) corrections are not applied to the altimetry time series.

The goal of the investigation is to explore the possible causes of the different ASLR rates obtained by the two gridded altimetry datasets near Trieste, to look for clues directly into the original along track data from the Jason missions, reprocessed with advanced and coastal specific re-tracking (ALES) and improved coastal processing (X-TRACK). We also want to ascertain the suitability of the new SLCCI-AT record in long term coastal sea level monitoring. We concentrated on the altimeter track 196 of the SLCCI-AT dataset, which first crosses Marano Lagoon and a 0.5 km wide sandbar before entering the Gulf of Trieste from north, near Grado, and then flies over Umag and the full extent of the Istria peninsula. The retrieval is particularly problematic in the gulf area due to the complex morphology of the land. Moreover, some data loss could be due to sea-to-land and land-to-sea crossings that might influence the behavior of the on-board tracker. Operational altimetry products do not provide data over this section of the Gulf of Trieste, while the SLCCI-AT dataset provides 71 points along track, most of which yield over 70% of valid data (blue box in Figure 8). The most improvement is near the Istrian peninsula with more than 90% of data recovered. The valid data percentages decrease abruptly over a distance ranging 5 km from the coast. The reduced performance over the lagoon and islet (almost all data have been rejected) is probably related to the data corruption in the land-sea-transition. Note that at 1 Hz, any coastal altimetry along track product would give no more than 3–4 points along this 24 km long stretch of track 196.

Figure 8.

J1 + J2 track 196 geographical settings. Left: Percentage of valid data along the track. Right: Correlation with TRIESTE tide gauge. Adapted from the coastal sea level project of the ESA climate change initiative (SL_CCI bridging phase) document “Part II: Validation Results” (http://www.esa-sealevel-cci.org/webfm_send/588).

The data accuracy can be assessed in more detail comparing the altimeter-derived 20-Hz SLA with corresponding tide gauge sea level measurements. It should be noted that the TRIESTE TG is located in the harbor, and therefore it does not measure exactly the same ocean dynamics as the altimeter flying offshore. Nonetheless, the Pearson’s linear correlation coefficient of most of the 71 points along the section of track 196 facing the Trieste harbor exceeds 0.9 (red box in Figure 8). The RMS difference between altimetry observations and tide gauge measurements of instantaneous sea level is almost constant along the track 196 section in the Gulf of Trieste, and around 10 cm (not shown).

From the time series of SLCCI-AT SLA at each data point of the track 196 facing Trieste, we have calculated the slopes of the fitting lines, gradually growing the confidence interval from 68–95%, and performed Mann-Kendall statistical significance tests [70, 71] modified for autocorrelated data [72] on all the 71 fitting lines. The Mann-Kendall test is commonly employed to detect monotonic trends in time series. The null hypothesis is that the data come from a population with independent realizations and are identically distributed. The alternative hypothesis is that the data follow a monotonic trend. In Figure 9 the results of such calculation are reported for a preliminary version of the SLCCI-AT dataset at 20 Hz. The black diamonds mark the acceptance or rejection of the null hypothesis following this scheme:

• 1 - the null hypothesis “the sample has no trend” is rejected.

• 0 - the null hypothesis “the sample has no trend” cannot be rejected.

Figure 9.

Slopes and slop errors of the lines fitting the time series of along track SLA at every data point of track 196 near Trieste. Also plotted the Mann-Kendall test results. Black diamonds: 1 - rejection of the null hypothesis (the sample has no trend); 0 - no rejection. Green line: p_value. Adapted from the coastal sea level project of the ESA climate change initiative (SL_CCI bridging phase) document “Part II: Validation Results” (http://www.esa-sealevel-cci.org/webfm_send/588).

Already with a 68% confidence interval the null hypothesis (the trend is not statistically significant) is rejected in less than 24% of the data points. With a 95% confidence interval only for four fitting lines out of 71 the null hypothesis is rejected. In both cases the errors associated to the slopes are higher than the slopes themselves.

A similar analysis replicated on the final version of the SLCCI-AT dataset, published at the end of the SLCCI project, gave better results. Figure 10 reports the representation of the statistical characteristics of the slopes derived from the last version of the data of the SLCCI-AT X-TRACK/ALES SLA 20 Hz, with 95% confidence interval. The left panel shows slopes and associated errors at every data point latitude (low latitudes are near Umag, high latitudes near Grado); different colors indicate the statistical significance of the Mann-Kendall test (blue: significant; orange: not significant). The right panel shows the box and whisker plots of the two distributions (left: not significant; right: significant). The number of statistically significant slopes is much higher in the final version of the dataset, even if the variability is still rather high and difficult to explain because of the limited spatial variability along the track. Slopes are higher towards north (Grado), and lower near Umag. Considering only the statistically significant slopes in the SLCCI-AT dataset, their sample mean and standard deviation result to be 3.40 ± 1.01 mm y⁻¹ (Feb-2002 – Jun-2016) which is not far from the trends we have found in the Adriatic Sea at all the tide gauges. We recalculated the altimetry trends near TRIESTE in the SLCCI and C3S gridded products. The altimetry ASLR trends found so far in the analysis are summarized in Table 7.

Figure 10.

SLCCI SLA 20 Hz. Left: Slopes and slope errors of the lines fitting every data point of the track 196 in the Gulf of Trieste. Blue: Statistically significant slopes according to the Mann-Kendall test. Orange: Slopes not significant. Right: Box and whisker plots for the statistically significant and non-significant slopes. Red: Median value. Box: Upper and lower quartiles. Whiskers: Highest and lowest observations. Adapted from the coastal sea level project of the ESA climate change initiative (SL_CCI bridging phase) document “Part II: Validation Results” (http://www.esa-sealevel-cci.org/webfm_send/588).

$\dot{g}$ TRIESTE 2002–2016
SLCCI-AT; (mm y^−1)	SLCCI; (mm y^−1)	C3S; (mm y^−1)
3.40 ± 1.01	3.66 ± 3.97	5.07 ± 3.64

Table 7.

Trends for Trieste in February 2002 – June 2016 from SLCCI project and C3S altimetry. Column 1: SLCCI-AT along track 20 Hz product. Column 2: SLCCI gridded product. Column 3: C3S gridded product.

The trends calculated with the SLCCI dataset (along track and gridded) are in good agreement, apart from the different errors affecting the two results, due to the different methods used to calculate them. The C3S trend is instead higher than the other two. We believe that the difference between the results is to be ascribed to the different methodologies used in the two products. In any case the difference between the SLCCI and the C3S results is not yet explained by this further analysis, and the Gulf of Trieste remains a controversial place for the derivation of climatologically relevant oceanic variables from altimetry, because of the proximity of the land and the geometry of the surrounding coastline, and the very short time coverage of the altimetric datasets.

8. Summary and prospects

The sea level is a key variable of the climate system. Tide gauges measuring sea level variability are in operation since the 1900s. Satellite-based observations of sea level changes are more recent. Nevertheless, they play a crucial role in understanding the future coastal sea level changes. Advance in the processing of satellite radar altimetry have expanded the utility of this data set for climate-related studies and extended the potential exploitation in the coastal zone. The joint usage of the two different measuring systems (in situ and satellite) has two challenges. First how the two data sets can be consistently and systematically used in synergy to address that objective of estimating robust coastal sea level trends. Second how using high-rate (i.e. 20 Hz) altimeter measurements with a coastal-oriented processing could improve the satellite-based trend estimates with respect to the standard (1 Hz) data, especially near coast.

In this chapter, a more robust inverse method (called LIPWC-COV) has been proposed and tested in the Northern Adriatic Sea, where GPS data are available to conduct a realistic assessment of uncertainties. The results show that the classical approach of estimating VLMs provides less accurate trends than the LIPWC-COV method, and with lower errors. Moreover, the LIPWC-COV has demonstrated to compare better than the classic method with GPS derived VLMs.

In this chapter, the experimental SLCCI data set (high resolution along track) coastal sea level product (developed within SLCCI project) has been also assessed in the Gulf of Trieste, as it was possible only at that site. The retrieval is particularly problematic in the gulf area due to the complex morphology of the land. The trends calculated with the gridded and along track datasets show some differences, probably due to the different methodologies used in the generation of the products.

This study offers a more consolidated and improved understanding of the sea level trend variability in the Northern Adriatic Sea. The next step is to extend the application of the new methodology to the Mediterranean Sea.

Acknowledgments

The authors want to thank the European Space Agency, which propelled the Climate Change Initiative to produce a climate quality record of sea level from satellite altimetry. Likewise, the authors express gratitude to the European Union for making available the satellite altimetry data record through the Copernicus Earth Observation Programme Services. The authors also want to acknowledge Centro Previsioni e Segnalazioni Maree (CPSM) and Alvise Papa for providing tide gauge data under the CNR-ISP/CPSM technical–scientific collaboration agreement 2019–2022; the CNES distribution service AVISO+ for the dynamic atmospheric correction data; the Italian Institute for Environmental Protection and Research (ISPRA) for providing tide gauge and GPS data; CNR-ISMAR and Fabio Raicich for the time series of the Trieste tide gauge; PSMSL, national Oceanographic Centre, Liverpool, UK, for making available the tide gauge global record. The authors thank the SONEL Data Centre and the Nevada Geodetic Laboratory (NGL) for processing and making available the GPS trend global record. Finally, we acknowledge CPSM for funding the “Vento da Satellite+” project on the development of satellite data in storm surge forecasting; ESA who funded the “HYDROCOASTAL” project (ESA contract No. 4000129872/20/I-DT), and the “SL_CCI Bridging Phase” project (ESA contract No. 4000109872/13/I-NB), from which some of the images in this chapter were adapted.