## Assessment of Hydrological Regime in the Long Run of the Danube River Using Nonlinear Times Series Analysis and Wavelet Transform

Srđan Kostić^{1}, Milan Stojković^{1}, Nebojša Vasović^{2}

^{1} Institute for the Development of Water Resources "Jaroslav Černi", Jaroslava Cernog st. 80, Belgrade, Serbia; E-mail:
This e-mail address is being protected from spambots. You need JavaScript enabled to view it

^{2} Faculty of Mining and Geology University of Belgrade, Đušina 7, 11000 Belgrade, Serbia

### Abstract

In present paper the authors perform complex research of the annual flows and the climatic time series. In particular, the annual flows of the Danube River, recorded at the Orsova station in the period 1840-2006, are analyzed by means of a nonlinear time series analysis and wavelet discrete transform, while the latter is also used for the analysis of climatic drivers, i.e. rainfall and temperature, observed at the Bratislava climatic station. Results of the performed research indicate that the original unprocessed flow time series exhibit the properties of random fluctuations, as a consequence of measurement noise and inherent complexity of the river flow itself. For this reason, authors invoke the method of wavelet discrete transform, with the idea to decompose the observed time series, in order to obtain more regular oscillations, which could directly be linked to the corresponding climatic drivers. The obtained results demonstrate a statistically insignificant negative trend in annual flow time series, while the approximation component of annual flows implies a strong negative trend. One should note that approximation components of annual rainfall and air temperature show the same behavior as the corresponding component in annual flows. The results of this study may be utilized for establishing convenient models to simulate flow time series based on the climatic records, by coupling the wavelet discrete transform and nonlinear regression or machine learning.

**Keywords:** wavelet transform, long-term flow behavior, nonlinear time series analysis, the Danube River.

Introduction

Analysis of long-term properties of hydrological regimes has an important role in water management systems, since it enables development of effective strategies to cope with the impact of extreme climate changes, by means of deriving reliable models for estimating the flow trend in the future. In particular, variability of annual flows could have a significant effect on hydroelectric power generation, since it predetermines the water availability in river basins (Sale et al. 2012; Shrestha et al. 2014; Kao et al. 2015; Dahal et al. 2016). In this paper the authors analyze the long-term hydrological properties of the Danube River using observed flows and climatological data at Orsova and Bratislava station, respectively.

There are many models for flow estimation that are based on long-term flow statistics (Pekarova et al., 2003, Pekarova and Pekar, 2006; Stojković et al., 2016) which preserve the characteristics of long-term hydrologic processes. Here, the authors use two different techniques for assessing the long-term characteristics. One is based on the wavelet transform, and the other is based on the application of nonlinear time series analysis. The principle idea is to characterize the observed time series, both from the hydrological viewpoint (using wavelet transform) and from a purely dynamical point of view (using techniques of nonlinear time series analysis).

Wavelet transform has recently been used largely in hydrological analyses. In particular, Tran et al. (2016) used wavelet transform in order to increase the accuracy of autoregressive time series forecasting. Apparently, a hybrid technique was developed as a combination of significant input variable selection using partial linear correlation and input data transformation and discrete wavelet transform for decomposing the input time series into low and high frequency components. On the other hand, Nourani et al. (2011) coupled the wavelet transform to the adaptive neural fuzzy inference system, and showed that the developed model is more appropriate (compared to the artificial neural network model) for long-term runoff discharges, since it uses the multi-scale time-series of rainfall and runoff data in the input layer. In this paper, the method of wavelet transform is engaged in order to decompose the observed time series into components suitable for linking with climate drivers.

