#### Abstract

We study the multifractal properties of water level with a high-frequency and massive time series using wavelet methods (estimation of Hurst exponents, multiscale diagram, and wavelet leaders for multifractal analysis (WLMF)) and multifractal detrended fluctuation analysis (MF-DFA). The dataset contains more than two million records from 10 observation sites at a northern China river. The multiscale behaviour is observed in this time series, which indicates the multifractality. This multifractality is detected via multiscale diagram. Then we focus on the multifractal analysis using MF-DFA and WLMF. The two methods give the same conclusion that at most sites the records satisfy the generalized binomial multifractal model, which is robust for different times (morning, afternoon, and evening). The variation in the detailed characteristic parameters of the multifractal model indicates that both human activities and tributaries influence the multifractality. Our work is useful for building simulation models of the water level of local rivers with many observation sites.

#### 1. Introduction

Long-range dependence (LRD) and multifractality are the inherent characteristics of many natural phenomena. Long-range dependence (LRD) was found in the flood process by Hurst [1] when he was studying the flow of the Nile River. Fractional Brown motion (FBM), which is introduced by Mandelbrot and Van Ness [2], is the commonly used model for studying long-range dependence. Long-range dependence is characterized by the Hurst parameter, which describes the intensity of long memory phenomena. Many studies [3–6] for estimating the Hurst parameter (or Hurst exponent) to judge the correlation effect have been made available in recent years. These studies reveal that, in many areas of applied sciences, such as climate change, stock markets, telecommunication network, and river flow, the phenomena appear as long-range behaviour, and the estimation of Hurst has been widely used in decisions and predictions.

The concept of multifractality was originally introduced by Mandelbrot [7] in the context of turbulence. The formalism of multifractality was first introduced using box-counting algorithms for the complexity of measures [8–10] and further generalized to functions using the structure function method [11]. Muzy et al. [12] first developed the multifractal formalism using wavelets in the discussion of turbulence data. Schertzer and Lovejoy [13] concluded that a multifractal field could be regarded as the result of a multiplicative cascade process. Until now, the multifractality and its formalism have been applied in hydrology and many other contexts [14–20].

Previous studies on the LRD and multifractality of hydrology (see, e.g., [14–20]) were based on relatively large sampling periods, such as 1 day, 1 month, or even 1 year, rarely on high-frequency hydrological data. Venugopal et al. [21–23] studied the multifractality of high-resolution temporal rainfall data using wavelet methods. We focus on studying the properties of high-frequency hydrological data and would like to further explore the typical LRD phenomena and multifractal behaviours in hydrology. Our records are the high-frequency water level records (every 6 minutes) of one river (a branch of the Haihe River in China) at its 10 observation sites from April 2011 to September 2013. We have not yet seen studies conducted for such short sampling intervals and such a wide range of observation sites (along the river from upstream to downstream).

It is known that the hydrologic data are often affected by trends or other nonstationarities, for example, due to the seasonal cycle or a change in climate, which may lead to an overestimation of the Hurst exponent. Furthermore, under the influence of the trends, the uncorrelated data may behave as long-term correlated data. Thus, the methods applied in this paper are based on wavelet analysis (estimation of Hurst exponents, multiscale diagram, and wavelet leaders for multifractal analysis (WLMF)) and multifractal detrended fluctuation analysis (MF-DFA). The wavelet method can eliminate some trends as a result of the vanishing moment property, and MF-DFA can avoid spurious detection of correlations that are artifacts of nonstationarity in the records [24]. The wavelet-based estimator of the Hurst exponent was proposed by Veitch and Abry [25] and has been improved in recent years [26, 27]. This estimate was further extended to the situation of random fields [28, 29]. In addition to the estimation, the wavelet-based method can also be used to obtain a multiscale diagram for detecting the multifractal behaviour [8, 9, 26, 30]. The multifractality has been studied in several methods since it was proposed; the primary and advanced methods are wavelet analysis [10, 31, 32] and multifractal detrended fluctuation analysis [24, 33]. The MF-DFA is based on the standard fluctuation analysis (FA) and is a multifractal generalization of the detrended fluctuation analysis (DFA) [24, 33]. Compared with FA, the method can be used to conduct nonstationary data for its detrending function. In recent years, the well-established MF-DFA has been widely used in hydrologic data (see Koscielny-Bunde et al. [14, 16]). Lashermes et al. introduced a new multifractal formalism based on wavelet leaders [34], which allows the multifractal spectrum to be obtained over its entire range. Since then, many studies on multifractality and its applications have been continuously conducted [35]. The analysis method, namely, wavelet leaders multifractal analysis (WLMF), is a discrete version of wavelet transform modulus maxima (WTMM) method [36, 37]. The DFA and WTMM methods are compared in many literatures (see, e.g., [16, 24, 38]). Comparing with WTMM, DFA is easy to implement and has smaller MSE than WTMM for very short time series. But the DFA is lack of deep theoretical concepts underlying the multifractal formalism [38].

