{'fraction_non_alphanumeric': 0.04075372480280456, 'fraction_numerical': 0.009430324276950043, 'mean_word_length': 4.597076424997548, 'pattern_counts': {'":': 0, '<': 7, '<?xml version=': 0, '>': 6, 'https://': 0, 'lorem ipsum': 0, 'www.': 2, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
The Fokker-Planck Approach to Complex Spatio-Temporal Disordered Systems

When the complete understanding of a complex system is not available, as, e.g., for systems considered in the real-world, we need a top-down approach to complexity. In this approach one may start with the desire to understand general multi-point statistics. Here such a general approach is presented and discussed based on examples from turbulence and sea waves. Our main idea is based on the cascade picture of turbulence, entangling fluctuations from large to small scales. Inspired by this cascade picture, we express the general multi-point statistics by the statistics of scale-dependent fluctuations of variables and relate it to a scale-dependent process, which finally is a stochastic cascade process. We show how to extract from empirical data a Fokker-Planck equation for this cascade process, which allows to generate surrogate data to forecast extreme events as well as to develop a non-equilibrium thermodynamics for the complex systems. For each cascade events an entropy production can be determined. These entropies fulfil accurately a rigorous law, namely the integral fluctuations theorem.


INTRODUCTION
For quite some time research on complex systems has been considered as a continuation of investigating nonlinear or chaotic dynamics. The main difference between these systems may be understood by realizing that nonlinear or chaotic systems are spatially homogeneous and, thus, are described by lowdimensional nonlinear differential equations, cf. (1,2), whereas complex systems possess spatial and temporal inhomogeneities. Due to the interdependence, relationships, or interactions between units of a complex system the understanding of the entire system is not attainable by simply understanding each part, or by the local features. Complex systems are in general composed of many interacting subunits, where nonlinearities play an important role, so that complex spatio-temporal structures emerge; see e.g. (3,4). A consequence of the interaction between the subsystems and the overall behaviour is that it is often difficult to achieve full comprehensive understanding of complex systems dynamics. Additionally, surprising new collective phenomena may emerge. Examples of emergent behaviours include short-and long-term climate changes, hurricanes, cascading failures, evolution, learning, and intelligence, to name just a few (5).
In this article we will take turbulent flows and sea waves as examples of complex systems. The main task for a good understanding of the appealing complexity of flow patterns, like shown in Fig. 1, is to characterise the clearly visible structures, as well as their large variability (6). Watching such turbulent flows, one recognises immediately the flow type by its overall structure, but at the same time one gets the impression that over the time exactly same patterns are never seen twice. This mutuality of order and stochasticity is one exciting aspect of flow patterns. The two examples, selected for this article, reflect also nicely this mutuality. One of the challenging problems of turbulence is the small scale structure and its deviation from Gaussian statistics, cf. (6,7,8). The anomalous statistics can be seen in connection with the millennium problem, defined by Clay Mathematic Institute, where it is asked for the local structure of a solution of turbulent flows described by the Navier-Stokes equation (9). An open question is whether there are special small scale coherent structures explaining the anomalous statistics. For sea waves coherent structures seem natural, but here we should note that we are not interested in the case of periodic wave structures but in the cases of the rough sea. Hereof a most prominent wave structure is the monster wave, also called freak or rogue wave. Still an open question remains whether these structures are part of the disordered wave state or somehow independent of it (cf. (10,11,12)). Turbulent-like features of waves are treated as wave turbulence. This leads to the open problem of what the basic features of such complex structures are. The questions can be formulated, whether there are some clear structures (coherent deterministic structures) that can be singled out and may serve as a skeleton to access the complexity. Alternatively, one may ask, whether such systems can be understood best by their stochasticity and statistics? Quite often one approaches such complex systems in a pragmatic way by either studying the structures or the stochasticity.
Another common approach to characterise such complex patterns is to apply linear correlation mea- Smoke visualisation of flow patterns of a turbulent wake behind a sphere with dimples. The overall structure of the wake flow is recognised right away by its large structures, highlighted by arrows. A closer look at this turbulent wake will show that flow patterns are never repeated exactly either in time or space, see blow-up. For many complex disordered systems the challenge is to understand this interplay between clear structures and stochasticity.
sures, as such patterns decay in time or, respectively, in space. Typically for complex systems, like for turbulence, no simple exponential decay of such correlations is found, consequently we are facing the problem of multi-scale correlated systems. Indeed, for such systems the low order statistics, i.e. one point statistics and two-point correlations, are not sufficient to grasp the observed complexity.
Our present work has to be seen in a clear difference to the works devoted on the analysis of nonlinear dynamics with and without noise in the framework of time series analysis. For time series anaylsis the proper embedding, Lyapunov-exponents, fractal dimensions, fixed points, stable and unstable limit cycles, reconstruction of dynamic equations etc. are of interest, cf. (13,5). Although methods of time series analysis will be used, this work focuses on the extended spatial disordered structures with the aim to get a comprehensive characterisation at arbitrarily many points, i.e. by a general N-point joint statistics. The development of such a method has been stimulated by research on turbulent flows (14) and may be considered as a top-down approach. Knowing such a general description, it should be possible to determine all statistical aspects of the system. Moreover, the common problem of structures versus statistics should be sorted out, as the general joint N-point statistics can grasp any sequential ordering of some patterns as multi-point structures. Also the mentioned multi-scale correlations and higher order statistics can be captured by general N-point statistics. The question is how complicated such joint statistics become. For the empirical estimation of N-point statistics the question whether sufficient data can be provided rises immediately.
For a general approach to N-point statistics, we propose a hierarchical ordering of the N-point statistics. This hierarchical ordering is in analogy to the common cascade picture of turbulence, which describes how structures on larger scales interact with structures on smaller scales in a hierarchical way, so that a downwards cascade from large to small scales is obtained. Inspired by the idea of a cascade we investigate how the scale-dependent structures will change with the scale at each location. Vividly interpreted, this can be taken as a zooming-in process of the complex structure. To get access to the high demanding multi-point statistics, we set the scale-dependency of the complex structure in the context of a stochastic process evolving in scale. The novelty is that we do not consider the common time-evolution of stochastic processes but an evolution in scale. In particular we show evidence that this zooming-in process can be approximated by a Markov process, i.e. that this process has no memory in its scale evolution. This approximation allows to derive a Fokker-Planck equation, evolving in scale, for the hierarchical ordering of the N-point statistics. The Fokker-Planck equation is not only a compact description of the whole complexity, but also enables to derive several other aspects ranging from scaling behaviour to thermodynamics, as outlined in this paper.
This work has its origin in a series of former works started in 1996 (15,14,16). Initially the idea of a scale-dependent process has been worked out without paying much attention to how to nest larger structures into smaller structures, and can be seen as a continuous formulation of the propagator description of the cascade (17,18). Reviews on this approach can be found in (5,19), where the stochastic processes in general and the difference between commonly known stochastic processes in time and the new processes in scale are worked out. In (5) applications and citations of stochastic processes in scale are given, which range from turbulent flows, financial data and surface roughness to earthquakes, cosmic background radiation and iEEG recordings from epilepsy patients. All these examples posse remarkable multi-scale features and the complexity seems to be related to a hierarchical ordering connected to cascade-like structures. At this stage, the expression of multi-point and multi-scale statistics was used more or less synonymously. With the attempt to reconstruct time series from the knowledge of the multi-scale processes (20) the meaning of the correct placement of the smaller scale structures within the larger ones became clear and the relation between mutli-scale and multi-point statistics has been worked out. In (21) this was done for financial market data and a short-time forecasting has been worked out. For turbulence data (22) it has been realised that an extended class of stochastic cascade processes, expressed by a family of Fokker-Planck equations, is needed. In the present work we work out in detail the multi-point approach and relate this to the stochastic cascade processes in scale. The technical details for the handling of empirical data will be given in the corresponding sections.
In our approach, we first consider the simplification of one-dimensional cuts of complex patterns (for turbulence this simplification is related to the Taylor's hypothesis of frozen turbulence). Thus, a quantity q(x) along an axis x is considered. Second, the hierarchical ordering is introduced by changing the scales r, thus we ask how the structure looks like on different scales r, where the changes go from large to small distances, as explained in the next section. Furthermore we show how the N-point statistics can be expressed by a joint multi-scale statistics. In Sec. 3 a three-point approach or three-point closure for the hierarchical multi-scale statistics is described, which finally opens up the possibility of projecting the general N-point statistics on the stochastic processes in the scale parameter r, ending in a scaledependent Fokker-Planck equation, see Sec. 4. In Sec. 5 special self-similar or fractal solutions of the stochastic cascade process are discussed. Two further consequences are deduced from this approach. On the one side, we show in Sec. 6 that surrogate data sets can be produced with the same statistical properties and patterns (original processes). On the other side, in Sec. 7, the stochastic approach is put in the context of non-equilibrium thermodynamics for the complex systems, relating complex structures with the general fluctuation theorem.

TAYLOR'S HYPOTHESIS OF FROZEN TURBULENCE
In 1938 G. I. Taylor introduced a hypothesis by which he deduced the spatial fluctuations of a turbulent velocity profile from the corresponding measurements of temporal fluctuations at a single point. This hypothesis, known as the Taylor frozen-flow hypothesis, relies on the existence of a mean flow u that translates the spatial structures past a stationary probe in a time smaller than the inherent evolution time of the fluctuations (23). q(x + dx) ≈ q(x − u dt), where we take the mean flow in x-direction. Therefore, analyses are taken to be equivalent regardless of whether they are taken as snapshot in space (see Fig. 1) or as time sequence of the structures passing over a sensor by the mean flow velocity. Here, we discuss only the spatial complexity.

MULTI-POINT STATISTICS EXPRESSED BY INCREMENT STATISTICS
One important basic aspect of this work is the connection between the general multi-point characterisation of a complex structure and its multi-scales properties. We start with a quite formal consideration and show how this connection can be worked out mathematically, the next sections show consequences and apllications. We consider the case that a complex structure is given as space and time dependence of a quantity q( x, t). For the example of a turbulent flow q( x, t) is given by a velocity field u( x, t). For a surface q( x, t) is the spatial pattern of the height h( x, t). As mentioned in the introduction, we simplify this system by assuming that its temporal and spatial structures are statistically similar, and in addition only one direction is of interest. We also assume that the characterising quantity is a scalar q(x). We are interested in multi-point statistics, i.e. the probability of finding a sequence of events q(xi) for several discrete locations xi with i = 0, · · · , N , which is given by the joint probability density function (jPDF) where we used the abbreviation qi := q(xi). Here W (q0, q1, · · · , qN )dq0 · · · dqN is the probability that the random variables q(x0), q(x1), · · · , q(xN ) belong to the intervals q0 ≤ q(x0) ≤ q0 + dq0, · · · , qN ≤ q(xN ) ≤ qN + dqN . Instead of this joint N + 1-point PDF W , one may be interested in the conditional probability of obtaining the value q at one selected point under the condition of the remaining events.
In the following we select the value of q at the last point xN as reference value. Therefore the conditional probability of finding q(xN ) for the given preceding data q(xi), with i = 0, . . . , N − 1, is given by the conditional probability density function (cPDF), which can also be taken as a transition probability, The multi-point probabilities can also be expressed in another way by considering the statistics of relative changes from one selected point. For xN as the point of reference we denote the distance ri := xN − xi. Therefore we introduce increments (other notations are common in the literature like δrq(xi), ∆q(xi), qr(xi), ...) for i = 0, . . . , N − 1, which quantify the differences of q over the distances or scales ri, as illustrated in (4) |J| denotes the determinate of Jacobin for transformation (q0, . . . , qN ) → ξ0, ξ2, . . . , ξN−1, qN and is unity. Here W (qN ) is the one-point probability density function (PDF) for the value qN . Based on the natural ordering x0 < x1 < . . . < xN the scales ri are ordered as ri > ri+1. Note that we have defined a scale evolution of ri running with the index i from large to small scale, which can be illustrated as a process where one zooms in to resolve smaller and smaller structures. With Equation (4) we have expressed the general N + 1-point statistics by the statistics of N increments ξi taken in a right-justified way from the point qN . If the statistics of the complex structure is homogeneous (or, for time dependencies, stationary) the probability W (ξ0, ξ2, . . . , ξN−1, qN ) does not depend on the location xN but on the values of the describing quantity, here qN . describing quantity: q(x i ) = q i relative distance: there are other common notations like qr δq.
x q(x)

Figure 2
Scheme of the hierarchical ordering of increments ξ i = q(x N ) − q(x N − r i ) and scales r i , in order to describe N -point statistics. Sorting the scales after their sizes, from large to small, one obtains a zooming-in process, to which a trajectory ξ i (·) can be assigned. This process then decribes the evolution of the increments ξ i when scales r i evolve from largest to smallest size. for the scale-dependent quantityξi(ri), which could be besides an increment also another scale-dependent quantity (see Sidebars on Wavelets). When increments ξi = ξ(xN , ri) are considered for multi-point statistics, the definition of the reference point xN is of special relevance. Besides the left-jusitified increment ξ l i (xN , ri) = q(xN + ri) − q(xN ) also a right-justified definition ξ r i (xN , ri) = q(xN ) − q(xN − ri) can be used, as we do here. For the multi-scale statistics also a centered version ξ c i (xN , ri) = q(xN + ri/2) − q(xN − ri/2) has been used, cf. (24,25), where further details on the relation between these different definitions are given. In this paper we restrict ourselves to the right-justified definition, as this is required for the multi-point reconstruction of data in section 6.

MULTI-POINT AND MULTI-SCALE STATISTICS
Note that from the multi-point statistics the multi-scale statistics can be derived, but as the reference value qN is not any more taken explicitly into account in the multi-scale statistics, one can not derive the multi-point statistics from the multi-scale statistics in general, see also (21,22).

CLOSURES OF MULTI-POINT STATISTICS
The introduction of the hierarchical ordering of the increments ξ0, ξ1, . . . , ξN−1, with ri > ri+1, (see Fig. 2) has been used so far just to reformulate the N-point statistics by increment statistics. As a next step we consider this hierarchical ordering of the increments as scale-dependent fluctuations of the quantity q that go from large to small scales or vice versa. Later on, this will become an essential aspect for working out a cascade idea for the description of complex systems. Before we come to this point some formal aspects of the joint probabilities have to be discussed. The jPDF of Equation 4 can be expressed by a product of conditional probabilities A tremendous simplification arises, if the multi-conditioned PDF only depends on the increment of the next larger scale p(ξi|ξi−1, ..., ξ0, qN ) = p(ξi|ξi−1, qN ).
We then obtain a much simpler form of Eq. (5), i.e.
Note that these simplified cPDF p(ξi|ξi−1, qN ) are three-point statistics, is a three-point closure of the the general (N + 1)-point jPDF. As a remark we mention, that possible closures are of central interest for the turbulence problem, and that in (26) such a three point closure is discussed for the Lundgren-Monin-Novikov Hierarchy, a description of turbulence by multi-point probabilities. Two further simplifications are given, first, if the cPDF are independent of the reference qN and second, if the cPDF are independent of larger increments The last conditions correspond to a two-point closure for which the general N + 1-point probability factorises completely to products of simple one increment or one scale probabilities as rely on the dependence on qN . Only for the case that the increment statistics are also independent of the reference value p(ξi|qN ) = W (ξi), a complete knowledge of the single increment PDFs W (ξi) for all scales characterise completely the multi-scale disorder of the considered complex structure, an aspect which is important for the fractal characterisation of complex structures, see also Sec. 5. We can conclude that a characterisation of a complex system, which is done only by the statistics of increments based on W (ξi), is a special two-point characterisation. With the knowledge of the PDF W (ξi) all higher order moments of ξ n i for all scales are known, but nothing is known on more than two-point correlations. Thus the question of how far the general multi-point problem can be reduced is of central importance for the proper characterisation of complex systems.
Having data from concrete complex systems, the validity of the simplifications can be tested. As we focus in this contribution on the reduction to three-point statistics, see Equation (6), the validity of the simplification can be seen by investigating p(ξi|ξi−1, . . . , ξ0). This is easily done by determining different increments for the same reference value qN . As an example the results from one turbulent data set is shown in Fig. 3. Here the quantity of the system is the local velocity in the direction of the mean flow, thus qi = ui. It can clearly be seen that p(ξ3|ξ2, ξ1) depends on ξ2, as the contour lines are not parallel to the ξ2-axis. Thus the simplification of Equation (9), namely the reduction to two-point statistics, does not hold. It can, however, be seen that the double conditioned PDF p(ξ3|ξ2, ξ1) shown as red contours, is similar to the single conditioned PDF p(ξ3|ξ2). This result is a good indication that Equation (6), the three-point closure, holds. If many data are available, conditions on further larger-scale increments can be investigated. The quality of how well this condition is fulfilled can be tested by statistical tests, see for example (14,27) and has been found for many data sets of turbulence (28,29,30) and other data like financial data (e.g. (21)) and surface heights, see (5). In (22) it is shown that this simplification holds also for p(ξ3|ξ2, ξ1, qN ) = p(ξ3|ξ2, qN ) and that for turbulence data p(ξ3|ξ2, ξ1, qN ) depends on the reference values qN , too. (Another way of showing Equation (6) for experimental data is given by (31)).
Single-conditional PDF p(ξ(r 3 ) | ξ(r 2 )) are presented as solid lines, double-conditional PDF p(ξ(r 3 ) | ξ(r 2 ), ξ(r 1 ) = 0) as symbols. σ∞ is the standard deviation of the velocity values u i , or respectively q i . The conditional PDF shows two points clearly. First, the dependency on the condition is seen in part (a) by the change with x-axis, this shows that p(ξ(r) | ξ(r )) = p(ξ(r)) or a two-point closure is not supported. Considering features of the increments for only one selected scale, like done by investigating structure functions, see Equation 27, are incomplete. Second, the conditional PDF do not depend on a second condition on an increment on an even larger scale r , p(ξ(r) | ξ(r ), ξ(r )) = p(ξ(r) | ξ(r )). This is taken as an indication that the r-evolution of the increments is memoryless, i.e., the knowledge of the value of the increment ξ(r) suffices to determine the next step by which the system evolves to ξ(r ), with r < r < r.

FOKKER-PLANCK EQUATION IN SCALE
So far the characterisation of complex disordered structures has been discussed by multi-point and multi-scale statistics, as well as the possible simplifications of three-and two-point closures. Next, the hierarchical ordering of the increments shown in Fig. 2 together with the simplification of a three-point closure, which is achieved by Equation (6), is put in the context of cascade processes, for which we work out a description by stochastic differential equations. Thus the aim is to grasp the whole complexity by some stochastic equations. The basic first idea is to look at the increment for a chosen location xN , i.e. ξ(xN , ri) as a quantity that changes with r and denote it by the increment trajectory ξ(·) ≡ ξ(xN , r), which describes the above mentioned zoom-in process. Equation (6) is now nothing else than that ξ(xN , ri) depends only on the increment of the next larger scale ξ(xN , ri−1). At the same time ξ(xN , ri) is independent of further increments on larger scales. This means that the evolution of ξ(·) has no memory, note, this is the definition of a Markov process in scale r. Thus this evolution has Markov properties and can be considered as a stochastic process evolving in r, or more precisely, based on our definition, evolving with decreasing r. (For readers who are not familiar with stochastic processes, we refer to the sidebar 4 and

WAVELETS, INCREMENTS, AND CORRELATIONS
The description of the complex structure by multi-increment statistics can be set in analogy to a wavelet analysis, cf. (33,17,34,35). Wavelets ψ a,b (x) are characterised by a scale a (width) and a location b. Using the difference of two Dirac functions δ(x), the scale a = ri and the loaction b = xN it is possible to define the following wavelet also sometimes called poor man's wavelet. Increments are nothing else than the coefficients of these wavelets As discussed in (17), we analyse the evolution of these coefficients of the wavelets with scales as stochastic processes. In principle the discussion in this paper can also be performed with general wavelets ψr,x(·) and their coefficientsξ(x, r). A main difference is that the increment statistics can be related directly to Npoint statistics. For instance the correlation functions are given by the second oder moment ξ 2 (x, r) = 2 q 2 − 2 q(x − r)q(x) . Consequently, higher-order and mixed-order correlations are directly related to higher-order moments of the increment PDF, such as ξ n (x, r) . the continuative literature (36,37,38).) For experimental or empirical data like from turbulence or from financial markets (40) the absence of memory gets lost for the smallest scales. Such a behaviour has already been proposed as a natural one by Einstein in his pioneering work in Brownian motion (41). This defines a lower bound scale, which we call Einstein-Markov length rEM and which can be determined by the validity of Equation (6) as xN−2 converges against xN−1 (30). Thus our consideration is valid for r > rEM and may be treated as a small scale cut-off, see e.g. (42). Note that rEM is more than a lower bound, but is also the finite step size that coarse-gains the whole r-evolution from largest to smallest scale. In this way the Markov process may be taken as a stochastic process modelling in a continuous manner the coarse-grained process.
The evolution of the cPDFs p(ξ|ξi, qN ) with r < ri describes the transition probability of ξi(ri) → ξ(r) for the given reference qN . An equation for the evolution of this transition probability with extending the difference between r and ri is given by the Kramers-Moyal expansion (more precisely the Kramers-Moyal forward expansion) (36) Note that if p(ξ|ξi, qN ) is independent on qN and if one multiplies the equation with W (ξ) and integrates over ξi, Equation (17) reduces to www.annualreviews.org • Fokker-Planck Approach
Therefore marginal PDF W (ξ1) and conditional PDF p(ξ k |ξ k−1 ) are sufficient to describe Markov processes. The probability distributions (marginal and conditional) of Markov processes satisfy a partial differential equation of order one in the scale and order infinity in the state variable ξ. The governing equation is known as Kramers-Moyal (KM) equation, see Equation 16.
In this respect the Pawula theorem states that there are only three possible cases in the KM expansion: i)the Kramers-Moyal expansion stops at n = 1 means that the processes are deterministic, ii) the KM expansion stops at n = 2, the resulting equation is the Fokker-Planck or Kolmogorov equation and describes diffusion processes and finally iii) the Kramers-Moyal expansion stops at n = ∞. Any truncation of expansion at finite order n > 2 would produce non-positive probability density W (ξ) (36).
For the case (ii) the KM expansion reduces to the Fokker-Planck equation, which means that the first and second KM coefficients D (1) (ξ, r) (drift coefficient) and D (2)  where noise η(r) is a zero mean white Gaussian with intensity 2, which means that η(r)η(r ) = 2δ(r − r ).
which leads to the pure description of the increment statistics of W (ξ). This is clearly less information than p(ξ|ξi), which we know from the discussion of the multi-point statistics in the previous section. Note that D (n) (ξ, r, qN ) and D (n) (ξ, r) are related by D (n) (ξ, r, qN )p(ξ|ξi, qN )W (qN )dqN = D (n) (ξ, r)p(ξ|ξi).
Looking at the evolution of the increments ξr with scale, the Markov property of Equation (6) means that only delta correlated noise acts on the trajectory. The reduction of the Kramers-Moyal expansion of Equation (16) goes along with the requirement that the involved noise in the stochastic process is not only delta correlated but has also Gaussian distribution. This is also called Langevin noise, for which a corresponding differential equation for a single event or path ξ(·) is given as

COMMENTS ON MARKOV PROCESSES direction of the process
The evolution of the scales from large to small values is considered, for this ourpose we used in Equation 16 a negative prefactor of −r. Having −r/∂r = −1/∂ln(r) shows that we have implicitly used a log scaling of the r-evolution, which is of advantage for complex structures with self-similar properties, see Chapter 5.

inverse direction of the process
If it is shown that the data fulfill the Markov conditions Equation 6 from large to small scales, the Markov condition is fulfilled also in the other direction, from small to large scales (39) p(ξi|ξi+1, ..., ξN−1, qN ) = p(ξi|ξi+1, qN ) where η(r) denotes zero-mean Gaussian white noise with a variance of 2, i.e. η(r)η(r ) = 2δ(r − r ).
Here we use the Itô interpretation. For Stratonovich and other descriptions see cf. (36,37). From the Langevin equation it is evident that D (1) describes the deterministic part of this equation and is called drift coefficient. The function D (2) (ξ, r), which is called diffusion coefficient, determines the amplitude of the noise. The case that D (2) changes with ξ is called multiplicative noise. An essential point for the stochastic description of the scale dependent increments is the knowledge of the Kramers-Moyal coefficients D (n) , which can be determined directly from the data as conditional moments (cf. (36,5)).
First let us define the n-th order moments for two increments in scales that are separated by δ The values of M (n) depend on the value of δ for some chosen or fixed values of ξ, r, qN . The Kramers-Moyal (KM) coefficients D (n) are given by The definition presented here for KM coefficients differ by the factor of r with common definition, which is due to our description of a stochastic process in scale, see Equation (19). It should be mentioned that this definition of the Kramers-Moyal coefficient can already be found in an early work of Kolmogorov in 1931 (43).  We see that the cPDFs, p(ξ |ξ, qN ) play again an important role, as for the lim δ→0 the differential equation (16) can be estimated from their knowledge. Note that this lim δ→0 can be considered as a fusion process of two increments ξ j → ξj, which is of interest for a field theoretical approach to such complex systems (5,44,26).
An important aspect of Equations 21,22 is that they can also be read as a concept to estimate the Kramers-Moyal coefficients directly from given data, as shown for some turbulence data in Fig. 4. That this is a very efficient method to analyse time series of noisy dynamical systems has been shown in (5). Experiences have shown that technically the limit lim δ→0 is best performed by investigating the moments M (n) (ξ, r, δ, qN ) with the help of the small step approximation in δ, for which M (n) (ξ, r, δ, qN ) = D (n) (ξ, r, qN )δ + O(δ 2 ), see Fig. 4b. If the given data do not allow this small step approximation, due to insufficient sampling rates or due to a too large Einstein-Markov length, respective corrections can be calculated (45,46). In Fig. 4b the sampling rate was sufficiently high. The deviation from a linear law for small values δ is due to the Einstein-Markov length. Other corrections arise, if additionally measurement or observation noise or another non-ideal noise contribution is given cf. (47,48,49,50).
In Fig. 5 the Kramers-Moyal coefficients are shown for a turbulent data set. The linear behaviour of the drift term D (1) and the quadratic behaviour of the diffusion term D (2) becomes clear. D (2) has an additional additive offset. The fourth order Kramers-Moyal coefficient can be taken as zero within the experimental precision. This indicates that for turbulence a Fokker-Planck equation can be used to describe the cascade process and thus the whole multi-point statistics. For the case that D (4) does not become zero, in principle infinitely many Kramers-Moyal coefficients have to be determined. A criterium to quantify the importance of higher order Kramers-Moyal coefficients has been worked out in (27), and in (29) it has been shown for turbulence that in contrast to the velocity fields the passive scalars require such higher order Kramers-Moyal coefficients. Another aspect of the results for the Kramers-Moyal coefficients shown in Fig. 5 is that the dependency on the reference point qN is only clearly present for D (1) . For positive qN the fixed point D (1) (ξ) = 0 is shifted to negative values, whereas for negative qN the fixed point is shifted to positive values. This result simply means that for positive qN the increments have the tendency to become more negative, what is in accordance with the boundedness of turbulent velocity data from stationary experiments.

Figure 5
First, second and fourth order Kramers-Moyal coefficients D (n) (ξ, r, q N ) calculated from data of a turbulent flow, using the procedure illustrated in Fig. 4. Three different reference values q N , labeled as q = −σq ± σq/6, q = 0 ± σq/6, q = σq ± σq/6, were chosen. To get sufficient data for the statistics, always a small but finite interval around the q N values was chosen. Here σq is the standard deviation of W (ξ). Note that only D (1) (ξ, r, q N ) shows clear dependence on q N . For further details see (22). The negative values for the fourth order Kramers-Moyal coefficient are due to the extrapolation procedure as indicated in Figure 4. Once these Kramers-Moyal coefficients are successfully estimated as in this figure, the corresponding Fokker-Planck equation is completely obtained.
calculated by the short time propagator (36,27) pST A sufficiently small step δ has to be used to stay in this limiting approximation. Based on this short time propagator any cPDF p(ξr|ξ r ) can be determined. Therefore the quality of the estimation of the coefficients D (1) and D (2) can be verified by comparing the cPDFs obtained from the data with those obtained by solving the Fokker-Planck equation with the estimated coefficients (51).

HIGHER-ORDER KRAMERS-MOYAL COEFFICIENTS
According to the Pawula theorem vanishing higher order Kramers-Moyal (KM) coefficients, specially fourth-order D (4) (ξ, r), guarantee that the process is statistically contineous and the Kramers-Moyal expansion Eq.(16) can be truncated after the second (diffusive) term cf. (36). For vanishing KM coefficients one can construct the Langevin equation with the computed drift function and diffusion coefficients form time series. Non-vanishing higher-order (n > 2) KM coefficients have been observed in various systems cf. (5), which indicates that the corresponding measured time series do not belong to the class of continuous diffusion processes and jump events should play a significant role in the underlying stochastic process (52).

SELF-SIMILARITY AND FRACTALS
With the previous Section 4 the derivation of a Fokker-Planck equation as a model of scale-dependent complexity has been completed. We now want to put this approach in the context of other analysis frameworks for complex systems, namely self-similarity and fractals. For complex structures the question is often posed, whether they possess self-similar structures, also called fractals. Particularly for the two examples of turbulence and sea waves, discussed here, the concept of self-similarity plays an important role, cf. (6,10). We start the discussion of self-similarity in a quite general way with the principles of scaling symmetries, from which we derive properties of the so-called structure functions ξ n . The structure functions are very often used for the characterisation of both turbulence and surface roughness like sea waves. We show how the concept of the stochastic cascade equations of the previous Section can describe these structure functions and their self-similiar structure.
Commonly the self-similarity is investigated by a local measure, which characterises the structure on the scale r at the location x. We denote the local measure again as ξ(x, r). Self-similarity means that in a certain range of r the quantities should have the same statistics. More precisely, the probability distribution of the quantity ξ takes the form with a universal function F (Q). The universality of F leads to the scaling behaviour Such a type of behaviour has been termed fractal scaling behaviour. The concept of fractals is widespread and many examples are known, like turbulence or surface roughness, just to mention two. The strict self-similarity expressed by Equation (24) is often just an idealised approximation. In fact the so-called multifractal behaviour is often more appropriate. Here the k-th order moments scale according to where the scaling indices ζ(k) are now not any more linear but a nonlinear function of the order k. Such a multifractal behaviour can formally be obtained by the assumption that the probability distribution W (ξ, r) has the following form This formula is based on the idea that the complex system is composed of subsets of different scaling indices α, whereW (α, r) gives a measure of the scaling indices α at a scale r (see e.g. (18,53) ). A shortcoming of the fractal and multifractal approach to complexity in scale is the fact that it only addresses the statistics of the measure ξ(x, r) at a single scale r. As we have derived above, one has to expect for a general N-point characterization dependencies of the measures ξ(x, r) and ξ(x, r ) from different scales, as well as dependencies in the value of a reference point q(x). The connection between the fractal and multifractal characterisation and the stochastic cascade description can be derived from the Kramers-Moyal expansion of Equation (16). The validity of the Markov property or, respectively, the three-point closure is assumed and the Equation (16) has been integrated over all values of qn, so that the dependency on the reference value is not any more taken into account The multiplying of this equation with ξ k and the partial integration over ξ, e.g. (27) leads to If the Kramers-Moyal coefficients have the form D (n) = dnξ n (where dq are constants) (54), scaling behaviour of (27) is guaranteed with Based on this formula it can be worked out which combinations of dn will result in the q-dependent function of ζq, which characterizes different multifractal models. For turbulence and increments (ξr) k is a common quantity to characterise different flow situations. (ξr) k is called the kth order structure function. If the Kramers-Moyal expansion truncates to a Fokker-Planck equation the structure functions can be found from the following equation From Fig. 4 it is obvious that the drift coefficient has a linear behaviour, D (1) (ξr) = d11 ξr, and the diffusion coefficient has a quadratic behaviour, D (2) (ξr) = d20 + d22 ξ 2 r . The dij may be r-dependent. The scaling index becomes now Due to the additive term in D (2) (ξr) a mixing of different structure functions with different orders takes place. For the case of d20 = 0, which is not supported by experimental data (see Fig. 3), and for constant values of d11 and d22, the so called Kolmogorov 1962 (K62) or lognormal model (55) is obtained with the intermittency parameter µ Thus the K62 scaling corresponds to d11 = − 3+µ 9 and d2 = µ 18 . The corresponding relation between stochastic processes and the other well known multifractal scaling models of turbulence has been worked out by (56, 26).

SURROGATE DATA AND FORECASTING
The Fokker-Planck equation (FPE) derived in Section 4, despite its compactness, achieves a comprehensive and powerful characterisation and description of a wide range of complex systems. However, the generation of surrogate data which fully obey a special FPE is not trivial. This section will develop an approach to this task.
Based on the relation of the general N -point statistics of a complex structure and a stochastic description by a Fokker-Planck equation we obtain the possibility to generate new data sets numerically or to forecast special events. Therefore we consider the case that preceding values q0, · · · , qN−1 fix the probability of a new value qN . Taking this row of values qi, i = 0, · · · , N − 1 as a sequence of events, the cPDF of Equation (2) can be seen as the prediction of the next event for the next time step (see e.g. (57)). Such a predictor can be expressed by the stochastic cascade process, using Equations (4,7) p(qN |qN−1, . . . , q0) = W (ξ0, ξ1, . . . , ξN−2, ξN−1, qN ) W (ξ0, ξ1, . . . , ξN−2, qN−1) Note that the increments ξi, defined in Equation (3), have to be taken from the reference value qN in the nominators and from reference value qN−1 in the denominators, respectively. Equation (37) enables to determine the probability of the new value qN based on the knowledge of the simple conditional PDFs p(ξi|ξj, qN ), which can either be calculated from the Fokker-Planck equation or which can be estimated directly from the data. As p(ξi|ξj, qN ) contains only knowledge of three values qi, qj and qN of the data, this is again a three-point closure of multi-point statistics.
The conditional probabilities p (qN |qN−1, . . . , q0) contain all relevant statistical information of the previous data points for a correct choice of the new value qN . Choosing now a random value from this distribution, the time series will be extended correctly by another point. Shifting the procedure by one step and repeating the same procedure may be used to generate new surrogate time series, which exhibits the correct joint probability density function for all considered scales. For technical reasons one should avoid zeros in cPDFs of the denominator of Eq. (37). The initial idea for reconstructing time-series following this procedure was developed in a similar way for fluid turbulence data, see (20), and has been used for turbulent data (22), for financial data (21), and also for sea waves (32). In Figure 6 we show two time series of wind speed measurements. In the upper panel the originally measured time series is shown, in the lower panel a time series obtained by the just mentioned reconstruction method. The color-coded left part represents the initial conditions of the first N values q0, ..., qN−1, used to start the reconstruction method. As this is a stochastic model, involving a deterministic as well as a random part, the two time series diverge quite fast. But the stochastic content in the sense of multi-point statistics is the same, which can be verified by reanalyising these surrogate data (22,21,32). Another interesting point is that apparently typical structures of a wave pattern could be reproduced by the stochastic method (58), thus it seems that the multi-point approach is capable to grasp the statistics as well as coherent structures. This will only work, if such structures are based on the special stochasticity and it will not work if special structures are added to a noisy background.
The method to reconstruct data sets with the conditional probabilities p (qN |qN−1, . . . , q0) can also be used for a short time forecast, as it was shown for financial data (21) and sea waves (32). In Fig.  7(a) typical time series of wave heights is shown. Note that the big wave at the end of the time series corresponds to a measured rogue wave. In the figure part (d) and (e) two selected conditional probabilities p(qN |q0, r0, . . . , qN−1, rN−1) are shown to illustrate our method. In addition to the conditional probabilities the single event probability p(qN ) = p(q) of all height values is shown (red curve). These figures show clearly how the conditional probabilities change with q0, r0, . . . , qN−1, rN−1 the values of the N wave heights seen before. There are cases when smaller qN values are expected in the next step, see Fig. 7(b), and there are cases when large qN values become highly likely, see Fig. 7(c). With this method a warning system for approaching large wave heights can be set up. The high quality of such a prediction was quantified according to the receiver operating characteristic curve (ROC) (32).

NON-EQUILIBRIUM THERMODYNAMICS OF COMPLEX HIERARCHICAL STRUCTURES
Up to this point we have outlined a statistical approach to characterise completely the disordered structure in one direction by multi-point statistics. This approach was a phenomenological one. For the case Time series of wind data as a real-world example of numerical data generation, cf. Sec. 6. The upper panel shows the measured data, while the lower panel of the figure presents reconstructed data using Equation (37). The data, which are used as initial condition for the reconstruction, are marked on the left sides of the panels It is well visible that the reconstructed data are not identical with the measured data it can be shown that they but follows the same statistics. Not only mean value and standard deviation are correctly reproduced, but also multi-point statistics and higher-order correlations.
of turbulence there are some works which show how such a Fokker-Planck equation can be related to the basic equations of fluid mechanics. In (59) the connection to the Navier-Stokes equation is shown, in (26) the Lundgren hierarchy was analysed based on a three-point closure or, respectively, Markov properties in scale. For the sea waves comparable results are not known to us. In this last Section we put our statistical approach in the context of non-equilibirum thermodynamics (for earlier approaches to non-equilibrium approaches see for example (60)). Based on the derived Fokker-Planck equations for the cascade process we can assign entropy values to each local structure of the complex systems. For these entropy values the validity of a fluctuation theorem, namely the integral fluctuation theorem, can be shown. This is a way for how the phenomenological stochastic approach can be linked to fundamental laws of physics, cf. (61).
In particular, the concept of stochastic thermodynamics is applied to turbulent flows (62, 63, 64) and sea waves (58). The novelty here is that concepts of non-equilibrium thermodynamics known to hold for microscopic systems are shown to be valid also for such macroscopic systems. These concepts enable to determine an entropy production of the cascade process. In particular for every individual trajectory ξ(·) = ξr; r = r0, ..., rN of the increments evolving from large to small scales, a total entropy production ∆Stot can be defined by The total entropy production is given by the sum of two contributions, ∆S med being the entropy variation due to the surrounding medium, which depends on the evolution of ξ(·) through the hierarchy of length scales r in the cascade. Here ∆Ssys is the entropy change of the system itself. In Equation ( Dealing with this thermodynamics (see also (62, 61)) one may interpret ∂ ξ ϕ(ξ(r)) in a less formal way as a force of the medium given by the "mean field" quantities D (1) and D (2) . The interaction of this force on the path "velocity" ∂rξr leads then to the entropy term ∆S med which represents an analogue of the work done by the medium on the single event ξ[·], which leads to a heat exchange with the bath (for more details of this analogy see (65)). The second entropy term ∆Ssys may be considered as an intrinsic contribution of the trajectory. The main point is that to each increment trajectory ξ[·] a value of ∆Stot can be determined like shown in Fig. 8. Thus "microscopic" entropy fluctuations can be determined, which may show positive and negative values. Most interestingly the negative entropy events seem to be related to extreme events in the increment statistics on the smallest scales, as can be seen in Fig. 8 and as reported in (62, 58). This again points in the direction that salient structures of a complex system can be a proper part of the multi-point statistics, somehow unifying the approach to complex systems by coherent structures or by statistical methods. If the complex structure is described correctly by the Fokker-Planck equation, the statistics of the entropy values should fulfill the integral fluctuation theorem (IFT) a fundamental entropy law of non-equilibrium thermodynamics, cf. (61). Here · · · N denotes the average

THERMODYNAMICAL INTERPRETATION OF TURBULENCE
Entropy values for cascade trajectories ξ(·) allow a thermodynamical interpretation. The potential ϕ(ξr) of equation 41 given by the drift and diffusion term D (1) (ξ) and D (2) (ξ) can be considered as the coupling of the trajectory or the subsystem to the bath, whereas the values of ξr and its probabilities are the intrinsic features of the individual trajectory. The connection with a possible thermodynamical interpretation becomes more clear if not only the velocity increments are considered but also the transferred energy r of the cascade. As shown in (66) where γ, G and m are positive vales which may depend on r. Note that now the increment process becomes purely additive, a well known effect (67,68). Such an in-stationary (r-dependent) Langevin equation can be interpreted in a thermodynamic way, following (65). Most interestingly r corresponds to the temperature. As r devolves its own fluctuations, the cascades can be considered as a mixture of temperatures. Note that after equation 43 these energy or temperature fluctuations increase in the cascade evolution to smaller scales. In this way the cascade pictures of Kolmogorov (55) and Castaing (18) are set in a new thermodynamic context.  over many different trajectories for the increments. In Fig. 9 the distribution of the entropy production www.annualreviews.org • Fokker-Planck Approach (a) Probability density function of the total entropy production ∆Stot obtained from turbulent data, after (69). The mean of the distribution is positive, which means that in average the entropy increases, in accordance with the second law of thermodynamics. For our system here we find quite big fluctuations of the entropy values and pronounced probabilities for negative entropy events exist. Such a distribution is a typical result for data like shown in Fig. 8. (b) Convergence of the exponential of the entropy production e −∆S tot for turbulence data, following the Integral Fluctuation Theorem (IFT) in Eq. 44. Shown is the evolution of e −∆S tot N as a function of the number N of trajectories ξ(·) and its convergence to the value 1. This result shows that the IFT is fulfilled within an accuracy of 1% and better, see (69). Most interestingly the IFT implies that the probabilities of the negative and positive entropy events are not arbitrary, but balanced by the IFT, which puts much weight on the negative events. In other words, events with negative entropy values must be accompanied by many events with positive entropy to fulfil the IFT law.
values for an experimental data set of a turbulent flow (27) is shown. Clearly the mentioned positive and negative entropy values can be seen. The mean value of this distribution ∆Stot is positive. By a weighting function e −∆S tot negative entropy values contribute much more and must be compensated by many large positive ∆Stot values so that the IFT Equation (44) is fulfilled. Thus the IFT is a relation which expresses the balance between the relative frequency of entropy-consuming (∆Stot < 0) and entropy-producing (∆Stot > 0) trajectories associated with the stochastic evolution of increment trajectories ξ(·) (individual stochastic trajectories).
In Fig. 9b the obtained values from Equation (44) for increasing numbers of trajectories are shown. The convergence to the absolute value of 1 is well obtained. Already after several thousand values the IFT gets fulfilled.
Another interesting result is obtained for sea waves (58). It has been found that extreme events, namely the rogue waves, are characterised by negative entropy values. Comparing different states of the sea waves it could be shown that the statistics W (∆Stot) change significantly from one state to the other, although for both cases the IFT was fulfilled in high quality. Coming back to the point that the IFT somehow balances the negative and positive entropy events. As the negative entropy events are correlated to large waves, like the rogue waves, we see that this non-equilibrium thermodynamics together with the stochastic cascade process grasps both, the statistics and the localised structure of the complex disordered system.

CONCLUSIONS
The leading topic of this work was the characterisation of complex systems and the question if an understanding of the complexity can be achieved by structures as basic elements or if higher order statistics are needed. For two examples, namely turbulence and sea waves, we showed how these two aspects of structure and statistics are interwoven. The description of the multi-point statistics by a stochastic process of a cascade, or, respectively, by a Fokker-Planck equation evolving in scale, allows to generate surrogate data sets as well as to determine entropy values for all data points. Most interestingly for our considered macro-systems large fluctuations of these entropy values are found, as it is up to now mainly discussed for mirco-systems cf. (61). We see how this concept of micro-systems can fruitfully be applied to our considered macro-systems. A key element is the Markov property in scale and the corresponding derivation of the Fokker-Planck equation. It is the Fokker-Planck equation which leads to a general law of non-equilibrium thermodynamics, namely, the integral fluctuation theorem, which is fulfilled for our data with an accuracy of 1 % and better. On the one hand the validity of this integral fluctuation theorem can be taken as an evidence of the consistency of our whole approach. On the other hand the integral fluctuation theorem expresses mathematically the balance between negative and positive entropy values. Each negative entropy values must be compensated by many positive entropy values to fulfil the integral fluctuation theorem. The exponential weight of the theorem means that each negative entropy values must be compensated by many positive entropy values. This seems to be also true for the structures of the complex systems, as we showed that the negative entropy values are connected with the large small-scale structures, which are the challenging properties for turbulence and waves. Thus we conclude that our work presents a new consistent approach to the mutuality of order and stochasticity in complex systems.

SUMMARY POINTS
1. N-point statistics representing an all-encompassing probabilistic approach to complex systems can be expressed by N-1 scale increment statistics. The increment statistics allows a hierarchical ordering. If the increment statistics only depend on increments of the neighbouring scale, a three point closure of the N-point statistics is achieved. 2. The three point closure of the N-point statistics is equivalent to a stochastic process in scale with Markov property, for which the process equations can be estimated via Kramers-Moyal coefficients from empirical or measured data. If Langevin noise is present a non-stationary Fokker-Planck equation for the cascade process in scale is obtained. 3. The knowledge of the non-stationary scale-dependent Fokker-Planck equation allows to generate numerically new data sets with the same N-point statistics and to forecast single events, most interestingly also extreme events. 4. Based on the non-stationary scale-dependent Fokker-Planck equation the entropy production of the cascade trajectory can be defined, for which the rigorous integral fluctuation theorem holds. Thus a connection with non-equilibrium thermodynamics is given, which balances the occurrence of negative and positive entropy events.

FUTURE ISSUES
1. The non-stationary scale-dependent stochastic description of N-point statistics can be generalised in a straight forward way to higher dimensional quantities q, like complex turbulent velocity fields, see (70). Here the Fokker-Planck equation depends on different variables of the vector field and D (1) becomes a vector, D (2) is a diffusion matrix. The problem how to extend this approach to two-or three dimensional spaces instead of the one dimensional cut, for which a hierarchical ordering is evident remains open. 2. It is a challenge to work out a meaningful non-equilibrium thermodynamics of these complex structures, relating it to quantities like energies of the systems. It should also be noted that there have already been different tries to set up thermodynamical approaches to complex systems. A relation between those would be important. 3. Often complex systems are also described by nonlinear partial differential equations. Is there a possibility to derive the non-stationary scale-dependent stochastic process equation directly from the partial differential equations? This would unify at least for these systems two different ways of description.