Methods of nonlinear time series analysis have been extensively applied in hydrological analyses, usually ''hidden'' under the term ''chaos theory''. Jayawardena and Lai (1994) were the first to apply the nonlinear time series analysis, and showed that the time delay embedding approach provides more accurate results in comparison with the traditional linear ARMA method. Sivakumar et al. (1999) used methods based on nonlinear time series analysis in order to reduce measurement noise in the real observed hydrological data. Moreover, Sivakumar et al. (2001) indicated a possible presence of deterministic chaos in the rainfall-runoff process. Zhou et al. (2002) used a complex approach including power spectrum analysis, phase space reconstruction and estimation of attractor dimension, in order to examine the possible chaotic dynamics of the flood series in the Huaihe River basin for the last 500 years. Results of their analysis indicated that the power spectrum of the analyzed time series is similar to that of the chaotic time series, while the attractor dimension is 4.66. Ng et al (2007) applied chaotic analytical techniques to daily hydrological series comprising of outliers. Obtained results indicated that the examined time series exhibits random-like fluctuations. Labat et al. (2016) employed a nonlinear dynamic method and the correlation dimension method to unique long, continuous, and high-resolution (30-min) streamflow data from two karstic watersheds in the Pyrénées Mountains of France. Results of the performed research implied the presence of deterministic chaos in the streamflow dynamics, with attractor dimension values below 3. This paper demonstrates how nonlinear time series analysis is used to examine the real observed time series (annual flows of the Danube River), so as to confirm or reject the possible regularity. It should be noted that apart from the application of nonlinear time series analysis in hydrology, this method is widely used in other areas of Earth Science (Kostić et al. 2013).

The structure of the paper is as follows: Applied methodology is described in Section "Methodology", while the analyzed data sets are described in Section "Data set". Section "Results" is devoted to the main results, while the final remarks with suggestions for further research are given in the final section "Conclusions".

Methodology

Assessment of the long-term changes in hydroclimatic time series is conducted by using the following procedure:

- Nonlinear time series analysis is engaged in order to characterize the observed flow time series;
- The recorded annual hydro-meteorological time series are decomposed into their approximation (a) and harmonic components - details (d) by means of the wavelet transform.
- Mann-Kendall test is applied for analysis of the recorded annual time series, including their decomposed components (approximation and details) so as to assess their monotonic trends.
- For each observed and decomposed time series the correlation coefficient is determined between the climatic factors (rainfall and air temperature) and the flows, so as to link the long-term changes in hydrological regime with simultaneous changes in climatic drivers.

Nonlinear time series analysis

Nonlinear time series analysis of the observed flow of the Danube River at Orsova station is performed through the following steps:

- Surrogate data analysis, which is performed by testing the three null hypotheses: (a) data are independent random numbers drawn from some fixed but unknown distribution; (b) data originate from a stationary linear stochastic process with Gaussian inputs and (c) data originate from a stationary Gaussian linear process that has been distorted by a monotonic, instantaneous, time-independent nonlinear function (Perc et al., 2008). Surrogates are generated using Matlab toolkit MATS (Kugiumtzis and Tsimpiris, 2010), while the zeroth-order prediction error is calculated based on the suggestion of Kantz and Schreiber (2004). Calculations in Matlab were conducted at the University of Belgrade Faculty of Mining and Geology.
- Determinism test, which is based on the assumption that if the time series originated from a deterministic system, the obtained vector field should consist solely of vectors that have unit length, indicating the average length of all directional vectors k to be equal to 1 (Kaplan and Glass, 1992). For calculating the optimal embedding delay, authors use the average mutual information method [Fraser and Swinney, 1986], which utilizes the first local minimum of mutual information as an optimal embedding delay. Minimal embedding dimension m is obtained using the procedure of false nearest neighbor detection suggested in (Kennel et al, 1992), based on the criterion that the normalized distance between the embedding coordinates of two presumably neighboring points is larger than a given threshold (in this case
*R*=10)._{tr}

Wavelet transform

The continuous wavelet transform is mainly used in order to reveal the characteristics of the hydrological time series under multi-temporal scales (Sang, 2013). It uses a mother wavelet ψ*(t)* across time series *f(t)* (Smith et al. 1998):

(3)

where

(4)

represents a family of functions commonly called wavelets. In eq. (3) *λ* is a scale parameter, *t* is a location parameter, and is the complex conjunction of *ψ(u)*.

Since hydrological time series are discrete signals, the discrete form of wavelet transform is given as follows:

(5)

where a_{0} and b_{0} are constants, integer *j* is a decomposition level, and k is a translation factor. For instance, the discrete wavelet transform is usually used by assigning a_{0}=2 and b_{0}=1 (Daubechies, 1992):