In this paper, we first apply the wavelet-based method to estimate the Hurst exponent with different sampling and wavelet scales and observe both the LRD phenomenon and multiscale behaviour in the high-frequency water level records. We think that the multiscale behaviour may induce multifractality. This multifractality is detected via multiscale diagram. Then we apply MF-DFA and WLMF to the records to study the multifractality in detail. The two methods give the same conclusion that the water level records at most sites (except some abnormal sites, where the water levels are seriously affected by human activities) satisfy the generalized binomial multifractal model, which is robust for different times (morning, afternoon, and evening). Another interesting result is that, by using MF-DFA, the values of vary from approximately 1 to less than 0.9; we conclude that this phenomenon occurs because there are four tributaries that flow into the river, changing the multifractality of the river. As a technical improvement in the estimation of the Hurst parameter , we show that the sampling intervals and choice of scales of the discrete wavelet transformation should be matched. Our work is useful for those who want to build simulation models of the water level of local rivers with many observation sites, because the LRD and multifractal properties of the simulation data can be checked by the methods used in our paper.

The logical structure and organization of the paper are described in Figure 1.

The remaining parts of the paper are organized as follows. We introduce the data source and data extraction method for analysis in Section 2. The methods employed in this paper are introduced briefly in Section 3. The detailed analysis and the results for the records are presented in Section 4. In Section 5, we compare the results obtained from Section 4 with the previous results. We conclude our work in Section 6.

#### 2. Data

The records that are collected in real time from one river (a branch of the Haihe River in China) represent the water level. The river is located in northern China, which has a significant continental monsoon climate. It is cold and dry in the winter and warm and humid in the summer. From upstream to downstream, 10 artificial canals have been built to control the water level. Next to these canals, 10 water level observation sites (Site 1–Site 10) have been set up along the river (see Figure 2). Four tributaries are located downstream of Site 7 that flow into this river. The data were obtained from these 10 sites and collected every 6 minutes over the past two and a half years from April 2011 to September 2013. Table 1 shows the sampling information at each site. Figure 3 shows the water level records at Site 1 denoted by and the cumulative sum of . is the mean of . The MF-DFA and WLMF methods require that the records should be a measure and its cumulative sum should be a function. The cumulative sum in Figure 3 displays a function form. This confirms the use of the MF-DFA and WLMF methods. Partial records for validation are available at GitHub: https://github.com/tongzhouzhao/water-level-records.

**(a)**

**(b)**

The water level of the river is always influenced by the weather conditions and human activities. Every year, during the traditional rainy season, the upstream canals are always filled with much rainwater. During the dry season, the upstream canals discharge the stored water. Both the storage and discharge of water should be adjusted at any time according to the growing condition of crops along the bank of the river and the actual demand of human activities. According to our survey, water storage and discharge at Site 2 and Site 7 are conducted frequently for crop growth and human activities.

The following three preprocessing steps were conducted before the multifractal analyses: () removing the outliers and () extracting the appropriate data. These preprocessing steps were conducted to exclude outliers, noise, and artifacts in the data. Firstly, we deleted the outliers caused by the complex environment where the data were acquired. For instance, at 12:42:00 7/18/2013, the water level measured at Site 8 was −1.89, which is obviously an outlier. Secondly, we adjust the artifact data records. The artificial control of the sluice may cause the water level to rise or fall in a short period of time. For example, at 13:12:00 7/5/2013, the record (Site 8) dropped from 27.3 to 0.28. The phenomenon of large fluctuations in a short time may be misdiagnosed, as there is noise in the dataset. Thirdly, we fill the missing data appropriately. There is a small portion of data missing during times when the data-collecting equipment did not work well. For instance, at Site 5, we received one day’s last record at 16:42:00 4/27/2011, which should have been followed by another record six minutes later at 16:48:00 on the same day. However, the next record time in our dataset is 8:48:00 4/29/2011, which means that almost one-day records are missing. A situation like this can occur at other sites. The missing record is directly relevant to the nearby records. For avoiding the effect of time in one day, we use the records at the same time of adjacent days to fill the missing data. The used method is linear interpolation. For example, at Site 5, the record at 16:18:00 4/28/2011 is missing. According to linear interpolation method, it can be filled with the mean of the records at 16:18:00 4/27/2011 and 16:18:00 4/29/2011. The preprocessing may ensure a uniformly spaced sampling interval.

