Discriminating chaotic and stochastic time series using permutation entropy and artificial neural networks.

Extracting relevant properties of empirical signals generated by nonlinear, stochastic, and high-dimensional systems is a challenge of complex systems research. Open questions are how to differentiate chaotic signals from stochastic ones, and how to quantify nonlinear and/or high-order temporal correlations. Here we propose a new technique to reliably address both problems. Our approach follows two steps: first, we train an artificial neural network (ANN) with flicker (colored) noise to predict the value of the parameter, [Formula: see text], that determines the strength of the correlation of the noise. To predict [Formula: see text] the ANN input features are a set of probabilities that are extracted from the time series by using symbolic ordinal analysis. Then, we input to the trained ANN the probabilities extracted from the time series of interest, and analyze the ANN output. We find that the [Formula: see text] value returned by the ANN is informative of the temporal correlations present in the time series. To distinguish between stochastic and chaotic signals, we exploit the fact that the difference between the permutation entropy (PE) of a given time series and the PE of flicker noise with the same [Formula: see text] parameter is small when the time series is stochastic, but it is large when the time series is chaotic. We validate our technique by analysing synthetic and empirical time series whose nature is well established. We also demonstrate the robustness of our approach with respect to the length of the time series and to the level of noise. We expect that our algorithm, which is freely available, will be very useful to the community.

Introduction

Chaotic and stochastic systems have been extensively studied and the fundamental difference between them is well known: in a chaotic system an initial condition always leads to the same final state, following a fixed rule, while in a stochastic system, an initial condition leads to a variety of possible final states, drawn from a probability distribution[1]. However, the signals generated by chaotic and stochastic systems are not always easy to distinguish and many methods have been proposed to differentiate chaotic and stochastic time series[2-10].

A related important problem is how to appropriately quantify the strength and length of the temporal correlations present in a time series[11-14]. The performance of these methods varies with the characteristics of the time series. As far as we know, no method works well with all data types, because methods have different limitations, in terms of the length of the time series, the level of noise, the stationarity or seasonality of the underlying process, the presence of linear or nonlinear correlations, etc. Moreover, any time series analysis method will return, at least, one number. Therefore, to obtain interpretable results, the values obtained from the analysis of the time series of interest need to be compared with those obtained from other “reference” time series, where we have previous knowledge of the underlying system that generates the data. Here we use as “reference” model a well-known stochastic process: flicker noise (FN).

A FN time series is characterized by a power spectrum , with being a quantifier of the correlations present in the signal[15]. Flicker noise has been extensively studied in diverse areas such as electronics[16,17], biology[18,19], physics[20,21], economy[22,23], meteorology[24], astrophysics[25], etc. Furthermore, related to this issue, many methods described in the literature are able to evaluate the time correlation quantification , such as the Hurst exponent [2,11-13,15,21,26].

In this paper, we propose a new methodology that simultaneously allows to distinguish chaotic from stochastic time series, and to quantify the strength of the correlations. Our algorithm, based on an Artificial Neural Network (ANN)[27], is easy to run and freely available[28]. We first train the ANN with flicker noise to predict the value of the parameter that was used to generate the noise. The input features to the ANN are probabilities extracted from the FN time series using ordinal analysis[29], a symbolic method widely used to identify patterns and nonlinear correlations in complex time series[30-32]. Each sequence of D data points (consecutive, or with a certain lag between them), is converted into a sequence of D relative values (smallest to largest), which defines an ordinal pattern. Then, the frequencies of occurrence of the different patterns in the time series define the set of ordinal probabilities, which in turn allow to calculate information-theoretic measures such as the permutation entropy (PE, described in ”Methods”). The PE has been extensively used in the literature, due to the fact that is straightforward to calculate, and it is robust to observational noise. Interdisciplinary applications have been discussed in Ref.[33] and, more recently, in a Special Issue[34].

After training the ANN with different FN time series, , generated with different values of , we input to the ANN ordinal probabilities extracted from the time series of interest, x, and analyze the output of the ANN, . We find that is informative of the temporal correlations present in the time series x. Moreover, by comparing the PE values of x and of (a FN time series generated with the value of returned by the ANN), we can differentiate between chaotic and stochastic signals: the PE values of x and are similar when x is mainly stochastic, but they differ when x is mainly deterministic. Therefore, the difference of the two PE values serves as a quantifier to distinguish between chaotic and stochastic signals. We use several datasets to validate this approach. We also analyze its robustness with respect to the length of the time series and noise contamination.