(6)

The decomposition of the recorded time series starts from the original series, and the result includes two types of functions called approximation and detail. The selection of wavelet and choice of decomposition level are two key factors in application of wavelet transform for the analysis of hydrological time series. Level of decomposition is determined by the periodicity of the corresponding component of the examined time series.

There are a number of wavelet families which are commonly used in the discrete wavelet transform. Some of them are given as follows: Daubechies (dbN), Coiflets (coifN), Symlets (symN), BiorSplines (biorM.N), ReverseBior (rbioM.N) and DMeyer (dmey), the latter of which is used in the present research.

Mann-Kendall test

Many trend analyses in hydrology are based on the non-parametric Mann-Kendall trend test (Douglas et al. 2000). This test represents a rank-based procedure which is not sensitive on extreme values in time series and does not require normal distribution of the examined time series (Partal, 2017).

In the first step, one needs to define Kendall's statistic *S*, which is computed in the following way (Kendall 1962):

(7)

where *Q* is the mean of the annual hydroclimatic time series at time steps *i* and *j*, n is the size of the time series, and *sign*(.) is equal to +1 if *Q _{i} > Q_{j}* and -1 otherwise. If

*S*> 0, the time series has a downward trend and if

*S*< 0, then there is an upward trend.

The standard normal variable zs characterized with zero mean and unit standard deviation is determined as follows (Douglas et al. 2000):

(8)

This standardized variable follows the standard normal distribution. Therefore, for *z _{s}* < 1.96, which corresponds to the significance level of 5%, the hypothesis

*H*on the absence of a monotonic trend in mean flows is not rejected. Otherwise, for

_{0}*z*> 1.96 the alternative hypothesis

_{s}*H*is that there is a trend in the series.

_{1}

Data set

Data set for the analysis is provided on the basis of cooperation between the Jaroslav Černi Institute from Belgrade (Serbia), the Institute of Hydrology SAS from Bratislava, (Slovakia) and the Global Runoff Data Centre within the German Federal Institute of Hydrology from Koblen (Germany).

In order to reveal the long-term changes in hydrological patterns, the authors choose the annual flows of the Danube River observed at the Orsova hydrological station (h.s.), since it is one of the largest observed series in a database. Furthermore, to associate such changes with meteorological drivers in the river basin, i.e. annual precipitation sums and mean annual air temperature, observations made at the Bratislava climatic station (c.s.) are utilized. Main data on the examined time series are given in Table 1.

Table 1: Hydro-meteorological time series used in the study within the available period, mean annual values and standard deviation.

Results

Nonlinear time series analysis

In order to test the null hypothesis that the data are independent random numbers drawn from some fixed but unknown distribution, 20 surrogates are generated by randomly shuffling the data (without repetition), thus yielding surrogates with exactly the same distribution yet independent construction. Then the zeroth-order prediction error is calculated for the original dataset (ε_{0}) and for each of the 20 generated surrogates (ε). It is clear from Figure 1 that ε_{0} is well within ε in all the cases, so the null hypothesis cannot be rejected.

Further confirmation of the results of the surrogate data testing could be obtained by a determinism test, through coarse-graining of the embedding space into 26 boxes in one dimension, of which only those boxes visited more than one time by the trajectory are included in the analysis. So as to be able to perform the determinism test, one firstly needs to embed the observed scalar time series into appropriate phase space. Regarding the optimal values of embedding parameters, the performed analysis showed that mutual information takes the first local minimum for s = 2 (Figure 2), while the fraction of the false nearest neighbors rises with the increase of embedding dimension, which could indicate random signature in the input data. Nevertheless, so as to be able to perform the determinism test, m = 3 is chosen as an optimal embedding dimension, since in this case, river flow is assumed to have three degrees of freedom (flow magnitude, rainfall and temperature).

Once the optimum embedding delay and dimension are determined, one could apply the determinism test, which, in the present case, resulted in the value of determinism factor k = 0.76 indicating possible randomness in the observed flow time series.