#### 3. Methods

##### 3.1. Wavelet-Based Estimation of the Hurst Exponent

Veitch et al. [25, 26] have proposed the wavelet-based method for many years. Recently, they improved it for non-Gaussian long-range dependent processes [27]. We use this method to estimate the Hurst exponent. The wavelet-based method can detrend the records with polynomial trend because of the property of vanishing moments of wavelet, can avoid the spurious result caused by trend, and is suitable for water level records.

Let , denote a second-order stationary process with a spectrum labelled , which satisfiesThe process is said to be long-range dependent.

Let be a wavelet function with a vanishing moment . The wavelet coefficients of are

Because of the stationary property of , we haveand considering the property of the spectrum density, (3) becomes

As a result,

Putwhere is the number of actually available at scale . can be estimated by . Then, the estimation of the Hurst parameter can be realized with a linear regression.

Because of the different variances of at different values of , the weighted least squares for this regression model are needed. The weight for this regression is the reciprocal of the variance of . Moreover,where . The term can be estimated for each octave by using the sample moment estimators of the fourth and second moments of .

Abry et al. [27] studied the statistical performance of this estimator for non-Gaussian long-range dependent processes obtained from point transformations of Gaussian processes. Their results suggest that the vanishing moment should satisfy . Furthermore, the linear regression should be performed over scales of , ranging from the finest to the coarsest available. To obtain the estimators at different scales, we do not comply with this suggestion in our experiments. It is also shown that the reduction of the dependence range of wavelet coefficients from long to short by increasing above 2, as observed for Gaussian processes, does not hold for non-Gaussian processes. This may lead to the fact that increasing above 2 does not result in any decrease in the variance of the estimator.

##### 3.2. Wavelet-Based Multiscale Diagram

Wavelet procedures provide us with a way to estimate the Hurst parameter. Sometimes, the Hurst index depends on time (i.e., multiple exponents), such as multifractional process and multifractal process. For the multifractal process, its sample path has local regularity:where denotes an element of the probability space [26]. Abry et al. [26] introduced a wavelet-based multiscale diagram for detecting multifractality of process. This method abandons the idea of estimating that changes with and studies the process statistically by Legendre multifractal spectrum, which can be obtained by calculating the structure function (or partition functions) from the sample path:and for the multifractal case, [26, 39]. Moreover, the Legendre transform of the function is the Legendre multifractal spectrum.

Wavelet processes are an effective tool for obtaining because of the fact that the increments of the sample path can be seen as simple examples of the wavelet coefficients. The wavelet coefficients are used to replace the increments in the partition function [9, 26]. Furthermore, the wavelet-based partition function exhibits power-law behaviour for small scales [8, 26, 30, 39]:According to the power-law, the statistics , for small . So can be estimated using linear regression [26].

If (linear form), the process is often called “monofractal,” whose multifractal spectrum is determined by one , such as self-similar process (e.g., FBM) and LRD process. Conversely, if is not linear, the process exhibits multiscale property. The multifractality of a process can be detected by measuring the deviation of from the linear form [26]. Let be the estimation of . The graph of against , together with the confidence intervals related to , is designated as “multiscale diagram (MD)” [26]. The confidence intervals are obtained according to Delbeke and Van Assche [40]. The detection of multifractality can be solved by examining whether there is a lack of alignment in the multiscale diagram.

##### 3.3. Multifractal Detrended Fluctuation Analysis

The wavelet-based multiscale diagram can detect whether the process is multifractal. We study the multifractality by using MF-DFA [24], which is based on the standard fluctuation analysis (FA).

In MF-DFA, one investigates the series:where is the length of and is the mean of .

Let . is divided into nonoverlapping segments with equal length . The fluctuation is calculated usingwhere is the best polynomial fit of in each of the segments .

