Numerical and experimental study of the effects of noise on the permutation entropy

We consider two dynamical systems: the one-dimensional logistic map and the three-dimensional Rössler system. In both cases, we study the effect of adding to the dynamical equations a Gaussian white noise, ${\xi }_{t}$ with $\langle {\xi }_{t}\rangle =0$ and temporal correlation $\langle {\xi }_{t}{\xi }_{{t}^{\prime }}\rangle ={\delta }_{t,t^{\prime} }$. We considered also the case of observational noise (not shown), where the dynamics is deterministic but the noise affects the observation, obtaining very similar results.

2.1.1. Logistic map

$$\begin{eqnarray}&&{x}_{t+1}=4{x}_{t}(1-{x}_{t})+\alpha {\xi }_{t},\end{eqnarray} \tag{ 1 }$$

where x_t is the state of the system at iteration t and α is the noise strength. In order to constrain the variable x_t in the interval $[0,1]$, the values of ${\xi }_{t}$ that would lead to ${x}_{t+1}\gt 1$ or ${x}_{t+1}\lt 0$ are simply discarded and redrawn. Thus, the noise ${\xi }_{t}$ is temporally uncorrelated, but not purely Gaussian due to this truncation effect. To investigate the variation of the permutation entropy with the noise strength, we computed the permutation entropy, for each value of α, from time series of length $N=1.2\times {10}^{7}$. We have also studied other nonlinear one-dimensional maps (Tent, Bernoulli and Quadratic) and obtained very similar results to those of the logistic map (results not shown).

2.1.2. Rössler system

The Rössler equations read

$$\begin{eqnarray}\dot{X} & = & -Y-Z+\alpha \xi (t),\\ \dot{Y} & = & X+{aY},\\ \dot{Z} & = & b+Z(X-c),\end{eqnarray} \tag{ 2 }$$

where $\{X,Y,Z\}$ are the states of the system at time t, α is the noise strength and $\{a,b,c\}$ are the local parameters set at {$0.1,0.1,18.0$}, respectively.

In order to apply the symbolic methods (ordinal patterns or blocks) we need to discretize the dynamics. Instead of employing temporal sampling [33], we introduce a Poincaré section [34] at X = 0, and analyze the time intervals between consecutive crossings of the Poincaré plane. The reason for this choice is that it is conceptually similar to how we discretize the experimental time series, as discussed in the next subsection. For each value of α, the permutation entropy is computed from time-series of $N=1.2\times {10}^{7}$ data points.

Experimental data was recorded from the output intensity of a semiconductor laser with optical feedback operating in the low-frequency fluctuations (LFFs) regime. In this regime, the laser intensity displays sudden and apparently random dropouts, followed by gradual recoveries. This spiking dynamics has received considerable attention because the intensity dropouts are induced by stochastic effects and deterministic nonlinearities. The optical feedback introduces a delay which renders the system in principle infinite dimensional. Therefore, the laser in the LFF regime generates complex fluctuations that, because of the stochastic and high-dimensional nature of the underlying dynamics, are interesting candidates to be investigated by means of complexity measures such as the permutation and block entropies.

The experimental setup is the same as in [32] and uses a 650 nm AlGaInP semiconductor laser (SONY SLD1137VS) with optical feedback. The feedback was given through a mirror placed 70 cm apart from the laser cavity, with a round trip of 4.7 ns. The feedback was controlled using a neutral density filter that can adjusts the light intensity injected into the laser. The laser has a solitary threshold current of ${I}_{\mathrm{th}}=28.4$ mA. The temperature and current of the laser were stabilized using a combi controller Thorlabs ITC501 with an accuracy of 0.01 C and 0.01 mA, respectively. The current used during the experiment was I = 29.3 mA and the temperature was set at T = 17 C. The neutral density filter was adjusted so that the threshold reduction due to feedback was about 7%. The signal was captured using a photo detector (Thorlabs DET210) connected to a FEMTO HSA-Y-2-40 amplifier and registered with a 1 GHz digital oscilloscope (Agilent Infiniium DSO9104A) with 0.2 ns of sampling. The intensity time series were acquired from the oscilloscope by a LabVIEW program that uses a threshold to detect the times when the intensity drops, and calculates the time intervals between successive threshold crossings (in the following, referred to as inter-dropout-intervals, IDIs). We recorded in this way time series of more than 10⁵ consecutive IDIs.