Figure 1: Surrogate data test for the hypothesis that the data are independent random numbers drawn from some fixed but unknown distribution. The red line denotes the zeroth-order prediction error for the original series and the black lines – zeroth-order prediction error for the surrogates, for n prediction units. The obtained results indicate that the error for the original dataset (ε_{0}) is within the error for surrogate data (ε), following that the null hypothesis could not be rejected.

Figure 2: Determination of the proper embedding delay – mutual information has the first minimum at minimum embedding delay *s* = 2, meaning that the values of the flow at dimensionless discrete time units, t + 2, t + 4, … (depending on the value of optimal embedding dimension) could be used for the embedding of the observed dataset in phase space.

Wavelet decomposition

The discrete wavelet decomposition is applied for the analysis of hydrological records of the Danube River at Orsova h.s. Such transformation enables one to examine the long-term behavior of the decomposed components at different time scales. In present research, the observed annual flows are decomposed into a long period approximation component (function) and five irregular and quasi-regular components (functions) commonly called details, namely d_{1}, d_{2}, d_{3}, d_{4} and d_{5} with characteristic periodicity of 2 years, 4 years, 8 years, 16 years and 32 years, respectively (Figure 3). The approximation component (a_{5}) is a residual part of annual flows at the fifth level of applied decomposition.

The same procedure is employed to detect changes in annual rainfall and air temperature, which are considered to be the most influential runoff drivers. the authors chose Bratislava c.s. since it is located in the Danube river basin and has a long record of observations as those at Orsova h.s. The results of decomposition for annual rainfall and air temperature are given in Figures 4 and 5, respectively.

Figure 3: Decomposition of annual flows at Orsova h.s. (the Danube River) into approximation component (a_{5}) and details (d_{1}, d_{2}, d_{3}, d_{4} and d_{5}) by applying the discrete wavelet transform.

Figure 4: Decomposition of annual precipitation at Bratislava c.s. into approximation component (a_{5}) and details (d_{1}, d_{2}, d_{3}, d_{4} and d_{5}) by applying the discrete wavelet transform.

Figure 5: Decomposition of annual temperature at Bratislava c.s. into approximation component (a_{5}) and details (d_{1}, d_{2}, d_{3}, d_{4} and d_{5}) by applying the discrete wavelet transform.

As it can be seen from Figure 3, the most influential detail is d_{1} with a 39% share of the total observed variance of annual flows at Orsova h.s. The lowest contribution to the observed variance comes from detail d_{5} (7.3%) with a periodicity of 32 years. Similar long-term behavior, as that detected for annual flows, is revealed for annual rainfall at Bratislava c.s (Figure 4). Actually, the greatest part in total variance comes from detail d_{1} (44%), while detail d_{5} (2.1%) possesses the lowest share in overall variance.

Figure 5 shows the observed and decomposed annual time series of air temperature at Bratislava c.s. The result suggests an upward trend in approximation component (a_{5}), while the rest of the decomposed time series visually implies only the long-term periodical behavior of annual rainfall.

Once the discrete wavelet decomposition has been applied for the analysis of hydro-meteorological time series, authors use the Mann-Kendal test to examine the trend in the observed annual flows (Q) at Orsova h.s, including the observed annual rainfall (R) and air temperature (T) at Bratislava c.s. Additionally, the trend test is also applied for the analysis of the approximation components and details of the examined time series. The results obtained are given in Table 2.

Table 2: Application of the Mann-Kendal test for the analysis of observed and decomposed (approximation component and details) annual flows (Q) at Orsova h.s. and annual rainfall (R) and air temperature (T) at Bratislava c.s.

The observed flows and rainfall have a decreasing tendency but trends are not statistically significant at significance level α = 5%. Nevertheless, in case of the observed annual air temperature results seem different - annual air temperature shows a statistically significant trend for a confidence interval of 95%.

One should note that approximation component a5 of flow time series exhibits a statistically significant downward trend at significance level α = 5%. This fact implies that negative tendency in annual flows could be clearly revealed by a strong negative trend in the approximation component. This could be at least partially explained by the effect of rainfall and air temperature, since a statistically significant negative trend in the approximation component of annual rainfalls, and a very strong positive trend in the approximation component of annual air temperature are also detected.