In FA, the fluctuation , the trend in the original series is not considered. In MF-DFA, the time series is detrended by subtracting . The -order polynomial used in the fitting procedure can eliminate the polynomial trend of order in the original series. This ensures the fluctuation analysis available for data affected by trends or other nonstationarities.

For a multifractal description of records, one must consider all moments :If is long-range dependent, satisfies [24]

For very large , that is, , becomes inaccurate due to the statistical errors produced by the small number of segments . So we will concentrate on the scales lower than . The minimum scale chosen is equal to to avoid spurious results produced by .

can be obtained by analyzing log-log plots versus . If , is the fluctuation exponent that is related to the correlation exponent by [24]. For a stationary time series, is the well-known Hurst exponent. The function is called the generalized Hurst exponent [24].

It is known that [14, 24], for , shows the scaling behaviour of the segments with large fluctuation, and, for , shows that with small fluctuation.

One can also get the Renyi scaling exponent by calculating . The way to describe the multifractality is the singularity spectrum (or multifractal spectrum), that is, the Legendre transform of [14, 16, 24]:where is singularity strength (or Hölder exponent). The strength of the multifractality of records can be characterized by .

If satisfies the following equation,which can be obtained from a generalized binomial multifractal model [16], , . In this case,

The three parameters , and can then be regarded as multifractal finger prints for records [14].

##### 3.4. Wavelet Leaders for Multifractal Analysis

For comparison with MF-DFA, we introduce the method of wavelet leaders for multifractal analysis (WLMF) proposed by Wendt et al. (see [41]).

Assume that is a wavelet function with compact time support and vanishing moment. Let be dyadic intervals and let denote the intervals . The wavelet leaders are defined asThe structure function is calculated by the averages of :It can been shown [41] that, under a mild uniform Hölder regularity condition on ,So the estimation for can be realized by a linear regression in versus :

The weights satisfy and .

To estimate the multifractal spectrum, Wendt et al. [41] suggest a parametric formulation that avoids the computation of the Legendre transform:where , , and .

We emphasize that the wavelet transform with vanishing moments can eliminate the polynomial trend of order in the original series.

#### 4. Computational Results

##### 4.1. Wavelet-Based Results

This section describes the wavelet-based estimation of Hurst parameters for the water level records of the 10 sites. We estimate the Hurst exponents at different scales and different sampling interval times. The results indicate that there is an interesting connection between the wavelet scale and the sampling period. For comparison, we repeat the procedure for a synthetic exponential long-range dependent process (ELRD) with length of 180000, which is close to the actual sampling data length. The vanishing moment of the used Daubechies wavelet is (db3). The simulation results for ELRD are listed in Tables 2 and 3. The values in brackets are errors bars. The errors bars are based on the weighted linear least square fit which is used to estimate and reveal 95% confidence interval.

The simulation results show that the simulations for both and are stable under changes in the interval length and estimation scale. One can see that the wavelet-based estimation of ELRD for the Hurst exponent is valid and reliable, and the selection of wavelet scale and sampling cycle is independent of the estimation.

Tables 4–6 display the results of the practical water level. There exists an expression ( is the element of the tables). For example, in Table 4,* H* = 0.82 (12 min, scales (1~5)) and* H* = 0.84 (6 min, scales (2~6)) are approximately equal;* H* = 0.70 (24 min, scales (1~5)),* H* = 0.68 (12 min, scales (2~6)) and* H* = 0.68 (6 min, scales (3~7)) are almost equal, and so forth. The results represent the consistent relationship between the wavelet scale and the sampling period. We present a mathematical proof in the Appendix.

We can also see that the estimations for the Hurst exponent are not stable under changes in the interval length and estimation scale, while the estimations along the secondary diagonal direction of the tables are still stable. We conclude that the massive hydrologic data have a multiscale consistent property, while the simulation data do not. This indicates that there may exist multifractality in the records.

##### 4.2. Multiscale Diagram for Hydrological Data

In this section, we use the multiscale diagram to test the multifractality that may exist in the records. The multiscale diagram is the graph of against , together with the confidence intervals related to . The multifractality can be detected by measuring the deviation of from the linear form [26].

We use a linear multiscale diagram to examine the alignment in the multiscale diagram, in which is drawn against . The multifractality can be determined by examining whether there is a horizontal alignment in the linear multiscale diagram. The horizontal alignment indicates that the Hurst parameter can control all of the scaling property. The value of can be received by estimating the level. The linear multiscale diagram and multiscale diagram are statistically equivalent.

