A Chaos Synchronization Diagnostic: Difference Time Series Peak Complexity (DTSPC)

1. Introduction

Certain dynamical systems have been known since Poincaré to be chaotic [1,2], i.e. they are highly sensitive to initial conditions and parameters. Classical models of deterministic chaos include the Lorenz system modeling weather [3] and the logistic map modeling population growth [4]. Dynamical systems can also display synchronization, when two or more different systems exhibit similar dynamics, either behaving identically as a limiting example, or more broadly possessing similarities in some key aspects, such as sharing the same frequencies or a constant phase lag [5]. Remarkably, though chaotic systems exhibit completely different behaviors with only slightly different initial conditions, adding a small amount of coupling can cause them to be synchronized, since periodicity is not required for synchronization [6]. More broadly, chaos synchronization occurs when two (or many) chaotic systems (either equivalent or nonequivalent) adjust a given property of their motion to a common behavior due to a coupling or to a forcing (periodical or noisy) [7]. Chaos synchronization has found applications in a wide variety of fields such as electrical circuits [8,9,10,11], power systems [12,13,14], medicine [15,16], and chemical reactors [17,18].

Historically, chaos synchronization has been classified into different categories including complete synchronization, practical synchronization, and phase synchronization[19]. Such categorization, however, only provides qualitative classification of chaos synchronization, so research in quantitative complexity measurement of chaos synchronization started to develop over the last two decades. Classic methods of complexity measurement include phase portraits, bifurcation diagram [20], and Lyapunov exponent [21] which require access to the dynamical equations. A growing body of work uses broader and more abstract techniques relying on entropy analysis such using Kolmogorov entropy to quantify the level of unpredictability in a chaotic system over time [22], spectral entropy to diagnose multi-stable fractional-order chaotic systems [23], and entanglement and relative entropy to investigate the chaotic behaviors in quantum systems [24]. One advantage of these entropy-based analyses is that it can often be applied to systems where only the dynamical data is available, in the absence of the actual equations, such as in experimental situations, and the results are typically robust to the presence of noise in the data or the dynamics.

