## 4 Time Series Analysis Method
### Spectral Analysis  
The sets of homogenised values \(T\) and \(H(T)\) are considered the only realisation of the random process, which is considered to be ergodic. Preliminary analyses of the observation series indicated that the process is not normal and that it contains a determining component. Drawing on the physical characteristics of the process, it can be expressed in terms of a function of time \(t\):  
\[H_{i}(T)=\sum_{k=0}^{K}A_{k}t_{i}^{k}+\sum_{j=1}^{J}R_{j}cos\Big{[}2\pi f_{j}t_ {i}-\ ^{0}\phi_{j}\Big{]}+s_{i}(t) \tag{1}\]where \(s_{i}(t)\) has a normal distribution with dispersion \(\sigma^{2}\) and a mean value of \(\mu_{s}(t)=0\), \(A_{k}\) are polynomial coefficients, \(f_{j}=1/P_{j}\) is the frequency, \(P_{j}\) the period (in days), \(R_{j}\) the amplitude, and \(\ {}^{0}\phi_{j}\) the phase. Our method sought to seek the optimal estimate of these parameters, preserving the ergodic conditions and the realisation of properties (1) of the process. We also introduced the notations:  
\[\sum_{k=0}^{K}A_{k}t_{i}^{k}=Q_{k},\,R_{j}cos\Big{[}2\pi f_{j}t_{i}-\ ^{0}\phi_{j} \Big{]}=Q_{j} \tag{2}\]  
In seeking the hidden periods, periodograms of the type:  
\[I(f)=(2\pi T)^{-1}\bigg{|}\int_{0}^{T}x(t)\ exp(-itf)dt\bigg{|}^{2}, \tag{3}\]  
are usually used as a criterion in which \(f\) is the frequency, and \(x\) the continuous realisation of the random process in the interval \(T\), \(i=\sqrt{-1}\). However, we used a slightly different procedure for seeking hidden periods, which was [PERSON]'s anharmonic analysis [24], as modified by [PERSON] [25] and [PERSON] [26].  
According to [PERSON] et al. [27], the optimum spectrum (normalised to the interval <0, 1>), using (2) reads:  
\[\sigma(x_{i})=1-\big{[}\varrho\big{(}x_{i}-Q_{k},\ Q_{j}\big{)}/\varrho(x_{i },\ Q_{k})\big{]} \tag{4}\]  
where:  
\[\varrho(F,G)=\sum_{t\in T}[F(t)-G(t)]^{2} \tag{5}\]  
defines the metric of the two functions \(F\) and \(G\) defined by \(T\). The normalised spectrum enables an unambiguous determination of the significance of individual periods according to the selected fixed criterion.  
Figure 1: Horizontal (E-W, N-S) and vertical coordinate differences at the EYPA station—the daily PPP solution smoothed by the anharmonic analysis. The coordinates are positive towards East, North, and Up and shifted for better visualisation. The time scale is in the Modified Julian Date (MJD).  
The coefficients \(A_{k}\) of polynomials \(Q_{k}\) are determined under the assumption of independence, so that:  
\[\varrho(x_{i},\;Q_{k})=min. \tag{6}\]  
and, similarly, coefficients \(R_{j}\) and \(\ {}^{0}\phi_{j}\) from \(Q_{j}\) are determined for every frequency \(f\in\Omega\), where \(\Omega\) is a set of reasonably optional frequencies, so that:  
\[\varrho(x_{i}-Q_{k},\;Q_{j})=min \tag{7}\]  
Coefficients \(R_{j},j=1,\ldots,J\), are therefore functions of frequency \(f_{j}\), and as regards to \(\sigma(x_{i})\), it is easy to prove that for the given realisation of the process:  
\[\sigma(x_{i}):\Omega\;\rightarrow\;\langle 0\;,\;1\rangle \tag{8}\]  
The separate maxima of function \(\sigma(x_{i},\;t)\), the spectrum peaks, determines the frequencies which make the principal contribution to the overall variation of the function \(x_{i}(t)-Q_{k}(t)\). Conditions (6) and (7) lead to the traditional method of least squares adjustment. The value of \(K\) is determined in (1) and (2) from the already known character of the function or the empirical computation.  
The described method entails finding the most significant wave first and filtering it by subtracting it from the original data. The difference is used as new input data and the next most significant wave is searched for. This is repeated until the maximum value of all found waves of the normalised spectrum is greater than the specified fixed criterium (e.g., 0.3).