This paper is organized as follows. In the main tex,t we present the results of the analysis of synthetic and empirical time series, which are described in “Datasets”. Typical examples of the time series analyzed are presented in Fig. 1. In ”Methods”, we describe the ordinal method and the implementation of the algorithm, schematically represented in Fig. 2.

Figure 1

Examples of time series analyzed, their probability density functions (PDFs) and power spectral densities (PSDs). (a,b) Time series generated by iteration of the map [Eq. (6) with ] and its PDF. (c,d) Uniformly distributed white noise and its PDF. We see that the PDF of the deterministic map is identical to the PDF of the noise. (e,f) PSD of the Schuster map, Eq. (8) with parameter , and of a Flicker noise with . We note that the PSD of the Schuster map has a long decay that is very similar to a decay of the noise. (g) PDF of m summed logistic maps (Eq. 7 with ), which approaches a Gaussian as m increases. (h–k) Examples of empirical time series analyzed: (h) a stride-to-stride of an adult walking in a slow velocity, interpreted as an stochastic process; (i) daily number of sunspots as a function of time (in years), where its fluctuations are interpreted as stochastic; (k) voltage across the capacitor of an inductor-less Chua electronic circuit, whose oscillations are known to be chaotic.

Figure 2

Schematic representation of the methodology. We compute the probabilities of the ordinal patterns and then use them as input features to the ANN. The ANN returns the temporal correlation coefficient . Then we compare the permutation entropy of the analyzed time series with the permutation entropy of a FN time series generated with the same value. Based on this comparison, we use Eq. (1) to classify the time series as chaotic or stochastic.

Results

Analysis of synthetic datasets

The main result is depicted in Fig. 3. Panel (a) shows the normalized permutation entropy (Eq. 4) vs. the time-correlation coefficient . The filled (empty) symbols correspond to different types of stochastic (chaotic) time series, and the solid black line corresponds to FN time series generated with , which is accurately evaluated by the ANN. For , FN has a uniform power spectrum, characteristic of an uncorrelated signal (white noise), with equal ordinal probabilities and, hence, . Otherwise, for , some ordinal patterns occur in the time series more often than others, and the ordinal probabilities are not all equal, which decreases the permutation entropy. These results are consistent with those that have been obtained by using different methodologies[7,10,26].

Figure 3

Temporal correlation and distinction of chaotic and stochastic synthetic signals. Panel (a) shows the normalized permutation entropy as a function of , evaluated through the ANN, for different time-series (stochastic and chaotic). The black solid line represents FN signals, which are used for training and testing the ANN. Filled symbols represent different stochastic signals, which are very close to the FN curve. Empty symbols represent chaotic signals, which are far away. The distance to the FN curve is measured by (Eq. 1) and it is shown in panel (b) higher values of indicate chaotic time series, while lower, vanishing values indicate stochastic ones.

In Fig. 3, we note that fBm signals have larger time-correlation ( closer to 2, a classic Brownian motion) than the other three stochastic systems . However, their permutation entropies are very close to those of the FN signals. The key observation is that stochastic time series all fall close to the FN curve, while chaotic ones do not, namely, map, Lorenz system, logistic map, skew tent map, and Schuster map. The distance to the FN curve thus serves to distinguish stochastic and chaotic time series. This is quantified bywhere is the permutation entropy of the analyzed time series and is the PE of a flicker noise time series generated with the value of returned by the ANN, . The results are presented in Fig. 3b, where we see that stochastic signals have , and deterministic signals have . To summarize this finding, Table 1 depicts and for ten representative systems.

Table 1

Results obtained from stochastic and deterministic time series: mean and standard deviation of the parameter and of (Eq. 1), calculated from 1000 time series generated with different initial conditions and noise seeds.

Stochastic process	α_e	Ω(α_e)
FN (α=0)	-0.008±0.013	0.00006±0.00005
fBm (H=0.5)	1.74±0.01	0.005±0.001
fGn (H=0.5)	-0.009±0.014	0.00006±0.00004
Cauchy	-0.008±0.004	0.00007±0.00005
Uniform	-0.008±0.003	0.00007±0.00006

Next, we study the applicability of our methodology to noise-contaminated signals. We analyze the signalwhere is a deterministic (chaotic) signal “contaminated” by a uniform white noise, , and controls the stochastic component of . For the signal is fully deterministic (fully stochastic).

Figure 4a shows as a function of for different chaotic signals. As expected, for , is high, but as grows, the level of stochasticity increases and decreases. At , the signal is strongly stochastic, as reflected by . For comparison, in Fig. 4a we also present results obtained by shuffling a chaotic time series. As expected, for all because temporal correlations are destroyed by the shuffling process.