It should also be emphasized that approximation component a_{5} of annual flow and annual rainfall shows significant linear correlation (r = 0.702), while less pronounced correlation is found in the case of the flow and air temperature approximation component. Separate analyses of the periodic wavelet components of the examined time series suggest that approximation components contain the predominant information on the observed trends compared with the details. The detected long-term changes in the approximation component of annual flows could be partially attributed to the effect of rainfall, whose approximation component exhibits qualitatively similar changes.

Conclusion

In this paper, the authors examined the observed data on the annual flow of the Danube River, recorded at the Orsova hydrological station, for the period 1840-2006, in order to establish a link with climatic drivers, i.e. rainfall and temperature, recorded at the Bratislava climatic station. The analysis consisted of three successive stages: (1) in the first stage, a nonlinear time series analysis was applied, so as to characterize the original recorded times series on the annual flow; (2) in the second stage, observed time series were decomposed using the wavelet transform, (3) in the final stage, changes in decomposed time series were linked to the corresponding changes in climatic drivers.

Results of nonlinear time series analysis indicate the presence of random fluctuations in the observed flow data, which could be the consequence of the measurement noise and the complexity of the recorded oscillations. Considering this, it was not possible to proceed with the modeling of the unprocessed original recorded data, but one needed to invoke methods for decomposing the time series into a convenient set of corresponding oscillations. Here, the authors deploy the method of wavelet transform. The annual records of flows, rainfall and air temperature at Orsova h.s. and Bratislava c.s. are decomposed into the approximation component and five details, whereby details d_{1}, d_{2}, d_{3}, d_{4} and d_{5} represent functions with a periodicity of 2 years, 4 years, 8 years, 16 years, and 32 years, respectively. Once the discrete wavelet decomposition has been applied for the analysis of hydro-meteorological time series, authors invoke the Mann-Kendal test in order to examine the monotonic trend in the observed time series and their approximation components and details. The approximation component a_{5} of the flow time series suggests a statistically significant downward trend at significance level α=5%. This fact implies that negative tendency in annual flows is clearly visible as a negative trend in the approximation component of the observed time series. Also, the behavior of annual flows is accompanied by a statistically significant negative trend in the approximation component of annual rainfall and a strong upward trend in the approximation component of annual air temperature.

Regarding the benefit of the results for future research, obtained results could be further used to create convenient models coupled with artificial neural networks or adaptive neural fuzzy inference system, where certain composing time series could be directly linked to corresponding discretization of the rainfall and temperature.

Acknowledgements

This research was partly supported by the Ministry of Education, Science and Technological Development of the Republic of Serbia, Contract No. TR 37005.

References

Dahal, V., Shakya, N.M. and R. Bhattarai. (2016). Estimating the impact of climate change on water availability in Bagmati Basin, Nepal. Environ Processes, 3, 1–17.

Daubechies, I. (1992) Ten Lectures on Wavelets. SIAM, Philadelphia, PA.

Douglas, E.M., Vogel, R.M. and C.N. Kroll. (2000). Trends in floods and low flows in the United States: impact of spatial correlation. Journal of Hydrology 240, 90-105.

Fraser, A. and H. Swinney. (1986). Independent coordinates for strange attractors from mutual information. Phys Rev A. 33, 1134–1140.

Jayawardena, A.W. and F. Lai. (1994). Analysis and prediction of chaos in rainfall and stream flow time series. Journal of Hydrology. 153, 23-52.

Kantz, H. and T. Schreiber. (2004). Nonlinear Time Series Analysis. Cambridge: Cambridge University Press. 388 p.

Kao, S.C., Sale, M.J., Ashfaq, M., Uria Martinez, R., Kaiser, D., Wei, Y. and N.S. Diffenbaugh. (2015). Projecting changes in annual hydropower generation using regional runoff data: an assessment of the United States federal hydropower plants. Energy, 80, 239–250.

Kaplan, D.T. and L. Glass. (1992). Direct test for determinism in a time series. Phys Rev Lett. 68, 427–430.