In this paper we explore the applicability of entropy to quantify chaos synchronization broadly. As we discuss below, our method uses `peak-to-peak’ dynamics which is present in almost all basic models of chaos, including the Lorenz system, the Duffing oscillator, the Chua cirsuit, and the Rössler attractor [25]. Specifically, we use the paradigmatic model of the basic coupled Lorenz systems, applying entropy analysis on ordinal patterns extracted from peak dynamics. We show that this entropy-based analysis in quantifying chaos synchronization not only captures overall similarities between two synchronized systems but also reveals a diversity of possible synchronization, non-monotonic changes in behavior with respect to linear change in parameters during transition from one regime to another, and complicated boundaries between different synchronization regimes.

Below, in Section 2, coupled identical Lorenz systems with dissimilar initial conditions are introduced, as well as coupled dissimilar Lorenz systems with identical initial conditions, the synchronization effects of which vary with the change in coupling parameters. Section 3 describes the specific methodology of the new synchronization diagnostic, focusing on converting the raw time-series information into ordinal patterns about the peak-to-peak dynamics before computing the Shannon entropy of the ordinal pattern populations obtained from the processed data. Results and findings about the coupled Lorenz systems after applying the new diagnostic are presented in Section 4 and are compared with the conventional methods, showing values of our new synchronization diagnostic. Finally, discussions and conclusions are presented in Section 5, revealing findings on other interesting behaviors and pointing out directions of further research.

2. Comparing Synchronization in Coupled Identical and Dissimilar Lorenz Systems

In this section we discuss various details about different kinds of synchronization behaviors with the example of coupled identical and dissimilar Lorenz systems, initially studied in weather forecasting. Under certain initial conditions and parameter values, the Lorenz system exhibits chaotic behavior. Before we move on to the more interesting dynamics in the coupled dissimilar Lorenz systems, we first examine the coupled identical systems.

2.1. Identical Lorenz Systems with Different Initial Conditions

Consider two Lorenz systems defined by

$$\begin{array}{cc}\hfill \dot{x_1}& ={\sigma }_1(x_2-x_1)\hfill \\ \hfill \dot{x_2}& =-x_1{x}_3+{\rho }_1{x}_1-x_2\hfill \\ \hfill \dot{x_3}& =x_1{x}_2-{\beta }_1{x}_3\hfill \end{array}$$

(1)

and

$$\begin{array}{cc}\hfill \dot{y_1}& ={\sigma }_2(y_2-y_1)\hfill \\ \hfill \dot{y_2}& =-y_1{y}_3+{\rho }_2{y}_1-y_2\hfill \\ \hfill \dot{y_3}& =y_1{y}_2-{\beta }_2{y}_3\hfill \end{array}$$

(2)

where ${\sigma }_1={\sigma }_2=10,{\rho }_1={\rho }_2=50,{\beta }_1={\beta }_2=8/3$. The initial conditions for the first Lorenz are $x_1=0,x_2=1,x_3=0$; the initial conditions for the second Lorenz are $y_1=1,y_2=1,y_3=1$.

Figure 1a,c show the 3D phase space diagrams of the two identical systems with different initial conditions. Although the two Lorenz systems have the same attractor, Figure 1e showing $x_1$ against $y_1$ makes it clear that the specific trajectories are quite different, demonstrated in time series data as well in Figure 1b,d and the difference time series of $x_1-y_1$ (Figure 1f), which behaves rather randomly. This sets our baseline to observe how coupling alters these differences.

We use the linear coupling model [7]

$$\begin{array}{cc}\hfill \dot{x_1}& ={\sigma }_1(x_2-x_1)+d_1(y_1-x_1)\hfill \\ \hfill \dot{x_2}& =-x_1{x}_3+{\rho }_1{x}_1-x_2+d_2(y_2-x_2)\hfill \\ \hfill \dot{x_3}& =x_1{x}_2-{\beta }_1{x}_3+d_3(x_1-y_1)\hfill \\ \hfill \dot{y_1}& ={\sigma }_2(y_2-y_1)+d_1(x_1-y_1)\hfill \\ \hfill \dot{y_2}& =-y_1{y}_3+{\rho }_2{y}_1-y_2+d_2(x_2-y_2)\hfill \\ \hfill \dot{y_3}& =y_1{y}_2-{\beta }_2{y}_3+d_3(x_3-y_3)\hfill \end{array}$$

(3)

where $d_1,d_2,d_3$ are coupling strengths, which indeed yields that the trajectories of the two chaotic Lorenz systems do overlap. The parameters and initial conditions are the same as above, which are ${\sigma }_1={\sigma }_2=10$, ${\rho }_1={\rho }_2=50$, ${\beta }_1={\beta }_2=8/3$, $x_1=0,x_2=1,x_3=0$ and $y_1=1,y_2=1,y_3=1$.

While Figure 2a,c look quite similar to the uncoupled systems (Figure 2a,c), the relationship between the two systems is now quite different, as visualized in the $x_1$ vs $y_1$ phase portrait. In Figure 2e, the values of $x_1$ and $y_1$ are nearly identical, suggesting that the trajectories of the two systems almost completely overlap. Moreover, as shown in Figure 2f, the difference between $x_1$ and $y_1$ gradually reduce to zero over time. In fact, $x_2-y_2$ and $x_3-y_3$ are also decreasing down to zero. This proves that the trajectories of the two coupled Lorenz systems indeed gradually overlap, so complete synchronization between two identical chaotic systems with different initial conditions is achieved through linear coupling. Note that though the systems synchronize, each of them separately remain chaotic, as can be seen from the attractors in Figure 2a,c. The situation, however, gets far more complicated in the case of coupled dissimilar Lorenz systems which possess identical initial conditions yet dissimilar parameters.

2.2. Dissimilar Lorenz Systems with Identical Initial Conditions

Above, we saw examples of complete chaos and complete synchronization after coupling. However, such clear distinctions do not completely cover the range of behaviors as we show below; even a simple division into practical synchronization or phase synchronization does not suffice to describe the complex behavior seen between complete chaos and complete synchronization. We can see this when we consider two slightly different Lorenz systems, one with parameters ${\sigma }_1=10,{\rho }_1=50,{\beta }_1=8/3$, the other with parameters ${\sigma }_2=10,{\rho }_2=40,{\beta }_2=8/3$. The initial conditions are still $x_1=0,x_2=1,x_3=0,y_1=1,y_2=1,y_3=1$. Each of the two systems are chaotic, and they are unsynchronized. On introducing linear coupling we see a wide range of different synchronization behaviors with more variety as shown in Figure 3 compared to the coupled identical Lorenz systems.

In particular we see in Figure 3 that as the coupling increases the range of $x_1-y_1$ decreases from roughly $\pm 4$ to roughly $\pm 2$. Notice in comparing Figs.Figure 3a–d with Figure 3e,f (at relatively different coupling strengths) we see that $x_1-y_1$ behaves similarly in phase space but displays more variety in the time series. Specifically notice that there are peaks and dips above and below the horizontal axis in all three cases, but the last two cases also exhibit ringing patterns (corresponding to entrainment at certain differences). Moreover, notice that as the coupling strengths increase from the second case to the third case, these ringing patterns change character – specifically oscillating for longer durations around locations further from the horizontal axis (zero difference). This change in ringing patterns as a function of coupling strengths is what distinguishes behaviors where phase-portraits could be a misleading indicator. We now quantify the behavior of the difference time series to more rigorously identify different synchronization patterns across changing coupling strengths, using an ordinal pattern analysis and a computed permutation entropy related to the complexity of the appearance of peaks in the difference time series which we term the Difference Time Series Peaks Complexity (DTSPC). This quantification allow us to map carefully, via a heat map of the DTSPC vs coupling strengths, different kinds of synchronization and transitions as well as the complex shape of the transition regimes. We now turn to describing the method.

3. Methods: Quantifying Synchronization

Across the various types of synchronization shown in Figure 3, although the variety of time series patterns are quite different, they all have the common element that the difference between $x_1$ and $y_1$ is centered around 0, i.e. the time series is centered on the horizontal axis $x_1-y_1=0$. We focus on behavior relative to this horizontal axis, specifically considering the behavior of the ringing patterns previously discussed in the context of the difference time series for $x_1-y_1=0$: We count the number of peaks between zero-crossings in the difference time series noting either the number of peaks or its negative depending on if $x_1-y_1>0$ between two zero-crossings or if $x_1-y_1<0$ respectively. This compresses the time series data into a list of the observed populations of number of consecutive peaks between zero-crossings in the original difference time series. This can be visualized as the bar chart shown in the third column of Figure 4, where the x-axis denotes the elements in the peaks list (i.e. the positive or negative number of peaks) and the y-axis denotes the frequency. We now present all this information in Figure 4 which adds bar charts for the peaks as described to the previous figure where again, as we look down the rows the coupling strength $d_2$ increases while $d_1$ stays unchanged, and the phase portraits become more linear.

These bar charts also indicate that there are various transitional behaviors in chaos synchronization. When the system is not synchronized, there are no ringing patterns but only single spikes, so there are usually very few instances of more than one consecutive peak, as shown in Figure 4c. As the transition begins, the range of possible values for numbers of consecutive peaks increases, whence the bar charts show a corresponding wide range of values on the x-axis. At the early stage of such transition, single peaks are still the majority, as shown in Figure 4f. We find that here the ringing pattern that exist often cross the horizontal axis such that only the last peak in a ringing pattern counts as a single `peak’ by our algorithm. Then, as the system get more synchronized with the increase in coupling strengths, the ringing patterns diverge from the horizontal axis, so the continuous peaks of all lengths are distributed more evenly, as shown in Figure 4o. All of these observations motivate us to compute the Shannon entropy of the distribution of peak populations, thus developing the Difference Time Series Peak Complexity (DTSPC) metric