Figure 4

Limits for identifying determinism. Panel (a) shows the influence of noise contamination: the quantifier , Eq. (1) is plotted as a function of the level of noise, , in Eq. (2). We see that as increases, decreases, but it remains, for high values of , different from the value obtained from shuffled data (black pentagons). Therefore, small values of can still reveal determinism in noise-contaminated signals. Panel (b) shows the effect of adding several independent chaotic signals: is plotted as a function of the number m of signals added. We see that as m increases, decreases, indicating that the deterministic nature of the signal can no longer be detected.

We expect that the addition of a sufficiently large number of independent chaotic signals gives a signal that is indistinguishable from a fully stochastic one. This is verified in Fig. 4b, where the horizontal axis represents the number, m, of independent chaotic signals added. Here a high value of is observed for (a single chaotic signal), but as m increases since the chaotic nature of added signals is no longer captured (examples of the PDFs of the time series obtained from the addition of m Logistic maps were presented in Fig. 1g, where we can observe a clear evolution towards a Gaussian shape).

To further explore the robustness of our methodology, we investigate the role of the length N of the analyzed time series in the evaluation of the quantifier (Eq. 1). Figure 5 shows as a function of N, where panel (a) depicts stochastic signals, and (b) chaotic ones. We see that even for , for all stochastic signals in panel (a) , which indicates that we can identify the stochastic nature of short signals. For the chaotic signals in panel (b), for (except for map), and for , for all signals, which demonstrates that our method can also detect determinism in short signals.

Figure 5

Robustness with respect to the time series length. We evaluate as a function of the time series length N for stochastic signals (a) and chaotic ones (b). In (a), we observe that even for , , which confirms the stochastic nature of the signal. In (b), even for , (except for map). Also, for all cases for , which classifies the signals as chaotic even with only 1000 datapoints. The error bars are the standard deviation over 1000 of samples. (c) Mean absolute error, , between the output of the ANN (), and the parameter used to generate the time series of flicker noise (depicted in color-code) as a function of the length of the time series, N.

As discussed in ”Methods”, the ANN was trained with flicker noise signals with data points. However, it is interesting to analyze how much data the trained ANN needs, in order to correctly predict the value of a flicker noise time series. To address this point, we generate FN time series and analyze the error of the ANN output, , as a function of the length of the time series, N, and of the value of used to generate the time series. The results are presented in panel (c) of Fig. 5 that displays the mean absolute error, . We see that as N increases, decreases. The error depends on both, and N, and tends to be larger for high due to non-stationarity and finite time sampling[10]. For FN time series longer than datapoints, the ANN returns a very accurate value of .

Analysis of empirical time series

Here we present the analysis of time series recorded under very different experimental conditions, as described in “Datasets”. Figure 6 displays the results in the plane (, ). The values obtained for the Chua circuit data and for the laser data confirm their chaotic nature[35,36] ( and respectively). For the star light intensity , confirming the stochastic nature of the signal[37]. For the number of sunspots, which is a well-known long-memory noisy time series, . In this case the value of obtained () confirms the results of Singh et al.[38] where a Hurst exponent close to 1 was found. Regarding the five time series of RR-intervals of healthy subjects, our algorithm identifies stochasticity () in all of them, which is consistent with findings of Ref.[9].

Figure 6

Analysis of empirical time series. Results obtained from each time series are presented in the plane (, ). Deterministic signals are the Chua circuit data (brown triangle) and the laser data (orange ‘X’ marker) that have . The other signals [the light intensity of a star (yellow dot), the number of sunspots (cyan diamond), the heart variability of healthy subjects (magenta thin diamond), and the different groups of human gait dynamics (green, blue, red, and black squares)] are stochastic and have . For the human gait, the inset depicts the as a function of the age of each subject. Consistent with[40], the correlation coefficient decrease with the ages.

The last empirical set analyzed reveals the nature of the dynamics of human gait: regardless of the age of the subjects, confirming the stochastic behavior discussed in[39]. In the inset we show that the value returned by the ANN decreases with the age, which is also in line with the results presented in[40], obtained with Detrended Fluctuation Analysis (see Fig. 6 of Ref.[40]). The authors interpret this variation as due to an age-related change in stride-to-stride dynamics, where the gait dynamics of young adults (healthy) appears to fluctuate randomly, but with less time-correlation in comparison to young children[40].

Discussion

