Layout 6 1 ANNALS OF GEOPHYSICS, 63, 6, GD656, 2020; doi:10.4401/ag-8403 Point wise Determination of Vertical Deformation at southern Coastal Areas of Caspian Sea by Combination of Satellite Altimetry and Tide Gauge Observations Mahmood Pirooznia1, Mehdi Raoofian Naeeni*,1, Taghi Yousefzadeh2 (1) Department of Geodesy, Faculty of Geodesy and Geomatics Engineering, K. N. Toosi University of Technology (2) Department of Geodesy, School of Surveying and Geospatial Engineering, College of Engineering, University of Tehran Article history: received December 7, 2019; accepted April 27, 2020 Abstract In this paper, a combination of satellite altimetry (SA) and tide-gauge (TG) observations is used to determine the pattern of vertical deformation at southern coasts of Caspian Sea. Satellite altimetry measures the variations of sea water surface with a respect to a reference ellipsoid, while tide gauge data gives the sea level relative to a particular ground station; consequently, the difference between these two observations will give the vertical motion at the ground station. However, altimeter satellites still have problems at shallow waters or coastal areas which leads to inaccurate and corrupted sea surface height at tide gauge stations. To minimize this imprecision, two strategies are applied for the selection of the suitable SA points near the TG stations. The First one is to choose those close passes to the TG stations, provided that the SA data is not corrupted (range observation is available). The second strategy is to extrapolate the SA points at TG stations using Spatio- Temporal Kriging method. In both point of views, use is made of re-tracking methods to improve the ranges of SA observations near the coast. Amplitudes of tidal components are removed from both SA and TG observations to derive the residual tide free signals. Afterward, a linear trend is fitted to the residual signals in the sense of least square; the difference between these two, supplies the vertical deformation of the ground station. The results illustrate that the average vertical crustal motion at Anzali, Noshahr and Neka are -1.387, -1.903 and -3.14 mm/year, respectively. Keywords: Vertical Crustal Motion (VCM); Satellite altimetry; Coastal TGs; Spatio-Temporal Kriging; Caspian Sea. 1. Introduction Caspian Sea is encircled by five countries, namely Islamic Republic of Iran, Russia, Republic of Azerbaijan, Turkmenistan, and Kazakhstan. It is the largest extant part of the old Tethys Sea; an ancient sea which had extended from North Pole to Indian Ocean in the first three eras of geology. Disconnection of Tethys occurred in the third period due to geological folding and resulted in the formation of famous Caucasus and Minor Asia mountain ranges; consequently, an uplift happened in European continent and Iran Plain. This in turn, eventuated in emergence of Caspian Sea. As shown in Figureure1, the Caspian Sea is situated north of the equator, from 36°33’ to 47°57’ north latitude and 46°43’ to 54°53’ east longitude. This geographical range covers an approximate length and width of 1030 and 240 km respectively [Arpe et al., 2000]. The structural tectonic features in the South Caspian Sea region are particularly important for geodynamical studies. These are mainly due to the active tectonics and structure of the South Caspian Sea and its role in the collision between Arabia and Eurasia [Jackson et al., 2002; Brunet et al., 2003]. The South Caspian depression is an autonomous tectonic unit that consists of: the southern part of the Caspian Sea; the West Turkmen depression adjacent to the latter; the northern lowland periphery of the Alborz system; the Lower Kura trough; and the Apsheron–Kobystan trough [Mammadov, 2008]. This interaction between tectonic plates makes the region subject to frequent and strong earthquake. The South Caspian Basin is surrounded by active earthquake belts on all sides. The southern Caspian Sea stands out as an aseismic block surrounded by belts of intense earthquake activity. These phenomena cause deformations in the crust of the Earth [Jackson et al., 2002]. A large number of geodynamic studies based on crustal deformations are done by using observations that typically derived from absolute gravimetry [Pagiatakis and Salib, 2003], global positioning system (GPS) [Larson and Van Dam, 2000; Johansson et al., 2002; Caccamise et al., 2005; Beavan et al., 2010], satellite laser ranging (SLR) [Noomen et al. 1996], Doppler orbitography and radiopositioning integrated by satellite (DORIS) [Soudarin et al., 1999], interferometric synthetic aperture radar – InSAR [Dalla et al., 2012], and very long baseline interferometry (VLBI) [Tesmer et al., 2009]. Most of the crustal deformation studies are based on Global Navigation Satellite System (GNSS) observations due to high levels of accuracy, moderate associated costs, and broad spatial coverage in a relation to SLR, VLBI and other techniques [Blewitt, 2007]. However, along the South of Caspian Sea, particularly coastal region, a dense, long term time series and homogenous GPS geodynamic network do not currently exist, and GPS geodynamic stations are located far from the coastal region. Consequently, studying geodynamic phenomena in south coastal area of Caspian Sea presents difficulties. Therefore, other alternatives are necessary. SA and TG records provide alternative and independent methods for estimating vertical crustal motion [Cazenave et al., 1999; Lin, 2000; Nerem and Mitchum, 2002; Fenoglio-Marc et al., 2004]. Due to strategic and economic importance of coastal areas and their impact on local people’s lives, dynamical changes of these areas and their sources have captured the interests of scientific communities in the recent years. Multitude of factors such as meteorological effects, oceanic influences, tide, climate change and vertical crustal motion control the dynamics of these regions [García et al., 2007; Pirouznia et al., 2017]. Among these factors, the considerable effect that vertical crustal motion (VCM) has on TG measurements, makes it a subject of interest in geophysical studies. For this reason, hitherto, so much effort has been made to estimate this motion [see Clark et al., 2002; Nerem and Mitchum, 2002]. Tide gauges measure the relative sea level changes, which is the variation of the mean sea surface with a respect to the solid Earth. To obtain such measurements, TG’s must be referenced to a benchmark on the crust of the Earth. However, in the most parts of the world, locations of these benchmarks are subjected to the Earth’s deformation processes such as earthquakes and their deformation cycles, tectonic plate motion, basin evolution, Glacial Isostatic Adjustment (GIA), and anthropogenic effects such as local subsidence due to water extraction, mining, oil or natural gas drilling. Therefore, the signals of TG records contain both sea level changes and VCMs [Kuo et al., 2008]. Observations of SA consist of (i) the precise satellite location in a relation to a reference ellipsoid, and (ii) the height of satellite over the surface measured using laser or radar pulses which are sent by satellite and bounced back from the target surface. Subtracting (ii) from (i) gives the sea surface height (SSH) relative to the reference ellipsoid at the given location and time. On the other hand, the TG measurement gives the sea level height relative to the solid ground on which the TG is fixed. Hence any vertical ground motion, regardless of its source, will be interpreted as an illusory change in sea level height observed by TG. As the land under the benchmark moves up and down, an apparent fall or rise will be seen in sea level which is due to the tectonic movement of sea floor. In other words, TG measures both pure sea level fluctuations which are variations in the water content and sea floor deformation which is due to tectonic activity. It is clear that the vertical ground motion can be obtained by subtracting these two data types, namely, altimetry and TG [Nerem and Mitchum, 2002; Garcia et al., 2007; Braitenberg et al., 2011]. Shum et al. [2002], used historic TG records and a decade of TOPEX/POSEIDON (T/P) and Geosat altimeter data in a stochastic adjustment to compute the estimation of vertical motion on the Great Lakes regions. In other paper, Mahmood Pirooznia et al. 2 Fenoglio-Marc et al. [2004] calculated the rate of VCM in Mediterranean Sea, combining T/P data and a set of TG records based on sea level height changes. Another research was performed by Kuo [2004] in which he introduced a new method, namely Gauss-Markov stochastic adjustment model, and used it to estimate the vertical crustal motions (VCMs) in areas near to open seas, semi-enclosed seas and lakes. Kuo et al. [2008] combined geocentric sea level measurements from high-precision T/P satellite altimetry observations and TG records at 25 stations distributed in coastal regions of the Baltic Sea and around the Great Lakes. Ostanciaux et al. [2012] computed VCM rates by combining TG records and local SA and then compared this database to previous studies that use geodetic techniques and TG records in order to evaluate the consistency of their results and previous ones. It should be pointed out that the quality of SA observations in coastal regions depends on the range between satellite orbit and a water surface, which is derived by measuring the time elapses between sending and receiving radar pulses from the satellite board and the surface of the Earth (the altimeter waveforms) [Fu and Cazenave, 2001]. However, near the coasts, altimeter waveforms may be corrupted due to variety of reasons, and the range observations are infected by some biases [Mathers et al., 2004]. To improve the precision of those altimeter observations with complex returning waveforms, a waveform re-tracking method is implemented to reprocess the original returning pulse. “Re-tracking” refers to the analyses of the waveforms, received in the satellite antenna as pulses to improve the range [Davis, 1995; Gomez-Enri et al., 2009]. Three common re-tracking algorithms such as Offset Center Of Gravity (OCOG), �-Parameter and threshold can be used for this purpose [see Rapley et al., 1987; Martin et al., 1983; Davis, 1997; Deng, 2003; Lee, 2008; Zhang, 2009]. Based on the author’s knowledge, the re-tracking correction has not been investigated in the studies of VCM using TG and SA observations and most of the papers ignored this correction. However, this effect may be very significant as it is later shown in this paper. In this study, missions of T/P, Jason-1 and Jason-2 altimeters are used for VCM’s determination at coastal parts of Caspian Sea (see Figure 1). These missions have 10-day repeat cycles. Since VCM is determined by difference between SA and TG observations, those passes of altimeter which cross near the TG station (adjacent to the coast) should be used. Due to the problems of SA in the shallow water of coastal regions, choosing a suitable SA points near the TG stations are very important. So in this study, two strategies are implemented (i) selecting the closest SA point to the TG stations on the condition that the SA data is not corrupted. (ii) Extrapolating SA points by Spatio- Temporal Kriging method to evaluate the altimeter observation at TG stations. In both cases, however, the re-tracking algorithms are used to improve the range observations near the coast. Tide effect should also be removed from SA and TG observations. In this study, to eliminate tide phenomena, it is modeled by considering significant tidal frequencies (empirical tide modeling), and we do not use global tide model such as FES2004 for tide correction. 3 Vertical crustal motion at Caspian Sea Figure 1. Satellites crossovers above the Caspian Sea (TP, Jason1 and Jason2). Mahmood Pirooznia et al. 4 Due to sample rate of SA observations (approximately 10 days for T/P, Jason-1 and Jason-2), we encounter with aliasing phenomenon for tide modeling (see section 2.4 for more explanation about aliasing phenomenon). So To deal with it, aliasing frequencies are used for tide modeling. Eventually the rate of tide free sea level height time series of SA and TG are determined which indicate the absolute and relative changes in sea levels, respectively. Difference of the above mentioned two rates gives the VCMs. These VCMs in turn, have been compared with VCMs determined by continuous GPS stations. 2. Data source and pre-analysis of observations In this section we concisely describe our data sources and some preprocessing strategies, including re-tracking of SA observations, selection of suitable SA point and tidal analysis, being used to prepare the two data sources for determination of VCM at Caspian Sea. 2.1 Data sources Tables 1 and 2 show the features of SA and TG data which are used in this study. Sensor Geophysical Data Record (SGDR) of SA is used in this study due to containing along-track waveform information that may be used for applying re-tracking corrections for the range observations. Tide gauge data in this study is provided from Ports and Maritime Organization (PMO) of Iran which were recorded by floating pressure tide gauges. This data is not available in international databases and no processing for quality control of which were applied. Table 1. SA data used in this study. Mission Cycles Period Source TOPEX/Poseidon 001 - 365 1992/09/23 - 2002/08/11 (NASA) Jason 1 001 - 259 2002/01/15 - 2009/01/16 (NASA,AVISO) Jason 2 001 - 226 2008/07/04 - 2014/08/21 (AVISO) Table 2. Location and time duration of TGs (coastal areas of Caspian Sea). Location Latitude Longitude Duration (time) Source Anzali 37.478° 49.4623° 21/3/2005-20/3/2014 Ports and Maritime Organization (PMO) of Iran Noshahr 36.6584° 51.5047° 21/3/2006-21/3/2014 Ports and Maritime Organization (PMO) of Iran Neka 36.8502° 53.3656° 1/1/2000-31/8/2012 Ports and Maritime Organization (PMO) of Iran 2.2 Re-tracking of SA observations The observable quantity in satellite altimetry is the range observation which is known as satellite range or satellite height above sea level. The range observation is contaminated by some source of errors which must be corrected via global or local modeling and can be summarized as follows: wet troposphere correction, dry troposphere correction, ionosphere correction, electromagnetic bias correction, ocean tide correction, inverse barometer correction, pole tide correction and center of gravity movement. They are better known as “geophysical corrections” [Benada, 1997; Pirooznia et al., 2016a; Soltanpour et al. 2017]. Also the re-tracking correction is used to correct range observation. The three common re-tracking methods, namely OCOG, �-Parameter, and Threshold are implemented. These corrections should be applied on satellite range for each repeat cycle to create the time series of sea height [Khajeh et al., 2014; Khaki et al., 2014]. Thus, corrected SSHs are derived from the following equations. (1) where C is the speed of light (299,792,458 ), 𝜏 the pulse duration (3.125 ns for T/P), G0 the nominal tracking gate (for T/P, Jason-1 and jason-2, it is 32.5, 31 and 31, respectively ), G the computed tracking gate (middle point of the leading edge which is computed from re-tracking algorithms), 𝛥𝑟 the re-tracking correction, h the altitude of the satellite from the reference ellipsoid, R the observed range, 𝛥𝑔 the geophysical corrections, and SSH is corrected sea surface height. As can be seen in computational analysis, the mean computed biases may be greater than 35 cm. Table (5) shows the correlation coefficient between TG and SA time series, for some re-tracking methods. It indicates that this correction can improve SA observations in coastal regions. 2.3 Selection of suitable SA points Due to the well-known problems of SA technique near the coastal regions and shallow waters, choosing a suitable SA point to compare with TG station is of crucial importance. For this purpose, two strategies are considered in this study. For the first strategy, we select those SA ground tracks that pass very close to the TG station provided that the range observation is available on them (see Table 3). It should be mentioned that the closer the satellite ground tracks get to the TG stations, the more accurate would be the results. In 2002, Nerem and Mitchum [2002] showed that as the satellite tracks get closer to the coastal TG station, their time series get more similar. Furthermore, increasing the number of cycles will result in the reduction of total error of modeling. Thus if there are more than one pass near the TG station, we consider all of them. For example, according to Table 3, for Anzali TG, three points are suitable points for comparison. This enables more reliable estimation of VCM. Choosing the pass of SA close to the TG station as it is done in the first approach, means that the selected SA points are still some distances apart from TG stations. Moreover it means that no rigorous criteria is applied to determine the true distance. How much the two passes can be distance apart? Moreover, there is no SA pass that exactly crosses the location of TG station. To deal with this issue, an extrapolation method is used to predict the value of SA observation at the location of TG station in order to further improve the accuracy of VCM. In this approach three point-wise SA sea level time series are created. Then, the Spatio-Temporal ordinary Kriging is used for predicting point-wise sea level time series at TG stations. Ordinary kriging is a geo-statistical interpolation method, which was originally developed for spatial data by Krige (1952). In this study, SA measurements along the Caspian Sea, which are scattered in space and time, are used to predict the sea level at any time and location along the Caspian Sea. Kriging creates a statistically unbiased estimator that is optimal with a respect to the mean squared prediction error. The predictor is a weighted average of the observed values, in which the weights depend on spatial or Spatio-Temporal dependence of observed and predicted locations. The dependency is expressed by the covariance C ((x1; t1); (x2; t2)) between the sea levels at two space-time points (x1; t1) and (x2; t2). The required covariance models reflect that the sea level, are estimated based on empirical covariance values between SA observations at 𝛥𝑟 = 𝑐 2 × 𝜏 × (𝐺 − 𝐺0) 𝑆𝑆𝐻= ℎ − (𝑅 + 𝛥𝑟 + 𝛥𝑔) 5 Vertical crustal motion at Caspian Sea various locations. It is obvious that the construction of an appropriate covariance function C ((x1; t1); (x2; t2)) to predict the sea level is of particular importance for SA data along the Caspian Sea. For more details about Spatio- Temporal kriging and empirical covariance function see Boergens, et al. [2016]. Figure 2 shows a typical example of an empirical and modeled covariance functions for cycle 10 of T/P satellite pass 16. The modeled covariance function is used for prediction. Mahmood Pirooznia et al. 6 Considerations Duration Position Name of the point Distance Km A point close to Anzali TG from pass 16 1992-2014 Lat:38.31° Lon:50.18° Anzali-16 112.7 A point close to Anzali TG from pass 209 1992-2014 Lat:38.70° Lon:49.26° Anzali-209 137.8 A point close to Anzali TG from cross over passes 16 and 209 1992-2014 Lat:39.20° Lon:49.59° Anzali-Crossover 192.7 A point close to Noshahr TG from pass 16 1992-2014 Lat:36.95° Lon:51.12° Noshahr-16 48.7 A point close to Noshahr TG from cross over passes 16 and 31 1992-2014 Lat:37.09° Lon:51.02° Noshahr-Crosso- ver 64.9 A point close to Neka TG from passes 107 1992-2014 Lat:37.01° Lon:53.8° Neka-107 42.5 A point close to Neka TG from cross over passes of 107and 92 1992-2014 Lat:37.09° Lon:53.85° Neka-Crossover 50.6 Table 3. Specification of closest points to TG stations derived from satellite ground track by the first strategy. Figure 2. Empirical and fitted (exponential) temporal (a) and spatial covariance functions (b). a) b) In table 4 the extrapolated SA points near each tide gauge station is shown. Figure 3 shows the position of TG stations, SA points, extrapolated SA points and available GPS stations in the southern coast of Caspian Sea. The GPS stations will be later used to determine the absolute vertical deformation at this region and in order to make a comparison with the deformation obtained by combination of SA and TG observations. 2.4 Tidal analysis Finally, the tidal effects should be removed from both data sources. Since, global tide models are limited in the number of tidal components (frequencies), and for comparison of these two data sources, it is necessary that remove the similar tidal components from both observational sources; instead of using global tide model, local tide model 7 Vertical crustal motion at Caspian Sea Table 4. Extrapolated SA points at the location of TG stations using Spatio-Temporal Kriging. Considerations Duration Position Name of the points A point close to Anzali TG 1992-2014 Lat:37.48 Lon:49.45 Anzali-01 A point close to Noshahr TG 1992-2014 Lat:36.63 Lon:51.50 Noshahr-01 A point close to Neka TG 1992-2014 Lat:36.84 Lon:53.35 Neka-01 Figure 3. Distribution of SA, TG, SA Kriging and GPS stations (Kriging stations are extrapolated from observed3 stations on the ground track of satellite altimetry). (The apparent distance between extrapolated SA point and TG 4 station is for better graphical representation. They are actually at the same location). should be conducted in order to calculate the tidal effects precisely. In this case, tide modeling is done by the following equation [Hou and Vaníček, 1994; Okwuashi and Ndehedehe, 2017]: ℎ(𝜙,𝜆,𝑡) = 𝑍₀(𝜑,𝜆) + 𝑎 × 𝑡 + ∑�� =₁ 𝑐�cos𝜔�𝑡 + 𝑠�sin𝜔�𝑡 (2) in which 𝑍₀(𝜑,𝜆) stands for Mean Sea Level of Caspian Sea, 𝜔� = 2𝜋𝑓� , (𝑖 = 1,2,..., 𝑛) is the angular frequency of tidal components, 𝑡 and ℎ(𝜙,𝜆,𝑡) denote the time of observation and sea level height respectively. We have two kinds of tide model, one for SA observation and the other for TG observation each done for all available time series of SA and TG. For the first model, ℎ(𝜙,𝜆,𝑡) stands for the sea level derived for SA measurement while for the second model, ℎ(𝜙,𝜆,𝑡) stands for the sea level derived for TG measurement. 𝑎 is the rate of sea level height and the parameters 𝑐� and 𝑠� denote Fourier coefficients which are the unknowns of the problem and must be estimated through least square optimization. After implementing this analysis, the phases Φ� and amplitudes 𝐴� of essential tidal components can be derived as (Pirooznia et al, 2016b). 𝐴� = �(𝑐²� + 𝑠²� ) 𝛷� = tan‒1(𝑐�⁄𝑠�) (3) To determine the optimum number of tidal frequencies, a spectral analysis is conducted which reveals that tide gauge observations in our case include 40 significant frequencies. It should be mentioned that the analysis of SA time series, suffers from aliasing phenomena. According to Shannon sampling theorem (William, 1977; Shannon, 1984), if 𝜐 is the highest present frequency in a signal, a sampling rate of at least 𝑓� = 2𝜐 (sample-second) would be adequate to reconstruct the whole signal in which, 𝑓� is the well-known Nyquist frequency. In other words, a signal can be fully reconstructed or recovered by sampling it at series of points which are 1/2𝜐 seconds apart. If there is any frequency content in the signal that is higher than this sampling frequency, it will be aliased to the lower frequency. In the case of sea tide, periods of SA observations (10 days) are longer than period of tidal components and thus the determination of tidal effects using this data is subjected to aliasing problem. To deal with this, aliasing period is defined by the following equation [Lindsley and Long, 2011]. (4) in which, 𝑓 is the frequencies of tidal constituent, 𝛥𝑡 is sampling period in days (orbit repeat period for satellite Jason2 is 𝛥𝑡 = 9.9156 𝑑𝑎𝑦𝑠) and [•] is used for a Fix function1 which rounds decimal numbers toward zero [Lindsley and Long, 2011]. It should be emphasized that the sea level time series from TG observations may be contaminated with outliers and it may suffer from offset and gap in the time series [see Pirooznia et al., 2019]. These problems have been resolved with a proper strategies and then the time series have been constructed. The main problem with Caspian Sea tide gauge stations is the improper collecting and compiling of them. In this study, the quality control of data is done according to steps of quality control which previously discussed in Pirooznia et al. [2019]. A quality control of sea level observations is important task in a detecting and removing the possible errors, such as outlier, offset that may contaminate the tide gauge data time series [Woodworth et al. 2012; IOC 1992]. These errors can greatly affect the results of time series analysis (see http://www.sonel.org/_Quality-control-of-the.html). The existence of outliers in the time series of observations can introduce some biases in the estimation of model parameters. Outliers may increase the confidence intervals for the model parameters, and thus, they may strongly influence modeling and predictions [Chen and Liu 1993a, 1993b]. Offset in the time series of tide gauge may occur during data acquisition process, due to some phenomena such as, boat hitting, tide gauge damage, earthquake motions and wrong recalibration of the instrument. Occurring of offsets in the observations can change the statistical parameter of signal and may affect on the trend computation [Gazeaux et al. 2013; Montillet et al. 2015; Lopez Bernal et al. 2016]. The first step in a dealing with these types of errors is modeling sea level time series. In this study, the functional model of the tide (Eq.) is used 1 The difference between Int and Fix is in negative numbers; for a number like -9.2 Int returns -10, whereas Fix returns -9. 2𝜋𝛥𝑡 2𝜋(𝑓𝛥𝑡-[𝑓𝛥𝑡 + 0.5]) 𝑇𝑎 = Mahmood Pirooznia et al. 8 http://www.sonel.org/_Quality-control-of-the.html to model sea level. Before that modeling, the low pass filtering (Savitzky–Golay filter) is used to reduce the noises in sea level time series [Shekhar, 2016]. Then, the time series is constructed (modeled) using significant frequencies and the 3-sigma threshold method is used to eliminate the possible outliers in the observations. Subsequently, for offset detection, noise assessment algorithm (least square variance component estimation) is used to find the effective noise component in the time series, and using the cofactor matrix of the noise as the weight matrix of observations, the observations are adjusted and the signal is processed to find the possible times of offsets in order to remove their corresponding observations from the time series. For more details see Pirooznia, et al. [2019]. For modeling sea level time series, the significant frequencies should be found. For this purpose, spectral analysis is performed on the sea level time series. To this end, angular frequencies are created using the period of observations (𝜔 = 2𝜋/𝑇). Assume that �𝜔₁,𝜔₂,...,𝜔�� is the set of real and positive angular frequencies of the observation vector that may be equal to tidal angular frequencies, or they can be constructed by period of observations. In this study, in order to find effective frequencies, the angular frequencies are determined based on the period of observations. To check the effective frequencies, the significant test is employed [Amiri-Simkooei and Asgari, 2011] and the effective frequencies which pass this test are selected. The result of spectral analysis on SA and TG of Caspian Sea shows that the annual and semiannual components have the highest effect on time series (see Figures 4 and 5). As mentioned previously, spectral analysis of SA time series due to sampling rate of observations faces with the aliasing phenomenon which creates the aliasing frequency. As an example, for the component of half daily S2, the main frequency is 2 cycles per day, whose aliased frequency is 0.01702 cycles per day corresponds to period of 58.74 days (see Figure 4). Figure 4 shows a typical power spectrum for SA time series of Anzali-crossover point, alias frequency of S2 component, annual (Sa) and semi-annual (Ssa) frequencies. The results of tide analysis indicate that the annual (Sa) and semiannual solar (Ssa) components on all of the SA and TG stations have the highest amplitude in comparison with the other components. The amplitude of annual components in three TG stations (Azali, Noshahr and Neka) are 16, 18 and 15 cm respectively, and for semiannual components are 2.8, 5.4 and 3.7 cm, respectively (Figure 5). 9 Vertical crustal motion at Caspian Sea Figure 5. Tide analysis of three tide gauge stations (Anzali, Noshahr and Neka). Figure 4. Power spectrum of Anzali-crossover point SA time series, alias frequency of S2 component, annual (Sa) and semi-annual (Ssa) frequencies. The tables which include amplitudes and phases of tidal components (both SA and TG data stations) are added in the appendix A. Figures 6-8 show the observed time series (SA and TG), tidal effect (SA tide and TG tide) and the corrected time series (show by residual signal in the Figures) in Anzali, Noshahr and Neka stations. In each Figure we show, (i) time series of TG observation (ii) time series of SA observations at selected SA points using two adopted approaches (Referring to Table 3, 4 and Figure 3 for more explanations). (iii) Corrected time series after removal of tide model. The composed tidal model appropriately estimates the general variations of sea level height time series. Judging from SA and TG time series ordinates (see Figures 6-8), it is evident that the Caspian Sea level has been decreasing during the last decades. According to Arpe et al. [2012], the Volga River as the main source of Caspian Sea water supply, has brought less water to the sea. Mahmood Pirooznia et al. 10 Figure 7. Left: Observed Sea level time series (purple line), tide model (blue line), and the Mean Sea Level (Red dashed line) at points Noshahr-16 SA, Noshahr-Crossover SA, Noshahr-01 SA and Noshahr TG station. Right: Residual (green line) at points Noshahr-16 SA, Noshahr-Crossover SA, Noshahr-01 SA and Noshahr TG station. Figure 6. Left: Observed Sea level time series (purple line), tide model (blue line), and the Mean Sea Level (Red dashed line) at points Anzali-16 SA, Anzali-209 SA, Anzali-Crossover SA, Anzali-01 SA and Anzali TG station. Right: Residual (green line) at points Anzali-16 SA, Anzali-209 SA, Anzali-Crossover SA, Anzali-01 SA and Anzali TG station. 3. Numerical results and discussion of VCM In this section, we present the numerical results of VCM at the coastal area of Caspian Sea. To do so, the SA and TG observations that are processed based on the steps of previous section are subtracted from each other. As mentioned previously, to improve the accuracy of SA time series, the re-tracking correction is implemented to range observation. Table 5 shows the correlation coefficient between the sea level that is resulted by Anzali, Noshahr and Neka tide gauge observations and SA observations (Anzali-16 point, Noshahr-16 point and Neka-107 point), without applying any re-tracking and while the re-tracking is performed. In this region the threshold re-tracking method, by 35% as threshold value, presents the highest correlation. That means, by using this re-tracking method the most accurate measurement of sea level will be achieved. It should be emphasized that the vertical datum of SA and TG observations are different; thus, to compute the correlation between SA and TG time series, after normalizing the time series of observations, the correlations are computed at the times when SA and TG observations overlap. Also due to the well-known problems of SA technique near the coastal regions and shallow waters, choosing a suitable SA point to compare with TG station is of crucial importance. For this purpose, two strategies are considered in this study. 11 Vertical crustal motion at Caspian Sea Figure 8. Left: Observed Sea level time series (purple line), tide model (blue line), and the Mean Sea Level (Red dashed line) at points Neka-107 SA, Neka-Crossover SA, Neka-01 SA and Neka TG station. Right: Residual (green line) at points Neka-107 SA, Neka-Crossover SA, Neka-01 SA and Neka TG station. Table 5. Correlation coefficient between TG and SA time series, for some re-tracking methods. Method Correlation coefficient Between Anzali TG and Anzali-16 SA point Correlation coefficient Between Noshahr TG and Noshahr-16 SA point Correlation coefficient Between Neka TG and Neka-107 SA point Not re-tracked 0.55 0.46 0.53 OCOG 0.77 0.68 0.62 �-Parameter 0.83 0.74 0.69 Threshold (15%) 0.72 0.68 0.73 Threshold (25%) 0.82 0.76 0.76 Threshold (35%) 0.85 0.82 0.81 Threshold (40%) 0.79 0.72 0.77 Threshold (45%) 0.78 0.70 0.71 Threshold (50%) 0.73 0.64 0.75 SA and TG respectively indicate the absolute and relative changes in sea levels. Difference of the abovementioned two rates gives the VCMs. Table. 6 lists the trends of both SA and TGs, their uncertainties and differences. The average vertical crustal motion estimated from SA-TG at Anzali, Noshahr and Neka are -1.387, -1.903 and -3.14 mm/year, respectively. Figure 9 shows Mean value of VCM at TG stations. Analyzing GPS time series at Geodynamic stations is an alternative approach to obtain VCM’s. This analysis supplies the auxiliary source to be compared with the method of this study. In this study, four geodynamic stations are used, namely, MABD, NKAD, RSHT and TKBN. Table 7 shows the approximate positions of GPS stations, time span of their observations and their distances to closest TG stations. These GPS time series were provided by national cartography center (NCC) of Iran. All processing steps were done in NCC. The processing of GPS raw data was done by GAMIT software, also it was integrated using GLOBK software [Herring et al., 2015]. According to recent report in NCC (www.ncc.org.ir), the first step in processing this data, is controlling the quality of rinex data. This operation is done with TEQC software. Then data processing on a daily basis is done using double differences method. Finally, the precise orbital data (SP3 files from International GNSS service (IGS)) with a delay of two weeks, mapping function equations for troposphere modeling, ocean and solid tide modeling and solar radiation model are applied. Moreover, the IGS stations are referenced with a respect to ITRF 2014 reference frame to minimize the calculation errors and increase the accuracy. Figure 10 shows vertical GPS time series at geodynamic stations. Mahmood Pirooznia et al. 12 Figure 9. Mean VCM estimated from SA-TG. Table 6. The rate of sea level changes via SA and TG’s and VCM estimation using SA-ALT. TG stations Rate (mm/year) SA points Rate (mm/year) VCM=SA-TG (mm/year) Anzali -26.40±0.65 Anzali-16 -27.78±1.12 -1.38±1.29 Anzali-209 -27.24±0.73 -0.84±0.97 Anzali-Crossover -27.96± 0.91 -1.56±1.11 Anzali-01 -28.17± 0.83 -1.77±1.05 Noshahr -30.22± 0.81 Noshahr-16 -32.85± 0.55 -2.63±0.97 Noshahr-Crossover -31.50± 0.76 -1.28±1.11 Noshahr-01 -32.02± 0.63 -1.8±1.02 Neka -21.10±0.48 Neka-107 -23.41± 1.03 -2.31±1.13 Neka-Cross over -26.11± 0.63 -5.01±0.79 Neka-01 -23.2± 0.85 -2.1±0.97 http://www.ncc.org.ir In order to compare VCMs from GPS with those from combination of SA and TG, it is necessary that GPS stations are located near to TG stations. Even though these GPS stations are not located sufficiently close to our TG stations (see table 7), both of them are positioned on the same faults (see Figure 11). 13 Vertical crustal motion at Caspian Sea Table 7. Positions of Geodynamic GPS stations and observation time span. Name of station Approximate position Duration Distance to the nearest TG stations MABD Lat:36.58 Lon:52.28 6/8/2005-18/7/2013 101.322 km to Neka TG station NKAD Lat:36.68 Lon:51.3 18/10/2005-08/9/2012 18.148 km to Noshahr TG station RSHT Lat:37.31 Lon:49.61 1/11/2005-18/7/2013 22.184 km to Anzali TG station TKBN Lat:36.78 Lon:50.91 16/8/2005-08/9/2012 54.540 km to Noshahr TG station Figure 10. Vertical GPS time series at MABD, NKAD, RSHT and TKBN stations. Figure 11. The location of GPS and TG stations in the dispersion of faults. The listed locations for TG and geodynamic stations shows that there is a notable distance between them. Due to long distances (few Kilometers) between geodynamic (GPS) and TG stations, comparing their results may look unsound. However, from tectonic point of view, these stations are situated on a common fault and their comparison would be meaningful. The resulted VCM’s from previous section are compared with the crustal motions deduced from GPS time series at geodynamic stations. The comparison shows an acceptable consistency (Table.8). Table 8 indicates that there is a good consistency between (SA-TG) and GPS vertical rates. The mean value of the trend differences between the closest GPS and the SA-TG vertical rates is 0.82, 0.21 and 2.71 mm/year for Anzali, Noshahr and Neka TG stations, respectively. At all stations, the rates of VCM are negative which indicates the subsidence phenomenon. Various factors may be responsible to this subsidence: tectonic changes due to existing faults, irregular groundwater extraction and etc. [Pacheco-Martínez et al., 2015; Carreón-Freyre et al., 2016]. According to the studies conducted by the Iranian ministry of energy (http://www.moe.gov.ir) using the data of piezometric wells, the rate of groundwater extraction has significantly increased in the recent years. This conclusion can be used as a confirmation to the negative rates of VCM in these areas. For instance, the highest rate of vertical change at MABD station might have originated from the extravagant amount of groundwater extraction which caused the crust to subside. For more investigation, the correlations between ground water level data of three wells near Anzali, Noshahr and Neka TG stations with GPS vertical time series of RSHT, NKAD and MABD stations are computed. These piezometric wells data were recorded from 23-Jul-2002 until 15-Dec-2014, from 28-Jan-2003 until 21-Dec-2013 and from 16-Aug-2004 until 29-Dec-2013 with sampling interval of approximately one month near Anzali, Noshahr and Neka TG stations, respectively. These wells are located at 37.36° north latitude and 49.86° east longitude, 36.63° north latitude and 51.50° east longitude, 36.72° north latitude and 53.26° east longitude, respectively. To compute the correlations between piezometric wells data and GPS time series, after normalizing the time series, the correlations are computed at the times when piezometric wells data and GPS observations overlap. Table 9 shows the computed correlations. Mahmood Pirooznia et al. 14 Table 8. The rate of VCM via SA and TG’s along with GPS vertical rate at geodynamic stations. Station SA - TG GPS Points Time span Rate (mm/year) Station Time span Rate (mm/year) Anzali Anzali-16 2005-2013 -0.82±0.72 RSHT 2005-2013 -0.05±0.46 Anzali-209 2005-2013 -1.06±0.84 Anzali-Crossover 2005-2013 -0.71±0.68 Anzali-01 2005-2013 -0.91±0.76 Noshahr Noshahr-16 2005-2013 -2.32±1.09 TKBN 2005-2013 -2.10±0.60 Noshahr-Crossover 2005-2013 -2.16±0.76 NKAD 2005-2013 -1.51±0.43 Noshahr-01 2005-2013 -2.45±1.35 Neka Neka-107 2005-2012 -5.12±1.07 MABD 2005-2012 -7.92±0.61Neka-Cross over 2005-2012 -5.64±0.86 Neka-01 2005-2012 -4.85±1.01 http://www.moe.gov.ir/ Figure 12 shows the groundwater level time series near Anzali, Noshahr and Neka TG stations. The rates of groundwater level time series are negative which may be the sign of groundwater extraction 4. Conclusion and recommendation In this paper, the pattern of vertical deformation in coastal parts of Caspian Sea is determined using a combination of satellite altimetry observations and tide gauge measurements. Satellite altimetry measures the absolute sea level, while tide gauge data computes the relative sea level. Thus the difference between these two observations can be used to determine the vertical motion at ground station. The result show that the sea level of Caspian Sea from 1992 to 2014 has had a decreasing trend. Moreover, the VCMs at the desired stations are negative that may indicate land subsidence. Typically, the numerical results illustrate that by considering SA crossover points, the vertical crustal motion at Anzali, Noshahr and Neka are -1.56, -1.28 and -5.01 mm/year, respectively. The innovations of this study are as follows: using re-tracking correction for SA data neat coastal TG stations that increases the correlation between SA and TG sea level time series, application of Spatio-Temporal kriging method by creating the empirical covariance functions to interpolate SA signal at the TG stations and implementing local and more accurate tide modeling instead of using global tide model such as FES2004. The future issues that may significantly improve the quality of the results are reported below The precision of orbit determination for satellite altimetry is 1-2cm. This leads to a precision of about 4-5cm in sea surface height. An increase in precision of orbit determination will result in a more accurate sea surface height data. 15 Vertical crustal motion at Caspian Sea Ground water level data near TG stations GPS stations Correlation coefficient Anzali RSHT 0.62 Noshahr NKAD 0.47 Neka MABD 0.43 Table 9. The correlation coefficient between ground water level data of three wells near TG stations and GPS. Figure 12. The underground water level time series near Anzali, Noshahr and Neka TG stations. TG stations should gather data in a unified vertical datum. This unification is of fundamental importance and simplifies the calculation process. In some researches the Mean Sea Level (MSL) which is determined from SA and TG data is used as an alternative for estimation of VCM, in this case there is no need to calculate the rate of change and VCM is calculated by subtraction of MSL from SA and TG data. This needs all TG stations to be geo-referenced. TG stations should be calibrated precisely. Furthermore, TG observations at nearby stations, should be given appropriate weightings and combined with each other in order to increase the accuracy. To integrate the elevations of different TG’s, stations must be fixed to a local benchmark; this benchmark then must be transferred to a national and finally an international datum. This procedure requires a modern geodetic technique. Various tectonic factors involved in vertical deformation and land subsidence due to groundwater extraction may be considered as causes of VCM at coastal regions of Caspian Sea. One of the drawback of this study is the lack of TG data for a long time. With long time data, it is possible to measure trends more accurately, and provide robust estimation of the trend. References Amiri-Simkooei, A. R. and J. Asgari. (2011). Harmonic Analysis of Total Electron Contents Time Series: Methodology and Results, Journal of GPS Solution, 16 (1): 77–88. doi:10.1007/s10291-011-0208x Arpe, K., L. Bengtsson, G. S. Golitsyn, I. Mokhov, V. A. Semenov and P. V. Sporyshev, (2000). Connection between Caspian Sea level variability and ENSO, Geophys. Res. Lett., 27(17), 2693-2696. doi:10.1029/1999gl002374. Arpe, K., S. A. Leroy, H. Lahijani, and V. Khan (2012). Impact of the European Russia drought in 2010 on the Caspian Sea level, Hydrol. Earth Sys. Sci., 16(1), 19-27. doi:10.5194/hess-16-19-2012 Beavan, J., P. Denys, M. Denham, B. Hager, T. Herring and P. Molnar. (2010). “Distribution of Present- Day Vertical Deformation across the Southern Alps, New Zealand, from 10 Years of GPS Data”, Geophys. Res. Lett., 37 (16): n/a-n/a. doi:10.1029/2010GL044165. Benada, J. R. (1997). “PO.DAAC Merged GDR (TOPEX/POSEIDON) Generation B User’s Handbook”, Version 2.0, see the JPL PO.DAAC Web page: http://podaac.jpl.nasa.gov, where this User’s Handbook can be viewed on-line., Jet Propulsion Laboratory, California Institute of Technology under Contract with the National Aeronautics and Space Administration. Blewitt, G. (2007). “GPS and Space-Based Geodetic Methods”. In Treatise on Geophysics, edited by Tom Herring, 3, 351–390. Cambridge, MA: Massachusetts Institute of Technology, doi:10.1016/B978-044452748-6.00058-4. Boergens, E., S. Buhl, D. Dettmering, C. Klüppelberg, and F. Seitz (2016). Combination of multi-mission altimetry data along the Mekong River with spatio-temporal krigin, J. Geod., 91(5), 519-534. doi:10.1007/s00190-016- 0980-z. Braitenberg, C., P. Mariani, L. Tunini, B. Grillo, B. and I. Nagy (2011). Vertical crustal motions from differential tide gauge observations and satellite altimetry in southern Italy,  J. Geodyn., 51(4), 233-244. doi:10.1016/j.jog.2010.09.003 Brunet, M., M. V. Korotaev, A. V. Ershov and A. M. Nikishin (2003). The South Caspian Basin: a review of its evolution from subsidence modelling, Sedimentary Geology, 156(1-4), 119-148. doi:10.1016/s0037-0738(02)00285-3  Caccamise, D. J., M. A. Merrifield, M. Bevis, J. Foster, Y. L. Firing, M. S. Schenewerk, F. W. Taylor and D. A. Thomas (2005), Sea level rise at Honolulu and Hilo, Hawaii: GPS estimates of differential land motion, Geophys. Res. Lett., 32, L03607, doi:10.1029/2004GL021380 Cazenave, A., K. Dominh, F. Ponchaut, L. Soudarin, J. F. Cretaux and C. Le Provost, (1999). Sea level changes from Topex-Poseidon altimetry and tide gauges, and vertical crustal motions from DORIS,  Geophys. Res. Lett., 26(14), 2077-2080. doi:10.1029/1999gl900472. Carreón-Freyre, D., M. Cerca, G. Ochoa-González, P. Teatini and F.R. Zuñiga (2016). Shearing along faults and stratigraphic joints controlled by land subsidence in the Valley of Queretaro, Mexico.  Hydrogeology Journal, 24(3), 657-674. doi:10.1007/s10040-016-1384-0. Chen, C. and L. Liu (1993a). Forecasting Time Series with Outliers.” Journal of Forecasting. 12 (1):13–35. doi:10.1002/for.3980120103. Chen, C. and L. Liu. (1993b). Joint Estimation of Model Parameters and Outlier Effects in Time Series.” J. Am. Stati. Ass., 88 (421):284–297, doi:10.1080/01621459.1993.10594321. Mahmood Pirooznia et al. 16 https://doi.org/10.1029/2004GL021380 Clark, J. A., P. E. Haidle and L. N. Cunningham (2002). Comparison of Satellite Altimetry to Tide Gauge Measurement of Sea Level: Predictions of Glacio-Isostatic Adjustment, J. Climate, 15(22), 3291-3300. doi:10.1175/1520- 0442(2002)015<3291:cosatt>2.0.co;2.  Dalla Via, G., M. Crosetto and B. Crippa (2012). Resolving Vertical and East-West Horizontal Motion from Differential Interferometric Synthetic Aperture Radar: The L’Aquila Earthquake. J. Geophys. Res.: Solid Earth 117 (2): 1– 14. doi:10.1029/2011JB008689. Davis, C. H. (1995). “Growth of the Greenland Ice Sheet: A Performance Assessment of Altimeter Retracking Algorithms”. IEEE Trans. GeoSci. Remote Sens. 33 (5): 1108–1116. doi:10.1109/36.469474. Davis C. H. (1997). “A robust threshold retracking algorithm for measuring ice-sheet surface elevation change from satellite radar altimeters”. IEEE Trans. Geosci. Remote Sens., vol. 35, no. 4, pp. 974–979. Deng X. (2003). “Improvement of geodetic parameter estimation in coastal regions from satellite radar altimetry”. Ph.D. dissertation, Curtin Univ. Technol., Perth, WA, Australia. Fu, L. L. and A. Cazenave. (2001). Satellite Altimetry and Earth Sciences: A Handbook of Techniques and Applications, 463. San Diego, CA: Academic Press. García, D., I. Vigo, B. F. Chao and M. C. Martínez (2007). Vertical Crustal Motion along the Mediterranean and Black Sea Coast Derived from Ocean Altimetry and Tide Gauge Data, Pageoph Topical Volumes, 851-863. doi:10.1007/978-3-7643-8417-3_13. Gazeaux, J., S. Williams, M. King, M. Bos, R. Dach, M. Deo and F. H. Webb. (2013). “Detecting Offsets in GPS Time Series: First Results from the Detection of Offsets in GPS Experiment”. J. Geophys. Res.: Solid Earth 118 (5):2397–2407. doi:10.1002/jgrb.50152. Gomez-Enri, J., S. Vignudelli, G. Quartly, C. Gommenginger, and J. Benveniste. (2009). “Bringing Satellite Radar Altimetry Closer to Shore”. Remote Sensing, SPIE Newsroom. doi:10.1117/2.1200908.1797. Herring, T. A., R. W. King, M. A. Floyd, S. C. McClusky (2015). Introduction to GAMIT/GLOBK, Release 10.6, Mass Inst. of Technol., Cambridge. Hou T. and P. Vaníček, (1994). Towards a real-time analysis of tides. International Hydrographic Review, LXXI (1), Monaco, pp.29-52. Intergovernmental Oceanographic Commission (IOC). 1992 Joint IAPSO-IOC workshop on sea-level measurements and quality control, Paris, 12–13 October 1992, Intergovernmental Oceanographic Commission, Workshop Report, No. 81. Jackson, J., K. Priestley, M. Allen and M. Berberian (2002). Active tectonics of the South Caspian Basin. Geophys. J. Int., 148(2), 214-245, doi:10.1046/j.1365-246x.2002.01005.x. Johansson, J. M., J. L. Davis, H. G. Scherneck, G. A. Milne, M. Vermeer, J. X. Mitrovica, R. A. Bennett, B. Jonsson, G. Elgered, P. Elósegui, H. Koivula, M. Poutanen, B. O. Rönnäng and I. I. Shapiro I. I., (2002). Continuous GPS measurements of postglacial adjustment in Fennoscandia 1. Geodetic results. J. Geophys. Res., 107(B8), doi:10.1029/2001jb000400 Khaki, M., E. Forootan and M. Sharifi (2014). Satellite radar altimetry waveform retracking over the Caspian Sea. Int. J. Remote Sens., 35(17), 6329-6356, doi:10.1080/01431161.2014.951741. Khajeh, S., Jazireeyan, I., & Ardalan, A. A. (2014). Applying Satellite Altimetry to Wetland Water Levels Monitoring (Case Study: Louisiana Wetland). IEEE Geosci. Remote Sens. Lett., 11(9), 1475-1478, doi:10.1109/lgrs.2013.2295831. Koch, K. (1999). Parameter Estimation in Linear Models, Parameter Estimation and Hypothesis Testing in Linear Models, 149-269, doi:10.1007/978-3-662-03976-2_4. Kuo, C., C. K. Shum, A. Braun, K. Cheng and Y. Yi (2008). Vertical Motion Determined Using Satellite Altimetry and Tide Gauges. Terrestrial, Atmospheric and Oceanic Sciences, 19(1-2), 21, doi:10.3319/tao.2008.19.1-2.21(sa). Kuo, C. Y. (2004). Vertical crustal motion determined by satellite altimetry and tide gauge data in Fennoscandia, Geophys. Res. Lett., 31(1), doi:10.1029/2003gl019106. Krige, D. G., (1952). A statistical analysis of some of the borehole values in the Orange Free State goldfield, J. Chemical, Metallurgical and Mining Society of South Africa, 47-70. Larson, K. M. and V. D. Tonie (2000). Measuring Postglacial Rebound with GPS and Absolute Gravity, Geophys. Res. Lett. 27 (23): 3925–3928, doi:10.1029/2000GL011946. Lee H. K. (2008). “Radar altimetry methods for solid earth geodynamics studies,” Ohio State Univ., Columbus, OH, USA, Rep. No. 489. 17 Vertical crustal motion at Caspian Sea Lin, J., (2000): Correction of tide gauge measurements to absolute sea level by vertical motion solutions. Master Thesis, the Ohio State University, Columbus, Ohio, USA, 66. Lindsley, R. D. and D. Long (2011). Fitting tidal constituents to altimeter data, UTAH SPACE GRANT CONSORTIUM 2011, Utah State University, Logan. Lopez Bernal, J., S. Cummins and A. Gasparrini (2016). “Interrupted Time Series Regression for the Evaluation of Public Health Interventions: A Tutorial”. International Journal of Epidemiology, doi:10.1093/ije/dyw098. Fenoglio-Marc, L., C. Dietz, and E. Groten (2004). Vertical Land Motion in the Mediterranean Sea from Altimetry and Tide Gauge Stations. Marine Geod., 27(3-4), 683-701, doi:10.1080/01490410490883441. Mammadov, P. (2008). The Subsidence Evolution of the South Caspian Basin. Caspian and Black Sea Geosciences Conference, doi:10.3997/2214-4609.20146089. Mathers L, P.A.M. Berry, J. A. Freeman (2004). Gravity, geoid and space missions. In: International association of geodesy symposia, Porto, August 30–September 3. Martin, T. V., H. Zwally, A. C. Brenner and R. A. Bindschadler (1983). Analysis and Retracking of Continental Ice Sheet Radar Altimeter Waveforms, J. Geophys. Res., 88 (C3): 1608, doi:10.1029/JC088iC03p01608. Montillet, J., S. D. Williams, A. Koulali and S. C. McClusky (2015). Estimation of Offsets in GPS Time Series and Application to the Detection of Earthquake Deformation in the Far-Field, Geophys. J. Int., 200 (2):1207–1221, doi:10.1093/gji/ggu473. Nerem, R. S. and G.T. Mitchum (2002). Estimates of vertical crustal motion derived from differences of TOPEX/POSEIDON and tide gauge sea level measurements,  Geophys. Res. Lett.,  29(19), 40-1-40-4, doi:10.1029/2002gl015037. Noomen, R., T. A. Springer, B. A. C. Ambrosius, K. Herzberger, D. C. Kuijper, G. J. Mets, B. Overgaauw and K. F. Wakker (1996). Crustal Deformations in the Mediterranean Area Computed from SLR and GPS Observations, J. Geodyn., 21 (1): 73–96, doi:10.1016/0264-3707(95)00015-1. Okwuashi, O. and C. Ndehedehe (2017). Tide modelling using support vector machine regression, J. Spatial Sci., 62:1–18, doi:10.1080/14498596.2016.1215272 Ostanciaux, É., L. Husson, G. Choblet, C. Robin and K. Pedoja (2012). Present-day trends of vertical ground motion along the coast lines, Earth-Sci. Rev., 110(1-4), 74-92, doi:10.1016/j.earscirev.2011.10.004. Pacheco-Martínez, J., E. Cabral-Cano, S. Wdowinski, M. Hernández-Marín, J. Ortiz-Lozano and M. Zermeño-de- León (2015). Application of InSAR and Gravimetry for Land Subsidence Hazard Zoning in Aguascalientes, Mexico, Remote Sensing, 7(12), 17035-17050, doi:10.3390/rs71215868. Pagiatakis, S. D. and P. Salib (2003). Historical relative gravity observations and the time rate of change of gravity due to postglacial rebound and other tectonic movements in Canada. J, Geophys. Res.: Solid Earth, 108(B9)., doi:10.1029/2001jb001676. Pirooznia, M., S. R. Emadi and M. N. Alamdari (2016a). The Time Series Spectral Analysis of Satellite Altimetry and Coastal Tide Gauges and Tide Modeling in the Coast of Caspian Sea, Open J. Marine Sci., 06(02), 258-269, doi:10.4236/ojms.2016.62021. Pirooznia, M., S. Rouhollah Emadi and M. Najafi Alamdari (2016b). Caspian Sea Tidal Modelling Using Coastal Tide Gauge Data, J. Geol. Res., 1-10, doi:10.1155/2016/6416917.  Pirooznia, M., M. Raoofian Naeeni and Y. Amerian (2019). A Comparative Study Between Least Square and Total Least Square Methods for Time–Series Analysis and Quality Control of Sea Level Observations. Marine Geod., 1-25, doi:10.1080/01490419.2018.1553806. Rapley, C. G., M. A. J. Guzkowska, W. Cudlip and I. M. Mason. (1987). An Exploratory Study of Inland Water and Land Altimetry Using Seasat Data, ESA Rep. 6483/85/NL/BI. Neuilly: European Space Agency. Shannon, C., (1984). Communication in the presence of noise, Proc. IEEE 72, 1192– 1201, http://dx.doi.org/10.1109/PROC.1984.12998. Shekhar, C. (2016). On simplified application of multidimensional Savitzky-Golay filters and differentiators, doi:10.1063/1.4940262. Shum C.K., C.Y. Kuo and J.X. Mitrovica (2002). Glacial isostatic adjustment in the Great Lakes region inferred by tide gauges and satellite altimetry, EOS Trans., Spring AGU Meeting 2002. Sitar, K., (2011). Determination of sea level trends and vertical land motions from satellite altimetry and TG observations at the Mediterranean coast of Turkey, a thesis submitted to the graduate school of natural and applied sciences of Middle East Technical University, 45. Mahmood Pirooznia et al. 18 http://dx.doi.org/10.1109/PROC.1984.12998 Soltanpour, A., M. Pirooznia, S. Aminjafari and P. Zareian (2017). Persian Gulf and Oman sea tide modeling using satellite altimetry and tide gauge data (TM-IR01), Marine Geores. Geotechn., 1-11. doi:10.1080/1064119x.2017.1366608. Soudarin, L., J. Crétaux and A. Cazenave. (1999). Vertical Crustal Motions from the DORIS Space-Geodesy System, Geophys. Res. Lett. 26: 1207–1210, doi:10.1029/1999GL900215. Tesmer, V., P. Steigenberger, M. Rothacher, J. Boehm and B. Meisel (2009). Annual Deformation Signals from Homogeneously Reprocessed VLBI and GPS Height Time Series, J. Geod., 83 (10): 973–988, doi:10.1007/s00190-009-0316-3. Vaníček, P. (1969). Approximate spectral analysis by least-squares fit, Astrophys. Space Sci., 4(4), 387-391, doi:10.1007/bf00651344. Wang Yu (2004). Ocean Tide Modeling in Southern Ocean, Department of civil and Environmental Engineering and Geodetic Science, Report No 471, The Ohio state university Columbus William, G. C. (1977). Sampling techniques, 3rd ed. Canada: John Willey & Sons Inc. Woodworth, P., C. Hughes, R. Bingham and T. Gruber (2012). Towards Worldwide Height System Unification Using Ocean Information, J. Geodetic Sci. 2 (4), doi:10.2478/v10156-012-0004-8. Zhang M. (2009). Satellite radar altimetry for inland hydrologic studies, Ohio State Univ., Columbus, OH, USA, Rep. No. 491. *CO R R E S PO N D I N G A U T H O R: Mehdi RAO O F I A N N A E E N I, Department of Geodesy, Faculty of Geodesy and Geomatics Engineering, K. N. Toosi University of Technology, e-mail: mraoofian@kntu.ac.ir © 2020 the Istituto Nazionale di Geofisica e Vulcanologia. All rights reserved 19 Vertical crustal motion at Caspian Sea