We compare two different methods to transform a time-series, $x(t)=\{x(1),x(2)$ ...$x(N)\}$, into a sequence of symbols, s(t): ordinal patterns and blocks.

In both cases, one needs to choose a dimension D for defining vectors made up of consecutive entries of the time series, i.e. $\{x(i),x(i+1),...,x(i+D-1)\}$. Ordinal patterns and blocks differ by the way in which the entries of these vectors are transformed into symbols. Ordinal patterns classify them according to the ranking (from the largest to the smallest value) of the D entries in the vectors. The total number of ordinal patterns of length D is then equal to the number of permutations, $D!$. For example, with D = 2 there are two ordinal patterns: $x({t}_{i})\gt x({t}_{i+1})$ corresponding to the ordinal pattern '01' and $x({t}_{i})\lt x({t}_{i+1})$ corresponding to the ordinal pattern '10'.

In cases when, due to finite resolution, equal values can occur, a small observational noise (10⁻⁸) is added to the time-series. The amplitude of the noise is sufficiently small to not modify the ordinal relations in the data set, except for entries of equal value. Equal values are very rare for all the time series we considered.

For blocks, the phase space is first divided into Q regions, associating a symbol to each region. Blocks represent all vectors in which each value of the time series correspond to the same symbol. For example, let us consider the time series $x(t)=\{0.1,0.6,0.7,0.3\}$, and partition the phase space into the two regions $[0,0.5)$ and $[0.5,1]$, associating to them the symbols 0 and 1 respectively. With D = 2, the blocks associated to the time series are $\{01,11,10\}$.

Then the estimation of the permutation entropy [1] or the block entropy [24, 25], is simply the entropy of the frequency p_i of the different patterns in the time series

$$\begin{eqnarray}&&{H}_{D}=-\displaystyle \sum _{i}^{M}{p}_{i}\mathrm{ln}{p}_{i},\end{eqnarray} \tag{ 3 }$$

where M is the number of possible patterns: for permutation entropy, $M=D!$, while for block entropy, $M={Q}^{D}$.

Figure 1(a) displays the permutation entropy, H_D, versus the dimension of the ordinal patterns, D, computed from time series of the logistic map at different noise strengths. It can be observed that H_D increases monotonically with D, regardless of the noise strength. As the noise increases, H_D approaches its maximum value, corresponding to equally probable ordinal patterns, ${H}_{\mathrm{max}}=\mathrm{ln}D!$ (solid black line). Note that at $\alpha =1$ (pentagons) the values of H_D is already very close to H_max. Figure 1(b) displays H_D as a function of the noise strength α. A clear transition from low-noise to high-noise can be observed, for a value of the noise strength approximately independent of D. The difference between the values of the entropies at low and high noise becomes more pronounced as D increases.

Zoom In Zoom Out Reset image size

To further investigate this transition, figure 2(a) displays the difference ${H}_{D}-{H}_{D-1}$ as a function of D, for various values of noise strength. As before, we indicate with a thin black line the noise-dominated limit in which all patterns are equiprobable, ${H}_{D}-{H}_{D-1}=(\mathrm{ln}D!-\mathrm{ln}(D-1)!)$. In the opposite limit of almost-deterministic, as D grows the expected value of ${H}_{D}-{H}_{D-1}$ is the Kolmogorov–Sinai entropy [26], which for a one-dimensional chaotic map is equal to the Lyapunov exponent λ. In the case of logistic map for a local parameter set at 4 one has $\lambda =\mathrm{ln}2$, indicated by the thick black line. As shown in detail in figure 2(c), we identify three possibilities:

• a almost-deterministic regime in which ${H}_{D}-{H}_{D-1}$ decreases for large D,
• a noise-dominated regime in which ${H}_{D}-{H}_{D-1}$ increases for large D,
• an intermediate regime in which ${H}_{D}-{H}_{D-1}$ remains nearly constant with D.

Zoom In Zoom Out Reset image size

In principle, this qualitative feature of the permutation entropy can be applied to experimental time series to assess whether the dynamics is dominated by noise or by the deterministic dynamics.