$$\zeta =-\sum _{i=1}^n{p}_i{log}_2{p}_i$$

(4)

to quantify the level of synchronization. Our investigations below show that $\zeta$ (DTSPC) indeed captures the variety in levels of synchronization of coupled dissimilar Lorenz systems, as shown in Figure 5 where we keep $d_1$ constant while varying $d_2$. Examining these line-graphs reveals that the level of synchronization varies with the change in coupling strength $d_2$, and comparing across different line graphs with different values for $d_1$ shows that the level of synchronization also changes with respect to $d_1$. Hence, we now turn to describing how $\zeta$ behaves for our two systems, with both $d_1$ and $d_2$ as changing parameters.

4. Results

To summarize the results of applying the DTSPC diagnostic on coupled identical and dissimilar Lorenz systems, we present in both cases a heat map of $\zeta (d_1,d_2)$, i.e. displaying the DTSPC as a function of coupling strength parameters. We find distinct regimes of different levels and types of synchronization via this heat map. The color coding for $\zeta$ ranges from warmer (dark red) to the colder (dark blue) where the warmer the color represents a smaller $\zeta$ and hence a less complex set of data from the synchronization dynamics.

Consider Figure 6 which shows a wide variety of $\zeta$ values thus indicating a wide variety of synchronization behaviors with changes in coupling parameters. In particular, we can identify regions or families across the map based on the change in color patterns in the heat map. Specifically, moving diagonally starting from the bottom left corner, we label the dominantly orange-red region as Region A; this is followed by distinct thin ridge of dark red we label as Transition Belt B. Beyond this Transition Belt B is Region C, where the color steadily changes from yellow to bluish-green. We note immediately that we can visually identify a speckle pattern suggesting that the transition need not be smooth as a function of coupling. We also see visually that these transitions are less sensitive to coupling (are slower) in Region A and C while they are more rapid in Transition Belt B.

