Entropy in Natural Time and the Associated Complexity Measures

Natural time is a new time domain introduced in 2001. The analysis of time series associated with a complex system in natural time may provide useful information and may reveal properties that are usually hidden when studying the system in conventional time. In this new time domain, an entropy has been defined, and complexity measures based on this entropy, as well as its value under time-reversal have been introduced and found applications in various complex systems. Here, we review these applications in the electric signals that precede rupture, e.g., earthquakes, in the analysis of electrocardiograms, as well as in global atmospheric phenomena, like the El Niño/La Niña Southern Oscillation.


Introduction
Non-equilibrium dynamics is the key for understanding self-organization [1].Self-organization is the main characteristic of complex systems that usually arise in physics, chemistry and mostly in biological sciences.The notion of self-organized criticality (SOC) introduced by Bak et al. [2] firstly provided a unifying picture [3][4][5] for the ubiquitous presence of power laws in various complex systems that definitely operate far from thermodynamic equilibrium.SOC reconciles non-equilibrium dynamics with the critical behavior obtained in equilibrium systems [6] without the necessity of a fine-tuned parameter, such as temperature or pressure in equilibrium.The mechanism behind this phenomenon is that the complex system self-organizes in a state for which power (or approximately power) laws appear naturally.Of course, not all complex systems or all self-organized systems give rise to power laws.A typical example of such a complex self-organized system is a living system [7].The physiological time series resulting from such a system do not exhibit huge, power law-like variations, but they rather vary within a well-defined range, which is necessary for life to be sustained.Yet, the (multi-)fractal behavior arises, and it is altered with disease or aging [7].Moreover, since physiological signals may contain both stochastic and deterministic components, the concept of entropy is also suitable for their study.
Within these lines of providing a general framework for the study of time series resulting from complex systems and focusing on the identification of forthcoming extreme events, our group [8][9][10] introduced natural time analysis (NTA) almost fifteen years ago.NTA enables the identification of novel dynamical features hidden behind the time series resulting from complex systems.The first applications of NTA have been reviewed by Varotsos et al. [11].In natural time, one can define a quantity reminiscent of the excessive "entropy" (see [12], pp.[26][27][28], called [8,13] entropy in natural time S. In the present work, we focus on the applications of natural time entropy to various complex systems in physics, geophysics, seismology, atmospheric sciences and cardiology.The paper is organized as follows: In Section 2, we present the necessary background of the analysis in natural time and introduce the entropy S, as well as the corresponding complexity measures in Sections 2.2 and 2.3, respectively.Section 3 focuses on the related applications to the electric and magnetic signals that precede rupture (see Section 3.1), the penetration of magnetic flux avalanches in type II superconductors (which is a typical SOC system; in Section 3.2), the Olami-Feder-Christensen (OFC) [14] earthquake model (Section 3.3), the preparatory processes before the occurrence of strong earthquakes (Section 3.4), the study of electrocardiograms (ECG), summarized in Section 3.5, and atmospheric physics in Section 3.6.These applications, as well as the future perspectives for new ones are discussed in Section 4. Section 5 accounts for various technical, though very important for real life, details on the aforementioned applications.Finally, the conclusions are presented in Section 6.

Analysis of Complex Time Series in Natural Time
In a time series of individual events resulting from a complex system, the natural time [8][9][10][11] associated with the k-th event is given by χ k = k/N, where N is the total number of events included in the time series.In NTA [8][9][10][11], the pair (χ k , Q k ) is studied where Q k is a quantity proportional to the energy emitted by the k-th event.This is made by considering the distribution: where p k is the normalized energy for the k-th event.Since the positive p k for k = 1, . . .N sum up to unity, they can be considered as probabilities (for more details, see [15]), and hence, the function [8][9][10][11]15]: may provide information about the distribution p k when ω → 0. Expanding Equation (2) around ω = 0, we obtain that: where κ 1 stands for the variance of natural time: with respect to the distribution p k .For a simple system, such as a classical gas in thermal equilibrium in a vessel where Q k might be considered as the energy of the particles that consecutively collide with the walls of the vessel, Q k are independent and identically-distributed random variables leading to p k → 1/N.This is the case of the so-called [11,13,16] "uniform" distribution leading to a value of κ 1 equal to κ u = 1/12 ≈ 0.083.For critical systems, where one can assume that p(χ) varies algebraically as a power law, Varotsos et al. [17] have shown that: for a variety of systems approaching criticality.Thus, κ 1 reaches the value of 0.070 for a critical system or 0.083 for a system exhibiting a stationary or quasi-periodic behavior.Figure 1 depicts how a time series from a variety of complex systems can be read in natural time.[18], (c) avalanches in 3D rice piles [19,20], (d) an ECG (for the so-called QT intervals), (e) a not obviously dichotomous electric (or magnetic) signal and (f) monthly Southern Oscillation Index (SOI) [21] values can be visualized in natural time.The meaning of the symbols, as well as the mathematical details for each case are discussed in Section 5.