We remind that this distinction can not be done for the block entropy, as in this case ${H}_{D}-{H}_{D-1}$ is necessarily a decreasing function of D (see e.g. [35, 36]). This fundamental difference between permutation entropy and block entropy can be appreciated by comparing the left and right panels of figure 2.

For the analysis of time series of the Rössler system, we considered the Poincaré plane X = 0, shown in figure 3(a), and analyzed the sequence of time-intervals between consecutive crossings. Figure 3(b) shows the difference ${H}_{D}-{H}_{D-1}$ versus D, for different values of α. The solid line indicates the expected value if all ordinal patterns were equally probable, ${H}_{D}-{H}_{D-1}=\mathrm{ln}D!-\mathrm{ln}(D-1)!$. Because of the high level of stochasticity, we calculate the confidence interval that is consistent with the null hypothesis of equally probable ordinal patterns: in figure 3(b) the gray region represents the expected value $\pm 3\sigma$, where σ is the standard deviation calculated for a hundred surrogated (shuffle) time-series.

Zoom In Zoom Out Reset image size

Before testing the method in experimental data we want to investigate how the choice of the Poincaré section influences the results. We consider a Poincaré section in the plane $Z=\beta$, as shown in figure 4(a), and varying β in the range $[0.05-26.7]$, for a fixed value of $\alpha =0$. In this case, to discretize the time series, we analyze the time values when the trajectory intersects the Poincaré section and Z grows.

Zoom In Zoom Out Reset image size

Figure 4(b) displays the difference ${H}_{D}-{H}_{D-1}$ versus D, for different values of β. We can see that the difference ${H}_{D}-{H}_{D-1}$ increases with β. This is due to the fact that, as β is increased, consecutive values in the time-series become increasingly uncorrelated, similarly to when increasing the noise strength. On the contrary, for the minimum value of β, the variation of ${H}_{D}-{H}_{D-1}$ with D is resemblant to the behavior under almost-deterministic observed in figure 3(b).

Next, we analyze experimental data from the laser output intensity, displayed in figure 5(a). To discretize the data we consider the thresholds indicated with horizontal lines in figure 5(a), and analyzed the time intervals between consecutive threshold-crossings [7, 9, 32]. Figure 5(b) displays the difference ${H}_{D}-{H}_{D-1}$ versus D, for different thresholds. Note that ${H}_{D}-{H}_{D-1}$ varies with the threshold in a similar way as in figure 4(b): as the threshold decreases, correlations between consecutive dropouts are lost.

Zoom In Zoom Out Reset image size

For all the thresholds, ${H}_{D}-{H}_{D-1}$ grows monotonically with D. The reason is that the empirical time series is very noisy and the 'almost-deterministic' regime is not seen, not even for the highest threshold. Nevertheless, the values of ${H}_{D}-{H}_{D-1}$ lie outside the gray region that indicates values consistent with equally probable ordinal patterns. This reveals that the sequence of intensity dropouts are not completely uncorrelated, and thus, this method can determine regularities also in very noisy data.

The time-delayed Lang–Kobayashi model [27] is a widely studied model of laser dynamics, which generates a high dimensional LFF dynamics. It is therefore of interest to compare time-series generated by this model with the experimental ones. The model equations and parameters are as in [32] (in our case, to model the experimental situation, no current modulation is considered). We use the same threshold value as with the experimental time series. We calculated data sets of more than 10⁶ consecutive time intervals. Results are in figure 6. One can observe a very good agreement with the experimental results.

Zoom In Zoom Out Reset image size

Finally, we consider the issue of the length of the time series. If the time series is too short, the statistics to compute the probabilities of patterns is insufficient, and the entropy is underestimated. Figure 7 displays the estimated value of the permutation entropy versus the length of the time series, for different dimensions of the ordinal patterns. Notice that the data requirements increases with D. As the number of possible patterns of length D is equal to $D!$, we analyzed time series of length $N=300\ {D}_{\mathrm{max}}$ for the simulations, where ${D}_{\mathrm{max}}$ is the maximum dimension considered D = 8; and for the experiment $N=10\ {D}_{\mathrm{max}}$ with D = 7. The vertical dotted line marks the length corresponding to this criterion and demonstrates that the permutation entropy is computed with sufficient statistics.

Zoom In Zoom Out Reset image size