The Impact of the Baltic Sea Non-tidal Loading on GNSS Station Coordinate Time Series: the Case of Latvia

The objective of this study is to discover the geodynamic processes of the Earth’s crust in the territory of Latvia occurred due to the effect of the Baltic Sea non-tidal loading, by way of using GNSS permanent station daily coordinate time series and tide gauge data to find correlations between two data sets for the period from 2012 up to 2018. For this study observations of 31 Latvian and 2 Estonian GNSS stations were used. Stations belong to the LatPos, EUPOS-Riga, EPN and EstPos networks. Station daily coordinate time series were computed using Bernese GNSS software v5.2 in a double-difference mode with 9 fiducial stations from International GNSS Service and EUREF Permanent GNSS Network. The analysis of obtained data significantly increases understanding of the Earth’s surface displacements occurring due to the loading effect in


Introduction
Non-tidal ocean loading is a geophysical factor resulting in the seasonal signals undermining GNSS (Global Navigation Satellite System) station coordinate time series. As GNSS observations have been essential in studies such as tectonic motion, sea-level rise and post-glacial rebound, correcting for the non-tidal ocean loading effects has thus been a critical prerequisite for estimating highly reliable secular rates of Earth's crust movements.
The detailed models and algorithms for eliminating the influence of global plate tectonics, solid Earth tides, polar tides, ocean and atmospheric tidal loading effects on station positions are given in the IERS (International Earth Rotation and Reference Systems Service) Conventions, and mostly ensured through the scientific postprocessing software packages. The correction for non-tidal ocean loading still needs further investigation (Geng and Williams, 2012).
The Baltic region is unique with the postglacial land uplift of Fennoscandia and the Baltic Sea itself. The Baltic Sea exhibits a number of remarkable phenomena. One of them is the sea level variations due to winds, complicated by the shape of gulfs and islands (Ruotsalainen et al., 2015;Lisitzin, 1959;Wübber and Krauss, 1979). In the case of tidal variations of the Baltic Sea the range is in the order of centimetres only, but under mentioned influence, the range of sea level variations can reach 3 m on the coasts of gulfs (Virtanen and Mäkinen, 2003). Such loading effect, occurred due to changing redistribution of the water masses, cannot be neglected especially for the coastal regions. For example, one-meter layer of water in the whole of the Baltic Sea will induce vertical deformation about -21 mm near Gotland Island. For the territory of west coast of Latvia, it will reach about -16 mm. But the horizontal maximums are about 4 mm . It should be taken into account, because the range of sea level variations at the Latvian coast can reach 2 m, what may cause larger deformations.
It has been confirmed that in GNSS time series 40% of the variation can be due to varying atmospheric, hydrological and ocean masses (Dong et al., 2002), which in turn makes the derived rate uncertainties larger than expected (Blewitt and Lavallée, 2002;Collilieux et al., 2012).
The main purpose of this study is to discover the influence of the Baltic Sea non-tidal loading on Latvian and two Estonian GNSS permanent station positions and to increase understanding of the geodynamic processes occurring due to this effect in the territory of Latvia.

GNSS station cumulative and annual horizontal velocities
For this study observations of 31 Latvian and 2 Estonian GNSS stations were used. Stations belong to the LatPos, EUPOS ® -Riga, International GNSS Service (IGS)/EUREF Permanent GNSS Network (EPN) and EstPos network. Station daily coordinate time series were computed for the period of 7 years: from 2012 up to 2018, in a double-difference (DD) mode using Bernese GNSS software v5.2 (Dach et al., 2015). Since 2015 the combined processing of GPS and GLONASS observations was applied (only GPS observations were used before). 9 EPN/IGS stations were used as fiducial stations: BOR1, GLSV, JOEN, LAMA, MAR6, ONSA, PULK, RIGA, VLNS. A cut-off elevation angle of 3° was selected. Station positions were corrected for the effect of solid Earth tides (Petit and Luzum, 2010) and the ocean tidal loading (FES2004 ocean tide model was used), and are without corrections of the atmospheric tidal loading. Coordinates obtained in IGb08 and IGS14 transformed to ETRF2000 using one-step transformation with 14 transformation parameters according to Boucher and Altamimi (2011). For the time series analysis Tsview software has been used with the standard settings. Written in Matlab, this software complements the GAMIT/GLOBK software package (Herring, 2003).
In Figure 1 the Latvian GNSS station cumulative horizontal velocities (in black) obtained for the period from 2012 up to 2017 (Varna et al, 2019;Haritonova, 2016a) are shown with the velocities from the NKG (Nordic Geodetic Commission) deformation model NKG_RF03vel (in yellow) (Milne et al., 2001;Lidberg, 2004). As can be seen, the orientations of the horizontal velocities are different comparing two velocity solutions.  Figure 2, were obtained also using Tsview software, but without removing of annual and semi-annual trends from the residuals.