Kendall, M.G. (1962) Rank Correlation Methods. 3rd ed. Hafner Publishing Company, New York.

Kennel, M., Brown, R. and H. Abarbanel. (1992). Determining embedding dimension for phase-space reconstruction using a geometrical construction. Phys Rev A. 45, 3403–3411.

Kostić, S., Vasović, N., Perc, M., Toljić, M. and D. Nikolić. (2013). Stochastic nature of earthquake ground motion. Physica A. 392, 4134–4145.

Kugiumtzis, D. and A. Tsimpiris. (2010). Measures of Analysis of Time Series (MATS):A MATLAB Toolkit for Computation of Multiple Measures on Time Series Data Bases. J Stat Soft. 33, 1–30.

Labat, D., Sivakumar, B. and A. Mangin. (2016). Evidence for deterministic chaos in long-term high-resolution karstic streamflow time series. Stochastic Environmental Research and Risk Assessment. 30, 2189–2196.

Ng, W.W., Panu, U.S. and W.C. Lennox. (2007). Chaos based Analytical techniques for daily extreme hydrological observations Journal of Hydrology. 342, 17– 41.

Nourani, V., Kisi, O. and M. Komasi. (2011). Two hybrid artificial intelligence approaches for modeling rainfall–runoff process. J. Hydrol. 402, 41–59.

Partal, T. (2017). Multi-annual analysis and trends of the Temperatures and Precipitations in West of Anatolia. Journal of Water and Climate Change. DOI: 10.2166/wcc.2017.109

Pekarova, P. and J. Pekar. (2006). Long-term discharge prediction for the Turnu Severin station (the Danube) using a linear autoregressive model. Hydrol Process, 20, 1217–1228

Pekarova, P., Miklanek, P. and J. Pekar. (2003). Spatial and temporal runoff oscillation analysis of the main rivers of the world during the 19th–20th centuries. J. Hydrol. 274, 62–79.

Perc, M., Green, A.K., Jane Dixon, C. and M. Marhl. (2008). Establishing the stochastic nature of intracellular calcium oscillations from experimental data. Biophys Chem. 132, 33–38.

Sale, M.J., Kao, S.C., Ashfaq, M., Kaiser, D.P., Martinez, R., Webb, C. and Y. Wei. (2012). Assessment of the Effects of Climate Change on Federal Hydropower. Technical Manual 2011/251. Oak Ridge National Laboratory, Oak Ridge.

Sang, Y.F. (2013). A review on the applications of wavelet transform in hydrology time series analysis. Atmospheric Research, 122, 8–15.

Shrestha, S., Khatiwada, M., Babel, M.S. and K. Parajuli. (2014). Impact of climate change on river flow and hydropower production in Kulekhani hydropower project of Nepal. Environ Processes, 1, 231–250.

Sivakumar, B., Phoon, K.-K., Liong, S.-Y. and C.-Y. Liaw. (1999). A systematic approach to noise reduction in chaotic hydrologicaltime series. Journal of Hydrology 219, 103–135.

Sivakumar, B., Berndtsson, R., Olsson, J. and K. Jinno. (2001). Evidence of chaos in the rainfall-runoff process. Hydrological Sciences Journal. 46, 131-145.

Smith, L.C., Turcotte, D.L. and B.L. Isacks. (1998). Stream flow characterization and feature detection using a discrete wavelet transform. Hydrological Processes, 12, 233-249.

Stojković, M., Plavšić, J. and S. Prohaska. (2016). Annual and seasonal discharge prediction in the middle Danube River basin based on a modified TIPS (Tendency, Intermittency, Periodicity, Stochasticity) methodology. J. Hydrol. Hydromech. 65, 165–174

Tran, H.D., Muttil, N. and B.J.C. Perera. (2016). Enhancing accuracy of autoregressive time series forecasting with input selection and wavelet transformation. Journal of Hydroinformatics 18, 791-802.

Zhou, Y., Ma, Z. and L. Wang. (2002). Chaotic dynamics of the flood series in the Huaihe River Basin for the last 500 years. Journal of Hydrology. 258, 100-110.