To determine the alignment of hydrological data, we impose the synthesized FBM () known as the “linear form” for comparison. The regressions for estimating for each of the 10 different are performed for small scales: according to the power-law behaviour of discussed in Section 3.2. For Site 1–Site 10, as well as the FBM, the multiscale diagrams are shown in Figure 4, and the values of for all 10 sites are also shown in Table 7 (the sampling interval is 6 minutes). We choose db3 as the wavelet of the wavelet transform, which is orthogonal to quadratic trends.

**(a)**

**(b)**

In Figure 4, the multiscale diagram in (a) shows a lack of alignment for the hydrologic data compared to the clear alignment of FBM for its self-similarity, which indicates that it has multifractal behaviour. Thus, the multifractal model may be relevant to the hydrological data. The linear multiscale diagram in (b) confirms this conclusion. The Hurst parameter of the synthesized FBM series can be estimated by the horizontal alignment, and its estimator is approximately 0.8. This is not the case for the hydrologic data. In summary, the hydrological data have multifractal behaviour and do not comply with the equation . The errors bars in Figure 4 are based on the weighted linear least square fit, which is used to estimate . They are calculated as . Under the assumption of normality, this reveals a 95% confidence interval. The detailed values of for other sites can be seen in Table 7.

##### 4.3. Multifractal Analysis and Results

In this section, we analyze the daily water level records at the 10 observation sites using MF-DFA analysis for the multifractal behaviour. In order to ensure that there are enough fluctuations in each segment so that MF-DFA can be applied, the high-frequency collected daily records are used. For sufficient samples, we get 100 sets of daily records randomly at different times of one day. The number of records of each sets is about six hundred. We also apply wavelet leaders for multifractal analysis (WLMF) based on the daily water level records for a comparison with MF-DFA.

*Statistical Tests.* It is well known that finite-size effects and lack of statistical convergence will introduce biases on the estimate of when applying the method of MF-DFA and WLMF. Also the robustness of the results obtained by different orders of polynomials in MF-DFA and different vanishing moments of wavelet in WLMF should be tested. In this part, we perform some statistical tests for these issues.

Following the test method reported by Khalil et al. [42], we consider the problem of statistical convergence of through the statistical convergence of and in the two methods. Recall that and can be transformed into a continuous integral over , means the modulus of or . means the probability density function (pdf) of . means the probability density function (pdf) of . The well-defined shape of the pdf’s plots indicates a proper convergence [42, 43]. Figures 5(b), 5(f), 6(b), and 6(f) show and computed at three different scales at Site 4 and Site 7. It can be seen that there are not enough statistics to give a well-defined shape of the pdf. However, the statistical convergence of and does not need so many statistics that can give a well-defined shape of the pdf. Then other statistical tests of convergence are conducted. We need to check the plots of or versus the number of datasets used (from 1 to 100, implicating or ). If or in plots is convergent with the increase of number of datasets used, the statistical convergence of or is confirmed. Figures 5(b), 5(f), 6(b), and 6(f) also show that the shapes of pdf are close to 0 with the decrease of scale. So we investigate the convergence for with a larger scale and with a smaller scale (see Figures 5(c), 5(d), 5(g), 5(h), 6(c), 6(d), 6(g), and 6(h) for Site 4 and Site 7). All the convergence plots except for Site 7 show that and are convergent for the range of from −10 to 10, the range of from to , and the range of from 3 to 6. From Figures 6(d) and 6(h), it can be seen that and at Site 7 are not convergent for with some small scales. This can be confirmed with the log-log plots in Figures 6(a) and 6(e).

Figures 5(a), 5(e), 6(a), and 6(e) show the log-log plots at Site 4 and Site 7 obtained by three different orders of polynomials in MF-DFA and three different vanishing moments of wavelet in WLMF. It can be seen that three different orders of polynomials in MF-DFA give the consistent results, as three different vanishing moments in WLMF do. The robustness of the results obtained by different orders of polynomials in MF-DFA and different vanishing moments of wavelet in WLMF is confirmed.