Latvian sea level data
Hydrological observations in Latvia are carried out by the Latvian Environment, Geology and Meteorology Centre (LEGMC). It is maintaining totally 9 coastal hydrologic stations on the Latvian coast, which provide continuous observations. For this study records of 7 Latvian tide gauges were used to perform correlation analysis: 5 of them are located on the coast of the Gulf of Riga (Kolka, Mērsrags, Lielupes grīva, Daugavgrīva and Salacgrīva), and 2on the coast washed by the open sea (Liepāja and Ventspils). Data are daily mean values of hourly records. The time series have a distinct increase in the amplitude during the winter season: the reasonwinds; this effect was pointed out by Ekman (1998) (Haritonova, 2016b). The short-term variations of the Baltic Sea (less than one month) are mainly internally driven variations, with maximum amplitudes in the far north and the far south, and a nodal line close to Stockholm in the middle. Thus short-term sea level variations are nearly eliminated at this site (Ekman, 2009). This is true for the Latvian coast as well, which is also located close to the nodal line of the Baltic Sea level variations. Concerning the Gulf of Riga, the range of sea level variations runs up to 1.50 m and higher. Figure 3 presents overall picture of dynamics of the sea level at the Latvian coast and used data for the correlation analysis.

Correlation between GNSS station and tide gauge observation data
For each observation year (2012-2018) and each GNSS station (33) correlation coefficients were obtained using GNSS station position residuals in North, East and Up components and sea level data of 7 tide gauges. To visualize results, the maximum values of correlation coefficients were used (see Figure 4). The results have a high probability of being statistically significant, as number of data pairs in any case is higher than 300 (during one year). For example, the correlations greater than 0.20 has less than 0.05 % probability of occurring randomly. Figure 6 presents statistics of 7 years for the Latvian and two Estonian GNSS stations in North, East and Up components, and Figure 7 displays correlation between observations of the year 2018 in North, East and Up components as well.
Both negative and positive correlations are observable in the Up component (see Figure 6 and 7). Mostly correlation is weak (correlation coefficient is less than |0.40|), except four stations: IRBE, KUL1/2, LIPJ and RUHN (see Figure 4). The observations at these sites are mostly moderately correlated. The maximum correlation coefficient in the Up component -0.60 is observable in the case of station RUHN in 2018, which is located on the Estonian Island Ruhnu in the Gulf of Riga. According to the data shown in Figure   The maximum correlation coefficients in the East component are observable in the case of station VAL1: -0.78, and IKLA: -0.80 in 2018. It is interesting that one site (station SLD1) shows inverse correlation for all years (see Figure 6 and 7).
In Figure 8 the results with strong negative correlation (-0.78) in the East component at station VAL1 are shown (year 2018). Green line presents GNSS station coordinate residuals. For better interpretation the sea level residuals (at tide gauge Ventspils), shown in black, are divided by factor 250 in the top graph; below the sea level data are shown with opposite sign, like as site displacement caused by the loading effect. The strong correlation between two data sets is obvious. Removing this signal, it would be possible to reduce the scatter of the GNSS time series.

Conclusions
The results presented in this paper show that the impact of the Baltic Sea non-tidal loading in the territory of Latvia is obvious. It can be assumed; this is due to the geographical location of the territory of Latvia.
The correlation between GNSS station daily coordinate time series and sea level data for the observation period from 2012 up to 2018 is well observable in North, East and Up components. Both negative and positive correlations are observable in the Up component. Sites located closer to the open sea display the highest negative correlations of all. The maximum correlation coefficient in the Up component -0.60 is observable in the case of station RUHN in 2018, which is located on the Island Ruhnu in the Gulf of Riga. The correlation is most prominent in the East componentit is observable for the most stations with highest values during 7-year observation period, and with maximum correlation coefficients -0.78 at VAL1 and -0.80 at IKLA in 2018.
Less pronounced correlation is observable in the North component, however, despite the fact that correlations are changing from year to year, even weak correlations at sites display regionally related results.
The unexplained fact is the inverse correlation in the East component at station SLD1 during the entire observation period comparing with the results of other sites.
Atmospheric loading is out of focus in this study. Suppose that the atmospheric surface pressure anomaly is near-uniform, i.e. it produces a near-uniform vertical deformation and near-zero horizontal deformations. According to this, the results obtained for the Up component could be affected to some extent by this impact.
The research will be continued possibly using modelled loading time series, where atmospheric loading will be taken into account.