## 4 Methods and Preliminary Analysis  
Analysis of real data usually comes with a number of challenges. Especially in Earth science, automatic data analysis should be preceded by careful validation of the raw data (e.g., missing data, disturbances, overdrive, interference, range of data). Any additional knowledge about the assumed behavior (e.g., seasonal, cyclic, random, deterministic) is valuable. The necessary part of the analysis and its results require a deep knowledge not only in mathematics and signal processing, but most of all in the specific field of Earth science. In this case, hydrogeology and meteorology are these two specific fields. The key challenge in this paper is to analyze signals related to groundwater level. The first task is related to the hydrological year (HY) self-similarity. Therefore, the signal has been segmented according to commonly accepted dates (black lines), and the winter and summer seasons have been marked in red (Figure 6). A key question relates the similarity of consecutive hydrological years in the investigated period. It can clearly be seen that in each consecutive year, the patterns of the groundwater level fluctuations differ from the previous year and that the time series of the groundwater levels for some HYs is similar and for some is substantially different (Figure 6). As also can be seen (Figure 6), the groundwater level data show a positive linear trend. Consequently, the next step was a detrending operation, which is commonly used as a preprocessing step, especially when looking for similarities between segmented data. In order to remove the trend, a linear dependence was identified and then subtracted from the whole time series of the groundwater level. If we denote the time series of the groundwater level as \(\{y_{t}\}_{t=1}^{m}\) (\(m\)-total number of observations), then the linear trend is defined by the following formula:  
\[\hat{y}_{t}=\hat{a}t+\hat{b},\ t=1,2,...,m. \tag{1}\]  
where \(\hat{a}\) is the estimator of the slope and \(\hat{b}\) is the constant coefficient in the linear regression for the time series \(y\). The \(\hat{a}\) and \(\hat{b}\) parameters are calculated by using the linear least-squares algorithm resulting in a detrended groundwater level curve (Figure 7). Figure 6: Time series of the groundwater level. The hydrological year is marked by the black line. Additionally, the summer/winter seasons are indicated with red broken lines. In order to find the relation between the HYs, we initially applied the Pearson correlation [52, 53]. However, because the interest is in the spatial dependence (dependence between time series corresponding to different HYs), the Pearson correlation matrix for segmented, detrended groundwater level signals was used. This matrix presents a Pearson correlation of time series \(\{z_{t}^{i}\}_{t=1}^{n}\) and \(\{z_{t}^{j}\}_{t=1}^{n}\) for all possible \(i,j=1,2,...,15\), \(n\): number of observations of the HY. There \(\{z_{t}^{i}\}_{t=1}^{n}\) is the vector of observations after the detrending process for corresponding \(i\)-th HY. The Pearson correlation coefficient \(\rho(z^{i},z^{j})\) of two time series \(\{z_{t}^{i}\}_{t=1}^{n}\) and \(\{z_{t}^{j}\}_{t=1}^{n}\) is a measure of their linear dependence and is defined as follows:  
\[\rho(z^{i},z^{j})=\frac{\sum_{t=1}^{n}(z_{t}^{i}-\overline{z}^{i})(z_{t}^{j}- \overline{z}^{j})}{\sqrt{\sum_{t=1}^{n}(z_{t}^{i}-\overline{z}^{i})^{2}}\sqrt {\sum_{t=1}^{n}(z_{t}^{j}-\overline{z}^{j})^{2}}}, \tag{2}\]  
where \(n\) is the sample size and \(\overline{z}^{i}\) and \(\overline{z}^{j}\) are the sample means of \(\{z_{t}^{i}\}_{t=1}^{n}\) and \(\{z_{t}^{j}\}_{t=1}^{n}\), respectively. The Pearson correlation coefficient \(\rho(z^{i},z^{j})\in[-1,1]\), which indicates a strong positive linear dependence when 1 and a strong negative linear dependence when -1.