Entropy in Natural Time
Apart from κ 1 , another useful quantity in NTA [11] is the entropy S given by [8,13,22]: where the brackets . . .≡ ∑ N k=1 . . .p k denote averages with respect to the distribution p k .The entropy S is a dynamic entropy that exhibits [23] positivity, concavity and Lesche [24,25] experimental stability.When Q k are independent and identically-distributed random variables, S reaches [22] the value S u ≡ ln 2 2 − 1 4 ≈ 0.0966 that corresponds to the aforementioned "uniform" distribution.
Moreover, upon reversing the time arrow and hence applying time reversal T , i.e., T p k = p N−k+1 , the value of S changes to a value S − ; see Figure 2. Hence, the entropy S in natural time does satisfy the condition to be "causal".As the concept of entropy is equally applicable to deterministic, as well as stochastic processes, the natural time entropies S and S − may provide a useful tool [11,22,[26][27][28][29] for the analysis of physiological time series (cf. the latter in most cases is due to processes involving both stochastic and deterministic components).It is worthwhile to mention that natural time entropy S has a particularly simple physical meaning since it may capture small trends existing in Q k .For example, when considering a small increasing trend (> 0) for p k versus k by studying the parametric family p(χ; ) ≡ 1 + (χ − 1/2) of continuous distributions for p k , one can show [27], using the definition of Equation (6), that: Expanding the last term of Equation ( 7) around = 0, we obtain that: Since for this family of continuous distributions, S − ( ) simply equals S(− ), we observe that an increasing trend in p(χ; ), i.e., > 0, corresponds to S − ( ) values higher than S( ).This result indicates that if we study the change of the entropy in natural time under time reversal: increasing or decreasing trends transform to negative or positive values of ∆S, respectively.Finally, since S( ) in Equation ( 7) is a nonlinear function of , we observe that the change of the entropy under time reversal ∆S is a nonlinear tool capturing alternations in the dynamics of the complex system.Using S and ∆S within a specified natural time window of length l, one can study time series resulting from complex systems by constructing the corresponding time series of S and ∆S obtained upon estimating these two quantities every l events as the natural time window of length l slides along the time series [30].This inspires the introduction of complexity measures that quantify the fluctuations of the entropy S and of the quantity ∆S upon changing the length scale, as well as the extent to which these quantities are affected when randomly shuffling the consecutive events.

Complexity Measures Based on the Entropy in Natural Time
Complexity measures quantifying the variability of S and ∆S upon changing either the (natural time window) length scale or shuffling the consecutive events randomly are classified into two categories depending on whether they make use of either S or the change ∆S of the entropy in natural time under time reversal [22,26,27,29,31].

Complexity Measures Based on S
As a measure of the natural time entropy fluctuations, we consider the standard deviation δS when we calculate the value of S for a number of consecutive events and study how S varies when sweeping this time window (each time by one pulse) through the whole time series.The following complexity measures based on δS have been suggested [22,26,31]: When the natural time window length changes from a short value, e.g., five events, to a shorter one, e.g., three events, the corresponding δS value also changes.This variation in the short(s) range is quantified by the measure λ s = δS 5 δS 3 , where the subscript in δS denotes the time window length chosen.If a longer (L) range, e.g., 60 events, changes to a short one, e.g., three pulses (events), the corresponding variation is quantified by another measure λ L = δS 60 δS 3 .Hence, one can in general quantify the effect of changing the natural time scale from three events to l events by considering the complexity measure: Considering now the data obtained by randomly shuffling the durations Q k of consecutive events, we can define the measure ν as: where δS and δS shu f denote the value of δS calculated when a natural time window of length l is sweeping through the original and the shuffled time series, respectively.The following ratios have been investigated [22]: ν s and ν L for the following natural time window lengths: 3-4 events from small scales and 50-70 events for the larger scales, respectively.It is clear that the quantity ν captures the extent to which the sequential order of events regulates the natural time entropy fluctuations.

Complexity Measures Based on ∆S
Complexity measures Λ and N can be defined [27,29], in a similar fashion with the measures λ and ν defined above, by using the time series of ∆S l obtained for a sliding natural time window l and considering its standard deviation σ[∆S l ] instead of δS.For example, we can define [27] the measure: or apply [29] the complexity measure: The measure N l quantifies the extend to which the ordering of the events contributes to the ∆S l values being equal to unity for a random process, whereas Λ l quantifies the variable effect of time reversal upon changing scale.

Applications of the Natural Time Entropy in Various Complex Systems
In this section, we will briefly present the various applications of S and ∆S in different complex systems.

Results for the Electric and Magnetic Signals that Precede Rupture
Seismic electric signals (SES) [12,[32][33][34][35] are low-frequency (≤1 Hz) variations of the electric field of the Earth that have been found in Greece [36][37][38][39][40] and Japan [41][42][43][44][45] to precede earthquakes.Many SES recorded within a relatively short time are termed SES activities and usually exhibit a dichotomous behavior; see Figure 1a.They are recorded at certain sites of the Earth's surface (sensitive sites) by using a multitude of measuring dipoles [33].SES exhibit the so-called selectivity phenomenon according to which SES related with earthquakes of a given earthquake mechanism from a given seismic area are recorded only at a given sensitive site(s) [33,35,[46][47][48][49][50].This enables the determination of the epicentral area based on the sensitive site that has recorded the SES.The experimental study of SES started [51] more than thirty five years ago in Greece, and the continuous collection of data enabled the construction of selectivity maps for several sensitive sites.When an SES is recorded at a sensitive site, the measured electric field variations are compared with previous SES recorded at the same site, and information on the epicentral area of the forthcoming earthquake can be obtained (for a recent review, see [11], Chapter 1).For SES registered at a given sensitive site and originating from a given seismic area (and earthquake mechanism), their amplitude expressed as ∆V/L (for a measuring dipole of length L with a given orientation) scales with the expected magnitude M of the impending earthquake according to the relation log 10 (∆V/L) = (0.34 − 0.37)M + β, where β is a constant depending on the dipole-seismic area and earthquake mechanism pair.The latter equation reflects that when the magnitude M increases by one-unit, the SES amplitude becomes almost double.All of the above properties can be understood [35,47,[52][53][54] on the basis of the SES transmission model suggested by Varotsos and Alexopoulos [12].The lead time of SES activities may vary [11,40] from three weeks to five and a half months, with an average around three months or so.For strong earthquakes, SES are accompanied by detectable magnetic field variations [55][56][57].Hence, it is of major importance to be able experimentally to distinguish the electric and magnetic field variations that precede the earthquake rupture from other similar-looking signals.
In this direction, Equation (5), which holds for critical systems, has been applied (e.g., see [9,13,16]) together with S and S − .In particular, both entropies S and S − have been experimentally found [13,58,59] to be smaller than S u for the electric signals that precede rupture: On the other hand, electrical noises that resemble in shape SES (they are of a dichotomous nature as shown in Figure 1a) violate in general one of the two conditions implied by the inequality (14).Table 4.6 of [11] presents an extensive list of recorded SES and noises for which Equation (5) and the above relations ( 14) can be used for a secure distinction between these two kinds of electric signals.Relations (14) have been also found to be valid [60] for long duration SES activities, such the one depicted in Figure 1e, although a more complicated analysis is necessary for this case (see Section 5).It is important to mention that S together with κ 1 have been also found to obey Equation ( 5) and Inequality (14) for the electric field variations that were recorded at Mexico [61] and preceded the magnitude M = 6.6 earthquake on 24 October 1993 and the M = 7.4 of the Guerrero earthquake on 14 September 1995.
As suggested by Varotsos [62], the experimental results on the electromagnetic (EM) precursory signals detected before earthquakes, as well as the ground-based ultra low frequency (ULF) magnetic field measurements (e.g., see [63][64][65][66] for such phenomena observed before the M = 9 Tohoku mega-earthquake on 11 March 2011) show that the main characteristics (e.g., amplitude, frequency range, lead time) of these precursory changes are similar to those discussed above for SES (which, as mentioned, for strong impending earthquakes are accompanied by detectable magnetic field variations).This observation has prompted the application [67] of the aforementioned criteria (Equation ( 5) and Relations ( 14)) to the pre-fracture EM emissions that preceded the 7 September 1999 Athens, Greece, M = 5.9 earthquake.The results show [67] that these criteria hold.The latter criteria have been also shown [68][69][70] to be valid for other EM emissions that preceded the strong earthquakes of 12 October 2013, 26 January 2014 and 3 February 2014 in Greece and the 2013 M = 6.3 earthquake at Kobe, Japan.

Results for the Penetration of Magnetic Flux Avalanches in Type II Superconductors
As first pointed out by de Gennes [71], when a type II superconductor is put in a slowly-ramped external magnetic field, magnetic vortices start to penetrate the sample from its edges.These vortices get pinned by crystallographic defects (e.g., dislocations), leading to the build-up of a flux gradient that is only marginally stable in a similar fashion as is the slope in a slowly-growing sand pile.Hence, it can happen that small changes in the applied field can result in large rearrangements of flux in the sample, known as flux avalanches [72][73][74].Thus, the system is very close to the original sand pile models suggested [2] as prototype for SOC.Natural time analysis of such magnetic flux avalanches have been performed [75,76] on the basis of the detailed measurements reported by Aegerter et al. [77] for thin films of YBa 2 Cu 3 O 7−x .It has been shown [76] that both S and S − satisfy the inequalities ( 14), although one cannot decide which of the two is larger than the other.
Moreover, these results, which correspond to a complex system at SOC, have been compared [76] with those obtained from 3D rice piles getting [19,20] progressively closer to SOC.The rice piles also resemble the prototype sand pile model of SOC, e.g., see [19,20,78,79], and the comparison has been made on the basis of the well-controlled experiments performed [19,20] on 3D rice piles.A typical example of such experiments is shown in Figure 1c.The results of Figure 2 of [76] indicate that both Equation (5), as well as Relation (14) hold in a fashion similar to that observed in the analysis of SES.This fact probably indicates that the system is entering a critical stage.Here, it is worthwhile to mention that similar results have been obtained [75,76] when analyzing the avalanches produced by a simple deterministic self-organized critical system introduced [80] to describe avalanches in stick-slip phenomena.Detailed simulations of this model have demonstrated [81] that SOC itself can spontaneously generate both critical avalanche statistics and long-range temporal correlations between avalanches in the presence of a slow temporarily-uniform external drive.It has been found [75,76] that Equation ( 5) holds together with S < S u as this system approaches SOC.When SOC is reached, the relations ( 14) hold only marginally [76].

Olami-Feder-Christensen Earthquake Model
The OFC earthquake model [14] is probably [82] the most studied non-conservative, supposedly, SOC model.It was conceived as a simplification of the Burridge-Knopoff spring-block model [83] for earthquakes by mapping it into a non-conservative cellular automaton.The OFC model introduced dissipation in the family of SOC systems.The spring-block model consists of a two-dimensional array of blocks arranged in a flat surface with friction.Each block is connected (by elastic springs) with its first neighbors and in the vertical direction to a driving plate that moves horizontally with velocity v.When the force acting on a block overcomes the static friction with the surface, the block slips (cf.[84,85] for important precursor dynamics before the transition to sliding).In the OFC model, the force on a block is stored in a site of a square lattice, and the static friction threshold is assumed to have the same value over all blocks.If energy input occurs in discrete steps instead of continuous and if thresholds are random, but not quenched, quasiperiodicity emerges combined with power-laws [82].The OFC model [14] runs as follows: we assign a continuous random variable z ij ∈ (0, 1) to each site of a square lattice, which represents the local "energy".Starting with a random initial configuration taken from a uniform distribution in the segment (0,1), the value z ij of all sites is simultaneously increased at a uniform loading rate until a site ij reaches the threshold value z thres = 1 (i.e., the loading ∆ f is such that z ij max + ∆ f = 1).This site then topples, which means that z ij is reset to zero, and an "energy" αz ij is passed to every nearest neighbor.If this causes a neighbor to exceed the threshold, the neighbor topples also, and the avalanche continues until all z kl < 1.Then, the uniform loading increase resumes.The number of topplings defines the size s of an avalanche or "earthquake".The quantity f = ∑ (∆ f ) represents the total increase of z ij due to the external force loading in each site.Since the loading rate is assumed uniform in time, f plays a role analogous to that of the conventional time T ≡ f .The coupling parameter α can take values from 0-0.25.Smaller α means more dissipation, and α = 0.25 corresponds to the conservative case.The parameter α is the only parameter of the model apart from the system size L, the edge length of the square lattice.The model can be supplemented by open boundary conditions in which the sites at the boundary distribute energy to the outer sites, which cannot topple; thus, energy is removed at the boundary.Another possibility is to use free boundary conditions: in this case, α varies locally: where n ij is the actual number of nearest neighbors of the site ij.For sites in the bulk n ij = 4; for sites at the edges n ij = 3; and for the four sites at the corners n ij = 2.The symbol K denotes the elastic constant of the upper leaf springs measured relative to that of the other springs between blocks [86].Obviously, the OFC model is non-conservative for K > 0 for which α ij < 0.25.
The criticality of the OFC model has been debated [87,88].Furthermore, the SOC behavior of the model is destroyed upon introducing some small changes in the rules of the model, e.g., replacing open boundary condition with periodic boundary conditions [89], introducing frozen noise in the local degree of dissipation [90] or in its threshold value [91], including lattice defects [92] (which should not be confused with the intrinsic lattice defects in solids [93], e.g., Schottky [94,95] or Frenkel [96,97] defects).Despite these findings, as well as others that show [98] that it is insufficient to account for certain aspects of the spatiotemporal clustering of seismicity, the OFC model appears to show many features found in real earthquakes.As far as earthquake predictability [99] or Omori's law [86,100] are concerned, the OFC models appear to be closer to reality than others [101].The predictability of the OFC model has been attributed to the occurrence of "foreshocks" (as well as "aftershocks") in the non-conservative case of the model [100].In addition, for certain values of the local degree "α" of dissipation (i.e., if "α" is chosen above 0.17), the OFC model exhibits an avalanche size distribution that agrees well [102] with the Gutenberg-Richter (GR) law [103], which states that the (cumulative) number of earthquakes with magnitude greater than or equal to M is given by: where b is a constant, which varies only slightly from region to region in the range 0.8 ≤ b ≤ 1.2 [104].
Considering that the seismic energy E released during an earthquake is related [105] to the magnitude through E ∼ 10 cm , where c is around 1.5, Equation ( 16) turns into: where γ = 1 + b/1.5.Hence, b ≈ 1 means that the exponent γ is around γ = 1.6-1.7.The above are some of the reasons why the OFC model is considered to be the prime example [106] for a supposedly SOC system for earthquakes, but the question of whether real earthquakes are described or not by SOC models of this type, or whether other kinds of mechanisms, e.g., [107][108][109], need to be involved, still remains unsolved [82,86,100,[110][111][112][113][114].
In the NTA of the OFC model [115], the size s was the quantity used as Q k in Equation ( 1), and the difference of the entropy in natural time under time reversal ∆S of Equation ( 9) has been evaluated for the time series of avalanches that preceded a large avalanche (e.g., one with s ≥ 1000 as in Figure 3).The values of ∆S have been stacked for a large number of such time series, and the results indicate (Figure 3; see also Figure 1 of [115] where a definition of ∆S opposite to that of Equation ( 9) has been used, or Figure 8.12 of [11]) a well-defined minimum of ∆S of Equation ( 9) exists that indicates that in the OFC model, there exists a time arrow as the system evolves towards a large avalanche.This can be interpreted as the existence of foreshocks or in general as an indication of the predictability of the OFC model.

Earthquakes
The most important power law for earthquakes is the GR law of Equation ( 16), which reflects that the seismic energy is emitted in the form of the power law of Equation (17).Both the entropy S and the entropy under time-reversal S − have been used for the determination of the occurrence time of an impending mainshock.This has been done [59,[116][117][118] in a variety of cases in Greece where SES have been identified.The procedure followed is that upon the recording of the SES, a candidate area for suffering a strong earthquake can be selected, and hence, the seismicity due to the small earthquakes occurring there can be studied in natural time as described in Section 5. Thus, upon the occurrence of each small earthquake, we can estimate the parameters κ 1 , S and S − .Experience has shown that a few days before the occurrence of the strong earthquake that gave rise to the recorded SES, Equation ( 5) as well as Condition ( 14) are valid.This has been independently verified before strong earthquakes in Greece [119,120] and Japan [44].As an example, we present in Figure 4 the analysis of the seismicity in the area N 38.6  36.0E 22.5 20.0 estimated by Varotsos et al. [121] based on an SES recorded on 14 January 2008 (see Figure 2(c) of [122]).This figure reveals that after the 56th event that occurred on 5 February 2008, both Equation ( 5), as well as Condition ( 14) have been continuously satisfied until the last event shown in this figure that occurred on 10 February 2008.This behavior has led to the public warning for an imminent strong earthquake in the aforementioned region that was announced on the 10 February 2008 issue of the Greek newspaper Ethnos (see the related discussion in Section 21.4 of Lazaridou-Varotsos [40], as well as [38,39]).Four days later, the strongest earthquake (M w (USGS) = 6.9) in Southern Greece for the last 30 years occurred at 10:09 UT on 14   [121] that it will suffer a strong earthquake based on a seismic electric signal (SES) recorded on 14 January 2008 (see Figure 2(c) of [122]), after discarding the two events (earthquakes) that were related to another SES activity [122].The black sticks correspond to the magnitude M L (ATH) (right scale) reported by the Geodynamical Institute of the National Observatory of Athens for each small earthquake.The horizontal red and blue lines correspond to 0.07 and S u , respectively.

Electrocardiograms
Physiological time series are best suited for NTA based on complexity measures based on the entropy S.This is so because usually, such time series give rise to κ 1 very close to κ u , as well as S close to S u due to the quasi-periodic behavior they exhibit.Hence, for the study of time series coming from ECG recordings, one should construct the corresponding time series of S or ∆S obtained for certain scales l of the length of the natural time sliding window.The ECG analyzed come from publicly-available databases [123] (for more details as well as the defitinitions of the RR, NN, QRS, and QT intervals that will be mentioned below see Section 5), and the following key results have been obtained [11,22,[26][27][28][29]: 1.When analyzing the RR and QRS intervals of the ECG in natural time, the subjects suffering from sudden cardiac death (SCD) violate [22] one or more of the four healthy limits related to ν s (RR), ν L (RR), ν s (QRS) and ν L (QRS). 2. When analyzing the QT intervals of the ECG in natural time δS l (QT) for l = 3-8 heartbeats, the SCD subjects exhibit [26] almost one order of magnitude larger fluctuations than those of the healthy ones.
3. The NTA of RR, QRS and QT intervals allows [26] the distinction of individuals into three categories: healthy, SCD and heart disease patients.Details on this ternary distinction can be found in [11,26].4. When employing ∆S in the analysis of the RR and NN intervals of long duration ECG of SCD patients, the fluctuations of ∆S 7 appear to maximize [27] during the last three hours before the ventricular fibrillation.5.When employing ∆S in the analysis of the RR and NN intervals of long-duration ECG, N 3 (RR), N 3 (NN) together with ∆S 7 (RR) and ∆S 7 (NN) may allow the separation [27] of individuals into four classes: healthy, SCD, congestive heart failure (CHF) and atrial fibrillation (AF) individuals.6.The analysis suggested in the previous point is also valid [28] for models of healthy or patient heart dynamics.

When analyzing the NN intervals of long-duration ECG recordings, the combination of σ[∆S 7 ]
with Λ 7 and Λ 49 enables (see Table 1 of [29]) the ternary distinction in healthy, SCD and CHF individuals.

Atmospheric Physics
Use of natural time entropy has been also made in atmospheric physics.

Ozone Hole Dynamics over Antarctica
Varotsos and Tzanis [124] applied S and S − in the time series of the maximum daily ozone hole area over Antarctica for each year during the period 1979-2009.The results obtained show that S for scales 3-7 years and S − for all scales (3-15 years) almost stabilizes during the last several years.On the other hand, characteristic features are found before the unprecedented event of the major sudden stratospheric warming and the subsequent break-up of the Antarctic ozone hole into two holes in September 2002.In particular, the following precursory changes have been identified [124]: First, for scales larger than eight years, S exhibits a gradual increase after around 1999.Second, from 2000 to 2001, S − shows an increase for all scales (3-15 years) except for the scale of 13 years.Third, the values of ∆S almost coincide at 2000 for the short scales 3-7 years and then decrease.Finally, the analysis in the natural time domain was also applied [124] on the eddy heat flux, which is proportional to the vertically propagating wave activity affecting the ozone hole over Antarctica, and the results drawn confirm those deduced from the ozone hole area diagnostics.
3.6.2.Forecasting the Intensity of El Niño/La Niña Southern Oscillation Varotsos et al. [125] have identified that when analyzing the monthly SOI values in natural time (see Section 5), the time series of the change of the entropy in natural time under time reversal estimated for periods of 20 months, i.e., ∆S 20 , may act as a predictor for the intensity of the El Niño/La Niña Southern Oscillation (ENSO) events (see also [126]).They suggested that when ∆S 20 exceeds a threshold ∆S thres , an El Niño event is expected to occur.In other words, the time series of ∆S 20 captures the dynamics of ENSO and exhibits clear precursory signals before strong El Niño events.Within this frame, the successful prediction that the 2015-2016 El Niño event would be rather a "moderate to strong" or even a "strong" event and not "one of the strongest on record", as that of 1997-1998, has been made [127].
More specifically, in Figure 5, we depict the ∆S 20 and Southern Oscillation Index (SOI) values' time series and indicate by a vertical arrow the data available at the time at which it was decided that the ongoing 2015-2016 El Niño event will not evolve into one of the strongest events ever recorded, i.e., October 2015 (the paper was initially submitted [127,128] on 13 November 2015).As is also shown in this figure, the subsequent values of SOI never left the green colored region that corresponds to the minimum values of SOI for "weak, weak to moderate, moderate, moderate to strong" El Niño events validating the aforementioned prediction.The latter has been based on the data presented in Figure 6 where the probability density function (PDF) of ∆S 20 is shown (for details on its determination, see [127]) and compared with the extreme ∆S 20 values measured during very strong El Niño events.An inspection of this figure together with Figure 5  The alarm is set on (black line) when ∆S 20 exceeds the threshold value ∆S thres = 0.0035 (red horizontal line).The selection of the threshold can be seen in Figure 2 of [127] and has been made as a compromise between the costs of false positive and false negative predictions.The two colored horizontal stripes represent the mean minimum negative values of SOI along with the one standard deviation bands for the two cases of "weak, weak to moderate, moderate, moderate to strong" (green band) and "strong, very strong" (yellow band) El Niño events [127].

Discussion and Perspectives
As we can see from the results presented in the previous section, natural time entropy S, as well as its value S − obtained under time-reversal is useful for identifying when a system approaches a critical stage, as is for example the case of electric signals before rupture (see Section 3.1), SOC systems (such as those studied in Section 3.2) or the seismicity before a strong earthquake (see Section 3.4).In all of the cases mentioned above, the validity of the inequalities ( 14) characterizes (together with other conditions, e.g., see [118,129]) the approach to criticality.
Apart from this application, the study of the time series of the difference ∆S of the entropy in natural time under time reversal may reveal useful information concerning the evolution of the system towards an extreme event.This has been discussed here in Section 3.3 for the OFC earthquake model and in the two cases in atmospheric physics of Section 3.6.Such applications may be very fruitful in the future as the kind of the complex system studied may vary.Hence, the use of ∆S time series for different complex systems is strongly encouraged.
Finally, for complex systems exhibiting quasi-periodic behavior, such as the case of heart rate variability (see Section 3.5), one should better study complexity measures based on the time series of S and ∆S such as those defined in Section 2.3.These can reveal the interplay between the various scales that are involved in the regulation of the behavior of the complex system (cf.regulation is a prerequisite of quasi-periodic behavior) and may be useful for the identification of severe situations (such as SCD), although the system may look to operate normally.Application of the results obtained in Section 3.5 to a wider range of individuals, as well as the possibilities of their practical use in the medical practice are currently under study.

Materials and Methods
For the analysis of dichotomous electric signals that precede rupture (e.g., SES activities [9,13,16]), we first convert the experimentally-recorded signal into its dichotomous representation (0 = no signal, 1 = signal), and then, we can easily obtain the power (see the middle panel of Figure 1a) and, hence, the energy emitted during each pulse, which is proportional to the duration of the pulse.As shown in the lower panel of Figure 1a, in NTA, we consider as Q k the duration of the k-th pulse.
In the case of earthquakes, instead of the earthquake magnitude M usually reported for each earthquake in conventional time (see the upper panel of Figure 1b), in NTA, we consider as Q k the seismic energy emitted during an earthquake, which is given by [105] E = 10 1.5M w +11.8 erg, where M w is the moment magnitude.For the case of the earthquakes shown in the upper panel of Figure 1b, appropriate conversion relations between M and M w suggested by Tanaka et al. [130] have been used to obtain the values of Q k depicted in the lower panel.In the latter panel, the quantity ε corresponds to the seismic energy emitted by an M = 3.5 earthquake.
For the case of 3D rice piles, e.g., see the upper panel of Figure 1c, which comes from the detailed experiments of Lőrincz [20] (see Figure 5.1 there), in natural time we consider as Q k the size of the k-th avalanche, as shown in the lower panel of Figure 1c.
In the case of the analysis of ECG, where traditionally the turning points are labeled by the letters P, Q, R, S and T (see the upper panel of Figure 1d), various durations (intervals) corresponding to the time elapsed between these points, such as the QT interval shown in the lower panel of Figure 1d, have been considered as Q k in NTA.For example, the RR interval which is the duration between two successive R turning points has been also considered as Q k in NTA.Other intervals that have been also used for the same purpose are the NN intervals, which are the durations between consecutive normal beats and intervals between pairs of normal beats surrounding an ectopic beat are discarded (see also [29] and references therein), and the QRS intervals that correspond to the duration of the QRS complex in a cardiac cycle of the ECG (see the upper panel of Figure 1d).Moreover, the ECG data analyzed in natural time so far (e.g., see [131]) come from various databases available at physiobank [123].As mentioned in Section 3.5, for the analysis of such ECG intervals in natural time, the time series of S and S − should be calculated when using a natural time window of length l.This window slides each time by one interval (event) in the ECG interval time series, so that when starting from the m 0 -th event, we have: S(m 0 , l) = χ ln χ w − χ w ln χ w (18) where χ ln χ w = ∑ l k=1 p k,w χ k,w ln χ k,w , χ w = ∑ l k=1 p k,w χ k,w with: and χ k,w = k/l.Similarly, S − (m 0 , l) is calculated by Equation ( 18) when p k,w of Equation ( 19) is substituted by: The ∆S l time series is obtained by the differences ∆S l (m 0 ) ≡ S(m 0 , l) − S − (m 0 , l), m 0 = 1, 2, . . ., N − l.
For a not (obviously) dichotomous signal, such as the electric signal depicted in the upper panel of Figure 1e that preceded [60] the M w 6.4 earthquake at 38.0 • N 21.5 • E on 8 June 2008, the analysis in natural time is made in two steps: First, we estimate the "instantaneous" power in the sense that we use data at the shortest time interval available from our (assumed even) sampling and construct the time series P i , e.g., see the middle panel of Figure 1e.Then, by selecting a threshold P thres (see the horizontal line in the same panel), we add every consecutive P i that is larger then P thres to estimate the energy emitted during an event.When P i falls below P thres , the event ends, and a new event starts when for the first time P i becomes larger than P thres again.This way, we obtain a representation of the signal in natural time as that shown in the lower panel of Figure 1e.By varying the threshold within a reasonable range, which can be determined by the distribution of P i , e.g., see Figure 6 of [60], we can obtain secure conclusions concerning the values of κ 1 , S and S − of the signal in natural time.
Finally, for data such as those usually collected in atmospheric physics involving hourly, daily or even monthly observations, which may acquire both positive and negative values, such as Troup [21] SOI, depicted in the upper panel of Figure 1f, the NTA can be performed [125,127] by adding an appropriate constant so that all values become non-negative, as shown for example in the lower panel of Figure 1f, and then proceed to the study of the resulting time series by ∆S l .

Conclusions
Here, we reviewed the applications of the natural time entropy S, its value under time reversal S − and the related complexity measures in a variety of complex systems.As an outlook for the future, one could propose the study of S and S − in systems that may approach criticality, the study of the ∆S time series as a general tool for the identification of precursory to extreme event phenomena and the use of complexity measures for the characterization of complex systems exhibiting quasi-periodic behavior.

Figure 1 .
Figure 1.How a time series of (a) dichotomous (e.g., zero or one) electric signals, (b) earthquakes[18], (c) avalanches in 3D rice piles[19,20], (d) an ECG (for the so-called QT intervals), (e) a not obviously dichotomous electric (or magnetic) signal and (f) monthly Southern Oscillation Index (SOI)[21] values can be visualized in natural time.The meaning of the symbols, as well as the mathematical details for each case are discussed in Section 5.

Figure 3 .
Figure 3. Results from averaging the last 1000 events before a large avalanche (s ≥ 1000, occurring at T 0 ) in the Olami-Feder-Christensen (OFC) model with L =100 and K = 2: the change ∆S (left scale) of the entropy in natural time under time reversal and the mean energy ζ ≡ ∑ L i=1 ∑ L j=1 z ij /L 2 (right scale) as a function of the "time" (T 0 − T) to the large avalanche.Note that ∆S minimizes before the occurrence time T 0 of the large avalanche and changes sign when ζ starts to increase.

Figure 4 .
Figure 4.The values (left scale) of κ 1 (red circles), S (blue lines) and S − (cyan lines) obtained from the study of seismicity within the region N 38.6 36.0E 22.5 20.0 , estimated on 1 February 2008 by Varotsos et al.[121] that it will suffer a strong earthquake based on a seismic electric signal (SES) recorded on 14 January 2008 (see Figure2(c) of[122]), after discarding the two events (earthquakes) that were related to another SES activity[122].The black sticks correspond to the magnitude M L (ATH) (right scale) reported by the Geodynamical Institute of the National Observatory of Athens for each small earthquake.The horizontal red and blue lines correspond to 0.07 and S u , respectively.

Figure 5 .
Figure 5.The entropy change ∆S 20 in natural time for the window length l = 20 months (red line, left scale) along with SOI monthly values (blue line, right scale) for the period January 2014-March 2016.The alarm is set on (black line) when ∆S 20 exceeds the threshold value ∆S thres = 0.0035 (red horizontal line).The selection of the threshold can be seen in Figure2of[127] and has been made as a compromise between the costs of false positive and false negative predictions.The two colored horizontal stripes represent the mean minimum negative values of SOI along with the one standard deviation bands for the two cases of "weak, weak to moderate, moderate, moderate to strong" (green band) and "strong, very strong" (yellow band) El Niño events[127].

Figure 6 .
Figure 6.PDF of ∆S 20 (blue curve, left scale) together with the corresponding histogram (red bars, left scale) obtained from the time series of ∆S 20 , which is also plotted versus time (green crosses, right scale) along the vertical axis.The arrows indicate when ∆S 20 exceeds 0.0205 and are labeled by the corresponding ongoing strong El Niño events.Taken from [127].
In the upper panel, the stronger pulses are emitted periodically, while in the lower panel, consecutively in the middle of the process.Although the Shannon entropies for both panels are equal, the entropies in natural time S are different in each panel.Moreover, when using the entropy in natural time under time reversal S − , we obtain values that are in general different from those of S.