*Selection of Parameters.* Before comparing MF-DFA and WLMF, we should select parameters to be used in the two methods. In MF-DFA, for detrending, following [16], we use a third-order polynomial to calculate the best polynomial fit of the profile in each segment, which can eliminate the second-order polynomial trend in the original series and is efficient for the hydrologic data in [16]. According to Section 3.3, the scales should be lower than , and the values of should be both positive and negative. For all of the sites, the number of daily water level records is approximately 600, and the range of time scale is selected from to . The upper scale is selected for sufficient statistics. The lower scale is selected for avoiding abnormal points caused by small scales (see Figure 7 for details). According to the statistical tests above and following [14, 16], we choose the range . In WLMF, compared with the results of MF-DFA, we appropriately choose 3 as the number of the vanishing moment, which can eliminate the second-order polynomial trend in the original series. The scale is chosen to be from 3 to 6.

*Multifractal Model.* After using MF-DFA and WLMF, we apply the least-squares method to fit the generalized binomial multifractal model (see (16)). Some other multifractal models may be available such as ansatz model with the functional formwith the three parameters , , and , which has been used successfully to describe the multifractal behaviour of rainfall and runoff records for [16, 44, 45]. However, the ansatz cannot be used for . The two-parameter binomial model is suitable for both positive and negative values of and has fewer parameters. So the two-parameter binomial model is relevant.

*Computational Results.* By using the results of the fit of the generalized binomial multifractal model, we calculate via (17). We can also calculate the fluctuation exponent with MF-DFA; see Table 8 for details. can be used to determine the long-range dependence because of its relationship with the correlation exponent . For a stationary time series, is the well-known Hurst exponent.

Figure 7 shows the log-scaling plots of the water level dynamics of all of the sites. It can be seen that except at Site 7, increases approximately with for and can be fitted well by a straight line. can be estimated by the slope of the regression line. Moreover, MF-DFA fails to estimate for at Site 7.

Figure 8 shows the comparison of the shape of assessed, respectively, by MF-DFA and WLMF, and the fit of the two-parameter binomial model. From the comparison of the shape of , although the values of estimated by the two methods are different, the shapes assessed by MF-DFA and WLMF at all sites (except Site 7) are similar and suggest the two-parameter binomial model. Especially at Site 1, Site 4, Site 6, and Site 10, the shapes of assessed by MF-DFA and WLMF can be consistent by up- and downtranslation. This indicates that the multifractality of the water level records assessed by MF-DFA is credible. Since the MF-DFA and WLMF fail to estimate for at Site 7, the fitting results at Site 7 are based on values of when . From Figure 8, it can be also seen that there are large errors bars for at Site 2, which also indicates poor estimation of for .

To further confirm the results of MF-DFA, we apply MF-DFA to the daily records at different times (morning, afternoon, and night) of a day. As seen in Figure 9, the shapes of at different times are similar. Then, through the mutual verification of the results at different times, we can further confirm the results of MF-DFA that are shown in Figure 8.

From Figure 8, we can also find some details of the results. As mentioned above, the estimation of for by MF-DFA and WLMF at Site 2 and Site 7 is poor, especially at Site 7. The adjacent sites are also affected by these two sites. It can be seen that at Site 6 and Site 8 near 0 deviate from their respective fitting curves. Excluding these sites, of most other sites can be fitted well with the two parameters (, ). The records of most sites satisfy the two-parameter binomial model (or satisfy (16)). For the interpretation of the special behaviour at Site 2 and Site 7, we note that the frequent water storage and discharge activities at these two sites lead to a different water level property when compared to other sites. Moreover, apart from Site 2, Site 7, and those that are seriously affected by these two sites (i.e., Site 1, Site 6, and Site 8), the two parameters and exhibit a relationship such that , whereas in the standard binomial multifractal model [24], the relationship is close to . However, for all of the sites downstream of Site 7, the sum of and is smaller than . The decrease in can also be seen in the figures for WLMF. From Figure 2, we think that the four tributaries downstream of Site 7 have a significant impact on the river; that is, they change the multifractality of this river’s water level.

From Table 8, we can see that all of the Hurst exponents are larger than . This indicates that the water level records of the river have the long-range dependence property, which is the same as the result given by the wavelet-based method. The values of are larger than , which indicates that there is strong multifractality at all sites.

#### 5. Discussion

The long-range dependence property and multifractality we detected in the water level records are identical to the previous results of hydrologic data (see, e.g., [14–16]), and our results are based on multiple observation sites and a high sampling frequency that we have not yet seen in other studies.

In the previous results (see, e.g., [14–16]), we would like to note that the lateral inflows have an effect on the runoffs (or water level). Our results of the change of seem to suggest that the contribution of lateral inflow (tributaries) leads to the change in the multifractality of the trunk stream’s water level.