While it is important to emphasize that we can construct this sort of heat-map given just the peak data without access to either the differential equations or the time-series themselves, to calibrate our ideas about what these various regimes represent, in Figure 7 we explore details of these different synchronization behaviors, showing examples from each region of the heat map. The first example is when $d_1=1.4,d_2=0.8$ in Region A, shown in Figure 7a-c. Region A represents unsynchronized chaos, where $x_1$ and $y_1$ follow completely random relationship as shown in the phase portrait (Figure 7b). Similarly, in the time series (Figure 7b) the difference between $x_1$ and $y_1$ appear to be random and chaotic. This results in its constant crossing of zeros, with few segments of the time series showing larger numbers of consecutive peaks. Thus, there are mainly 0 or 1 consecutive peaks according to our criterion, indeed as shown in Figure 7c which leads to the low $\zeta$ value. In Transition Belt B dynamics change can change quickly with parameters, as shown in Figure 7 d-f when $d_1=3.0,d_2=1.2$. The phase portrait (Figure 7d) shows that the relationship between $x_1$ and $y_1$ becomes mainly linear, and in time series (Figure 7e), accordingly, the difference between $x_1$ and $y_1$ gradually quiets down to zero, where there remains some bursts at the beginning in the first 20 time steps. This time series frequently crosses the horizontal axes such that the count for consecutive peaks remains even more concentrated at 0 or 1 than in Region A, such that the DTSPC $\zeta$ for systems in the Transition Belt B is even lower. The third example with $d_1=5.0,d_2=5.0$ is chosen from Region C and represents complete synchronization. In the phase portrait (Figure 7g) the relationship between $x_1$ and $y_1$ is more strictly linear, and the difference time series (Figure 7h) dies down to zero abruptly almost right after the two Lorenz systems get coupled. However, there could be nuanced and extremely frequent ringing patterns above or below the zero line (red horizontal line), as the number of consecutive peaks can range from $-600$ to 600 as show in Figure 7i. Hence, this analysis reveals overall that as a function of coupling parameters, the dynamical synchronization as quantified by $\zeta$ starts low, and then with the increase in coupling strengths it experiences first a decrease then an increase, which matches the general trend of the heat map (Figure 6). In summary we see that coupling in general increases synchronization however, the transition can be abrupt and non-monotonic in this case, and that the DSTPC captures this appropriately.

In Figure 8 we show a heat map for the coupled dissimilar Lorenz systems. This heat map is over the same range of coupling strengths as in Figure 6 above for coupled identical Lorenz systems but shows only one region of warm coloring in this one, indicating that the information complexity of the peaks distribution is in general less than the identical case above. In this case, we see that the heat map also starts from orange in the bottom-left corner indicated as Region $A^{\prime }$, transitions through the dark red Belt $B^{\prime }$, beyond which is Region $C^{\prime }$, where the entropy changes gradually. Beyond this Region $C^{\prime }$ is the orange-yellow Transition Belt $D^{\prime }$ separating it from Region $E^{\prime }$ above, which is bright yellow.

As before, we calibrate our understanding of this peak data information in Figure 9 where we see an example from each of the five regions in the heat map. As in the identical systems, Region $A^{\prime }$ represents chaotic dynamics for any measure of difference between the two systems, and we can see this for example in Figure 9b displaying the $(x_1,y_1)$ phase-space trajectories which wander chaotically between $\pm 30$. Because of this chaos, the observable time series frequently crossing zero and observations of number of consecutive peaks between zero-crossings are therefore mainly centered at 0 and 1, as shown in Figure 9c resulting in a low $\zeta$. In Transition Belt $B^{\prime }$ signals a rapid onset of synchronization as a function of coupling. In the phase-portrait (Figure 9d) we see that the in-phase component of the trajectories – that is, trajectory lines being traced from the bottom-left to top-right diagonal – increases, while the out-of-phase component – that is, trajectory lines being traced from the bottom-right to top-left – decreases. This also means that $x_1$ and $y_1$ are more directly correlated, and in general, becoming more linear. We see this in the difference time series (Figure 9e) which abruptly decreases to lie between $\pm 20$ most of the time, with the observable more symmetric around zero. This leads to the data in the bar chart (Figure 9f) becoming even more centered at 0 and 1 consecutive peaks, resulting in a decrease in $\zeta$.