We have proposed a new time series analysis technique that allows to differentiate stochastic from chaotic signals, and also, to quantify temporal correlations. We have demonstrated the methodology by using synthetic and empirical signals.

Our method is based on locating a time series in a two dimensional plane determined by the permutation entropy and the value of a temporal correlation coefficient, , returned by a machine learning algorithm. In this plane, stochastic signals are very close to a curve defined by Flicker noise, while chaotic signals are located far from this curve. We have used this fact to define a quantifier, , that is the distance to the FN curve. serves to distinguish stochastic and chaotic time series, and it can be used to analyze time series, even when they are very short (with time series of 100 datapoints we found that or , if the time series is stochastic or chaotic respectively, Fig. 5). We also found that small values of can be used to identify underlying determinism in noise-contaminated signals, and in signals that result from the addition of a number of independent chaotic signals (Fig. 4). We have also used our algorithm to analyze six empirical datasets, and obtained results that are consistent with prior knowledge of the data (Fig. 6). Taken together, these results show that the proposed methodology allows answering the questions of how to quantify stochasticity and temporal correlations in a time series.

Our algorithm is fast, easy-to-use, and freely available[28]. Thus, we believe that it will be a valuable tool for the scientific community working on time series analysis. Existing methods have limitations in terms of the characteristics of the data (length of the time series, level of noise, etc.). A limitation of our algorithm lies in the analysis of noise-contaminated periodic signals, because their temporal structure may not be distinguished from the temporal structure of a fully stochastic signal with a large value. Future work will be directed at trying to overcome this limitation. Here, for a “proof-of-concept” demonstration, we have used a well-known machine learning algorithm (a feed-forward ANN), a rather simple training procedure, and a popular entropy measure (the permutation entropy). We have not tried to optimize the performance of the algorithm. We expect that different machine learning algorithms, training procedures, and entropy measures can give different performances, depending on the characteristics of the data analyzed. Therefore, the methodology proposed here has a high degree of flexibility, which can allow to optimize performance for the analysis of particular types of data.

Methods

Ordinal analysis and permutation entropy

Ordinal analysis allows the identification of patterns and nonlinear correlations in complex time series[29]. For each sequence of D data points in the time-series (consecutive, or with a certain lag between them), their values are replaced by their relative amplitudes, ordered from 0 to . For instance, a sequence in the time series transforms into the ordinal pattern “0123”, while transforms into “0312”. As an example, Fig. 7 shows the ordinal patterns formed with consecutive values..

Figure 7

Schematic illustration of 24 ordinal patterns that can be defined from consecutive data values in a time series.

We evaluate the frequency of occurrence of each word, defined as the ordinal probability with , where i represents each possible word. Then, we evaluate the Shannon entropy, known as permutation entropy[29]:The permutation entropy varies from if the j-th state (while ) to if . The normalized permutation entropy used in this work is given by:In this work, to calculate the ordinal patterns, we have used the algorithm proposed by Parlitz and coworkers[32]. We have used and no lag, i.e., values overlap in the definition of two consecutive ordinal patterns. Therefore, we use as features to the ANN (see below) the probabilities of the ordinal patterns. For a robust estimation of these probabilities, a time series of length is needed. However, as we show in Fig. 5, the algorithm returns meaningful values even for time series that are much shorter.

Artificial neural network

Deep learning is part of a broader family of machine learning methods based on artificial neural networks (ANNs)[41]. In this work, we use the deep learning framework Keras[42] to compile and train an ANN. Since we want to regress the information of the features into a real value (classical scalar regression problem[43]) an appropriate option is a feed-forward ANN. The ANN is trained to evaluate the time correlation coefficient considering as features the 720 probabilities of the ordinal patterns. We connect the input layer to a single dense layer with 64 output units connected to a final layer, regressing all the information of the inputs into a real number. Other combinations were tested with different numbers of units () and layers. We found that a single layer with 64 units was sufficient to accurately predict the value. The ANN parameters and the compilation setup are given in Table 2. As explained in the discussion we have used the feed-forward ANN as a simple option for a “proof-of-concept” demonstration. Other deep learning/machine learning methods or a different compilation setup may give better results depending on the type of data that is analyzed.

Table 2

Compilation setup and parameters of the feed forward ANN.

Compilation setup
Optimizer	Adam
Loss function	Mean square error
Metrics	Mean absolute error

The training stage of the ANN is performed using a dataset of 50,000 flicker noise time series with points, where the parameter of each time series is randomly chosen in (see “Datasets” for details). We separate the dataset into two sets: the training dataset (), and the test dataset (). To quantify the error between the output of the ANN and the target, , we use the mean absolute error:where L is the number of samples. The training stage is concluded and then the parameters of the ANN are fixed. After that, we apply the ANN to the test dataset, and the error in the evaluation of regarding is .