In light of the results obtained by MF-DFA, the values of are between 0.12 and 0.40 for all of the sites, and the values of are between 0.54 and 0.70. From the previous results of 41 rivers’ runoffs studied by Koscielny-Bunde et al. [14], the values of and of the Würm river are in the range of our values. We would like to note that the Würm river is located in southern Germany (near Munich), which has a temperate continental climate. The river is also within a localized area, and its climate and area are similar to those of the analyzed river. The consistency of the two rivers in parameters validate the reliability of the obtained multifractal model.

As previously mentioned, the spurious multifractality may be induced by intricate nonstationarities, which often exist in the hydrologic data, for example, due to the seasonal cycle or a change in climate. We consider that the obtained multifractality is reliable. First the method we used are based on wavelet analysis and MF-DFA which can eliminate some trends and avoid spurious detection of correlations. Second the multifractality is confirmed consistently via multiple methods (multiscale diagram, MF-DFA, and WLMF). The MF-DFA and WLMF also give the same multifractal model, and this multifractal model is robust in different location and different time. Last the results obtained from the analyzed river are consistent with the previous study of actual example.

#### 6. Conclusions

In this paper, we proceed in an orderly manner to study the LRD and multifractal properties of high-frequency water level records of a river in northern China by using wavelet analysis and MF-DFA.

We observe the long-range dependence property and multiscale behaviour in the time series, which indicates the multifractality. This multifractality is detected via multiscale diagram. Then we focus on the multifractality of records using MF-DFA and WLMF. At most sites, the generalized Hurst exponent can be described with the function and the multifractal models are robust for different times (morning, afternoon, and evening). The validity of multifractality is (partially) supported with the same model obtained from the different methods (MF-DFA and WLMF) and by comparing previous study of actual example. The results also show that the poor estimation of for at Site 2 and Site 7 and the decreased change in indicates human activities and tributaries have a great influence on the multifractality of the river.

We obtain and list the values of , , , and for all available records, including some abnormal records. At all available sites, MF-DFA yields values of are larger than 0.5, whereas WLMF yields values are larger than , which shows a strong multifractality in the data.

As a technical improvement in the estimation of the Hurst parameter, we show that the sampling intervals and choice of scales of the discrete wavelet transformation should be consistent (Tables 4–6).

In contrast to the previous studies, our studies are based on a short sampling time, which we have not seen before. This study contributes a new insight into previous research applications, and it may also be the extension of the study of multifractality. Further investigations are needed, especially for the massive data cases, for example, comparing the fitting of that is obtained from the records partitioned for each year. As mentioned above, the human activities (frequent water storage and discharge activities) at Site 2 and Site 7 lead to poor estimation of for and the change in . This is the key point for the future work on detection and prediction. A real time multifractal model may be built for , , and . The detection and prediction of human activities and climate change in real time may be conducted by the change of , , and . Besides, it would be also interesting to consider the water level dynamics as a multifractional process, study at the different sites with different time lags, and provide predictive models of the water level records.

#### Appendix

#### A. Proof of the Consistent Relationship between the Wavelet Scale and the Sampling Period in Section 4.1

As we know that the process with sampling period is equal to the process with sampling period , therefore process with sampling period series is equal to the process with sampling period . For convenience, let .

Assuming is the Hurst exponent of time series for scale , then Let ; we have Furthermore, According to the estimation method, Under the discrete condition, the degree of the approximation of is different for the different time sequence by the same scale of sampling frequency. Then we only have More, .

#### B. The Related Materials Download Addresses

The multiscale diagram is implemented in Matlab by Veitch and available from his homepage at http://crin.eng.uts.edu.au/~darryl/index.html.

The multifractal detrended fluctuation analysis (MF-DFA) is implemented in Matlab by Ihlen and available at http://www.ntnu.edu/inm/geri/software.

The wavelet leaders for multifractal analysis (WLMF) is implemented in Matlab by Wendt and available from his homepage at https://www.irit.fr/~Herwig.Wendt/index.html.

#### Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

#### Acknowledgments

This work is supported in part by the National Basic Research Program of China (973 Program, Grant no. 2013CB910200), National Science Foundation of China (no. 61103136), the open project of the Hubei Province Key Laboratory of Intelligent Robot (no. 201007), and the Fundamental Research Funds for the Central Universities (no. JBK170164, no. JBK120509, and no. JBK140507).