In Region $C^{\prime }$ the differences continue to decrease, with the in-phase component of the phase portrait (Figure 9g) continuing to become sharper and linear. In the difference time series (Figure 9h) there start to appear some ringing patterns above and below the horizontal axis, though not yet very obvious. The difference between $x_1$ and $y_1$ further decrease to roughly between $\pm 10$. Transition Belt $D^{\prime }$ is where practical synchronization becomes rather apparent. As shown clearly in the phase portrait (Figure 9j), $x_1$ and $y_1$ follow a linear relationship. Accordingly, the time series for the difference between these two observables (Figure 9k) reduces to between $\pm 3$, and clear ringing patterns start to appear that lie strictly above and below the x-axis. However before this clear separation above and below the x-axis we also see some small-scale dynamics in this difference variable which presents as `signal fuzziness’ and a `stickiness’ when this difference crosses zero resulting from the slight entanglement of the two oscillators’ trajectories when they are exchanging positions in phase space during synchronization, and often resulting in a small single peak for our criterion. This results in many 0 and 1 consecutive peaks (Figure 9l), so the $\zeta$ remains low. Finally, Region $E^{\prime }$ is where the two systems are more synchronized, as the ringing patterns appear further apart from each other, shown in Figure 9n, with fewer crossings at zero. This indicates that there is less stagnation during a phase change, so the number of 0 and 1 consecutive peaks reduces. Thus, the data is more spread out in the bar chart, as shown in Figure 9o and $\zeta$ consequently experiences a significant increase. Thus, to summarize, $\zeta$ starting off low in Region $A^{\prime }$ then decreases and increases with the increase in coupling strengths as before for identical systems. However, the transition to synchronization now carries more regions of distinct dynamical regimes: Transition Belt $B^{\prime }$, Region $C^{\prime }$, and Transition Belt $D^{\prime }$. Within them, $\zeta$ first decreases a little, then increases, and slightly decreases again, each part corresponding to different dynamics of the system as we have just seen. This shows again the value of using $\zeta$ to quantitatively identify and classify synchronization behaviors.

We also note that both heat maps (as shown in Figure 6 and Figure 8) show that the synchronization is not locally smooth as a function of coupling strength, given the speckle pattern or occasional grids visibly brighter or darker than their surroundings. Even when we zoom in to some of these areas, such nuanced patterns still exist, as shown in Figure 10 below. While this sensitive dependence on the parameter of the dynamical complexity of the synchronization is perhaps not surprising, it is useful to see how our DTSPC technique based on Shannon entropy can capture this sensitivity.

5. Discussion and Conclusions

In this paper we computed the Shannon entropy of peak population distributions in the chaotic coupled Lorenz systems with peak-to-peak dynamics, introducing a Difference Time Series Peaks Complexity (DTSPC) diagnostic for chaos synchronization, which digitizes a characteristic ringing pattern that changes as a function of parameter in the dynamics of the difference variable for the paradigmatic cases of coupled identical and dissimilar Lorenz systems. We find that the DTSPC entropic analysis of ordinalized data indeed compresses useful information in being able to distinguish and quantify via a single metric the differences in system behavior during the transition from chaos to synchronization. We are thus able to more systematically study the its variation across certain set of parameter values and initial conditions, hence considering both the range of behaviors and the parametric dependence via a single heat map. In particular, in this test case, we can distinguish a variety of synchronization behaviors, clarifying that though system behaviors can be macroscopically classified into non-synchronization (or chaos) or complete synchronization, such rough classification omits detailed behaviors including non-monotonic parameter dependence. We again emphasize that these behaviors are not visible when examining phase portraits and can be constructed using substantially less access to to the complete time-series; these results are also robust to minor noisy perturbations.

In conclusion, the success of our Difference Time Series Peaks Complexity (DTSPC) diagnostic in creating digitized representation of diverse synchronization behaviors from phase space diagrams and time series demonstrates the effectiveness of using Shannon entropy to quantify chaos in systems with peak-to-peak dynamics, specifically utilizing the statistical information behind peak distribution. For further development of this research, one can test and evaluate the applicability of entropy and its variations on chaotic systems without peak-to-peak dynamics or of higher dimensions. Additionally, the complexity of zoomed-in heat maps shown in Figure 10 suggests a complex boundary between different behaviors as a function of parameter; understanding how these boundaries are shaped and the possible existence of such patterns across other systems during chaos synchronization are both potential directions for further research.