Datasets

Stochastic systems

In this paper, we use three types of stochastic signals: flicker noise (FN), fractional Brownian motion (fBm) and fractional Gaussian noise (fGn). Flicker noise (FN) or colored noise time series are used for the training of the Artificial Neural Network. They are generated with the open Python library colorednoise.py[44,45]. fBm and fGn time series are generated with the Python library fmb.py[46]. Both time series depend on a Hurst index [47]. For the fBm corresponds to the classical Brownian motion. If () the time-series is positively (negatively) correlated. For fGn characterizes a white noise[47]; if () the fGn time series exhibits long-memory (short-memory). The Hurst index is related to the of flicker noise: for a fBm stochastic process, and ; for fGn, and [47].

Chaotic systems

In this paper, we analyze time series generated by five chaotic systems:

The generalized Bernoulli chaotic map, also known as map, described by:where controls the dynamical characteristic of the map. Throughout this paper, we use , which leads to a chaotic signal[1].

The well-known logistic map[1]:we use to obtain a chaotic signal[1].

The Schuster map[48], which exhibits intermittent signals with a power spectrum . It is defined as:where we use .

The skew tend map, which is defined asHere, we use in order to obtain a chaotic signal[3].

The also well-known Lorenz system, defined as:with parameters , , and , which lead to a chaotic motion[49]. For analyze the time series of consecutive maxima of the x variable.

Empirical datasets

We test our methodology with empirical datasets recorded from diverse chaotic or stochastic systems. Additional information of the datasets can be found in Table 3. These are:

Table 3

Characteristics of empirical datasets.

Dataset	Number of samples	Mean length 〈N〉	Availability
E-I	1	6000	[51]
E-II	1	750	[52]
E-III	17	1600	[53]
E-IV	3	15,000	[57–59]
E-V	5	100,000	[63]
E-VI	60	800	[65,66]

Dataset E-I: Data recorded from an inductorless Chua’s circuit constructed as[50]. The circuit was set up and the data was kindly sent to us by Vandertone Santos Machado[51]. The voltages across the capacitors depict chaotic oscillations. To detect this chaoticity we compute the maxima values of the first capacitor .

Dataset E-II: Fluctuations in a chaotic laser data approximately described by three coupled nonlinear differential equations[36]. To detect the chaoticity of the laser, we compute the maxima values of the time series. The data is available in[52].

Dataset E-III: Light intensity of a variable dwarf star (PG1599-035)[36] with 17 time-series (segments). These variations may be caused by an intrinsic change in emitted light (superposition of multiple independent spherical harmonics[36]), or by an object partially blocking the brightness as seen from Earth. The fluctuations in the intensity of the star have been observed to result in a noisy signal[36]. The data is open and freely available in[53].

Dataset E-IV: Three time-series of the sunspots numbers for the period of 1976–2013[38], the daily sunspots numbers depicts a noisy “pseudo-sinusoidal” behavior. It is accepted that magnetic cycles in the Sun are generated by a solar dynamo produced through nonlinear interactions between solar plasmas and magnetic fields[54,55]. However, the fluctuations in the period in the cycles is still difficult to understand[56]. This type of data has been analyzed in[38], where its stochastic fluctuations depict a Hurst exponent , meaning the data carries memory. The data can be found at[57-59].

Dataset E-V: Five RR-interval time-series from healthy subjects. Each time series have RR intervals (the signals were recorded using continuous ambulatory electrocardiograms during 24 h). It still a debate if the heart rate variability is chaotic or stochastic[9]. While some studies suggest that heart rate variability is a stochastic process[9,60,61]. Much chaos-detection analysis has been identified as a chaotic signal[9,62]. The dataset is open and freely available in[63].

Dataset E-VI: Fractal dynamics of the human gait as well as the maturation of the gait dynamics. The stride interval variability can exhibit randomly fluctuations with long-range power-law correlations, as observed in[39]. Moreover, this time-correlation tends to decrease in older children[39,40]. The analyzed dataset is then separated into 3 groups, related to the subjects’ ages. Group No. 1 has data for subjects with 3- to 5-years old; Group No. 2 has data for subjects with 6- to 8-years old; Group No. 3 has data for subjects with 10 - to 13-years old; Group No. 4 has data for subjects with 18- to 29-years old[39]. The data is open and freely available in[64-66].