Analytical representation of Gaussian processes in the A−T plane
AAnalytical representation of Gaussian processes in the
A − T plane ∗ Mariusz Tarnopolski † Astronomical Observatory, Jagiellonian University, Orla 171, PL-30-244 Krak´ow, Poland (Dated: January 1, 2020)Closed-form expressions, parametrized by the Hurst exponent H and the length n of a timeseries, are derived for paths of fractional Brownian motion (fBm) and fractional Gaussian noise(fGn) in the A − T plane, composed of the fraction of turning points T and the Abbe value A .The exact formula for A fBm is expressed via Riemann ζ and Hurwitz ζ functions. A very accurateapproximation, yielding a simple exponential form, is obtained. Finite-size effects, introduced bythe deviation of fGn’s variance from unity, and asymptotic cases are discussed. Expressions for T forfBm, fGn, and differentiated fGn are also presented. The same methodology, valid for any Gaussianprocess, is applied to autoregressive moving average processes, for which regions of availability ofthe A − T plane are derived and given in analytic form. Locations in the
A − T plane of somereal-world examples as well as generated data are discussed for illustration.
I. INTRODUCTION
The characterization and classification of time series[1–4] is an important task in a variety of fields. Severalmethods have been developed for tasks such as detectingchaos [5], measuring complexity via entropies [6–8], esti-mating the Hurst exponent [9–11], distinguishing chaoticfrom stochastic processes based on graph theory [12], andmore. A connection between chaos and long range depen-dence was also established for low-dimensional chaoticmaps [13].The
A − T plane [14], spanned by the fraction of turn-ing points T and the Abbe value A , was initially in-troduced to provide a fast and simple estimate of theHurst exponent H , as tight relations of both A and T with H were discovered for fractional Brownian motion(fBm), fractional Gaussian noise (fGn), and differenti-ated fGn (DfGn). While A ( H ) and T ( H ) strongly over-lap for H ∈ (0 ,
1) for different processes, in the jointspace ( A , T ) the fBm and fGn intersect only at the pointcorresponding to white noise, i.e. (1 , / H based on the empirical relation A ( H ) and computed us-ing a wavelet method [15, 16] were consistent with eachother.The discriminative power of the A − T plane wasdemonstrated in a multiscale scheme [17] by employingcoarse-grained sequences, i.e. dividing the time seriesinto nonoverlapping segments of length τ and calculatingthe mean in each segment. This produces smoothed se-quences, and the evolution of A and T with varying tem-poral scale τ allowed to separate (i) developed, emergingand frontier stock markets; (ii) healthy and epileptic pa-tients based on their EEG recordings—moreover, for the ∗ Dedicated to my Mother for her 65th birthday. † [email protected] first time it was possible to distinguish healthy patientswith closed and open eyes; and (iii) patients with andwithout cardiac diseases based on the heart rate vari-abilty. Finally, it is also possible to differentiate betweenchaotic and stochastic processes based on the differentbehavior of paths, parametrized by τ , in the A − T plane[18]. The
A − T plane is therefore a powerful tool withseveral possible applications.The aim of this paper is to derive analytical descrip-tions of fBm and fGn, as well as autoregressive mov-ing average (ARMA) processes, in the
A − T plane. InSect. II the known facts about turning points are reca-pitulated. In Sect. III the exact and approximated ex-pressions for A are derived, for the first time, for fBmand fGn. Their A − T plane’s representation is depictedin Sect. IV. ARMA processes are discussed in Sect. V.Various applications, ranging from pure mathematics, tobiology and astrophysics, to financial markets, are brieflyoutlined in Sect. VI. Summary and concluding remarkswith future prospects are gathered in Sect. VII.
II. FRACTION OF TURNING POINTS, T A. Theory
Consider three values x t , x t + d , x t +2 d of a time series { x t } . For d = 1 the points are consecutive. Assume thereare no ties between the neighboring points, which forcontinuous processes, or empirical data with decent res-olution, should not be an issue (see also [19, 20]). Threepoints can be arranged in six ways, identified by an orderpattern π p (Fig. 1). If the smallest value among the threeis given an index 1, and the largest an index 3, then e.g.the relation for one of the four possible turning points, x t < x t +2 d < x t + d , is described by a pattern π p = 132.Denote the probability of encountering a pattern π p by p π p . Then the following theorem holds [21]: Theorem 1
For a Gaussian process X t with stationaryincrements, p = p = α/ , and the other patternsyield probability (1 − α ) / . a r X i v : . [ phy s i c s . d a t a - a n ] D ec
123 321 132 231 213 312
FIG. 1. Order patterns for three points.
The probability p ( d ) for a given delay d is given by p ( d ) = 1 π arcsin (cid:114) ρ ( d )2 , (1)where the correlation coefficient is ρ ( d ) = E [( X d − X ) ( X d − X d )] E (cid:104) ( X d − X ) (cid:105) . (2)The probability T of encountering a turning point amongthree consecutive points (i.e., for d = 1, the case consid-ered hereinafter), as per Theorem 1, is then T = 1 − p (1) . (3)Note that T ∈ [0 , T = 2 / B. T for fBm, fGn, and DfGn With the following properties of fBm B ( t ), fGn G ( t ),and DfGn Y ( t ), one can utilize the methodology fromSect. II A to calculate T for them as a function of H : E (cid:2) B ( t ) (cid:3) = t H , (4a) E [ B ( t ) B ( s )] = 12 (cid:0) t H + s H − | t − s | H (cid:1) , (4b) E (cid:2) G ( t ) (cid:3) =1 , (4c) E [ G ( t ) G ( s )] = 12 (cid:0) | t − s − | H + | t − s + 1 | H (cid:1) − | t − s | H , (4d) E (cid:2) Y ( t ) (cid:3) =4 − H , (4e) E [ Y ( t ) Y ( s )] = − | t − s | H + 2 (cid:0) | t − s − | H + | t − s + 1 | H (cid:1) − (cid:0) | t − s − | H + | t − s + 2 | H (cid:1) . (4f)Eq. (4c) and (4d) can be obtained from Eq. (4a) and(4b) by substituting G ( t ) = B ( t + 1) − B ( t ); likewise,Eq. (4e) and (4f) follow from Eq. (4c) and (4d) by using Y ( t ) = G ( t + 1) − G ( t ). [22] H DfGnfGnfBm
FIG. 2. Fraction of turning points T for fBm, fGn, and DfGn(red lines). The simulations were performed for n = 2 (blackpoints). The dispersion in case of fBm increases for H (cid:38) . T ). For fBm, one has T fBm = 1 − π arcsin (cid:0) H − (cid:1) , (5)which is roughly equal to Tµ T from [14], where T is thenumber of turning points in a time series of length n [2, 23], and µ T = ( n −
2) is the expected value for whitenoise [24]. In general, E [ T ] = ( n − T . The plot ofEq. (5) is shown in Fig. 2, together with data simulatedas in [14] (scaled herein from T /µ T to T ). For an fGn: T fGn = 1 − π arcsin (cid:115) H − H +1 − H − , (6)and for the increments of fGn, i.e. DfGn: T DfGn = 1 − π arcsin (cid:115) (2 H + 2) − (2 · H ) H − · H +1 + 15 , (7)both of which are also shown in Fig. 2. III. ABBE VALUE, A The Abbe value of a time series { x i } ni =1 is defined ashalf the ratio of the mean-square successive difference tothe variance [14, 26–29]: A = n − n − (cid:80) i =1 ( x i +1 − x i ) n n (cid:80) i =1 ( x i − ¯ x ) . (8)It quantifies the smoothness (raggedness) of a time se-ries by comparing the sum of the squared differences be-tween two successive measurements with the variance ofthe whole time series. It decreases to zero for time se-ries displaying a high degree of smoothness, while thenormalization factor ensures that A tends to unity for awhite noise process [30]. It was proposed as a test forrandomness [31–33]. It is straightforward to show thatfor a 2-periodic time series, { a, b, a, b, . . . } , A = 2, whilefor a 3-periodic one, { a, b, c, a, b, c, . . . } , A = 3 / { x i } ni =1 to be a realization of length n of anfBm, B Hn , with Hurst exponent H . Thence, its incre-ments { g i } n − i =1 ≡ { x i +1 − x i } n − i =1 form an fGn, G Hn − ,with the same H . One can then express A as A fBm ( H, n ) = 12 var (cid:0) G Hn − (cid:1) var ( B Hn ) , (9)where the dependence on H and n is highlighted. Simi-larly, for an fGn A fGn ( H, n ) = 12 var (cid:0) Y Hn − (cid:1) var ( G Hn ) , (10)where Y Hn − is the corresponding DfGn, i.e. the incre-ments of fGn, { y i } n − i =1 ≡ { g i +1 − g i } n − i =1 . The goal is tocalculate the variances of B Hn , G Hn , and Y Hn . A. Variance of fGn
We start with var (cid:0) G Hn (cid:1) since it occurs in both Eq. (9)and (10). The same methodology that was used in [34]to calculate var (cid:0) B Hn (cid:1) is employed, i.e.var (cid:0) G Hn (cid:1) = nn − E (cid:104)(cid:0) G Hn − E (cid:2) G Hn (cid:3)(cid:1) (cid:105) = 1 n − E n − (cid:88) j =0 G j − n − (cid:80) k =0 G k n . (11)Developing the sum for a few values of n , and utilizngthe variance and covariance from Eq. (4c) and (4d), onecan observe a pattern emerging, that leads to a formula:var (cid:0) G Hn (cid:1) = n − n H − n − . (12)It yields lim n →∞ var (cid:0) G Hn (cid:1) → , (13) i.e. it asymptotically approaches Eq. (4c). However, forfinite n , lim H → var (cid:0) G Hn (cid:1) = 0. The plot of Eq. (12)for n = 2 is shown in Fig. 3. For values H (cid:38) . H v a r ( G n H ) log n = FIG. 3. Deviation of var (cid:0) G Hn (cid:1) from unity. B. Variance and Abbe value of fBm
The variance of the discrete, finite length B Hn is givenin [34] as var (cid:0) B Hn (cid:1) = 1 n ( n − n − (cid:88) i =1 ( n − i ) i H . (14)For H → i H →
1; then n − (cid:80) i =1 ( n − i ) = n ( n − , sovar (cid:0) B n (cid:1) = , hence, per Eq. (9) and given Eq. (12), A fBm ( H = 0 , n ) = n/ ( n − Math-ematica one can calculate the sum in Eq. (14) to be[35] n − (cid:88) i =1 ( n − i ) i H = ζ ( − H − , n ) − nζ ( − H, n )+ nζ ( − H ) − ζ ( − H − , (15)where ζ ( s ) is the Riemann ζ , and ζ ( s, n ) is the Hurwitz ζ [36]. Hence, taking into account Eq. (12), one can givea closed-form formula: A fBm ( H, n ) = n ( n − var (cid:0) G Hn − (cid:1) ζ ( − H − , n ) − nζ ( − H, n ) + nζ ( − H ) − ζ ( − H − , (16) n ( a ) H = H ( b ) n = FIG. 4. The ratio from Eq. (18) at (a) H = 0 for varying n ,and (b) its dependence on H for a set n = 100. which, since ζ (0 , n ) = − n , ζ ( − , n ) = (cid:0) − + n − n (cid:1) , ζ (0) = − , and ζ ( −
1) = − , yields A fBm ( H = 0 , n ) = n/ ( n − A fBm ( H, n ), first observe that since at H = 0 and for n (cid:29) ζ terms and the Riemann ζ terms is big, i.e. ζ ( − , n ) − nζ (0 , n ) | nζ (0) − ζ ( − | = − n | − n | (cid:29) , (17)and that for H ∈ (0 ,
1) this ratio is monotonically in-creasing (Fig. 4), therefore ζ ( − H − , n ) − nζ ( − H, n ) | nζ ( − H ) − ζ ( − H − | (cid:29) , (18)hence nζ ( − H ) − ζ ( − H −
1) is a negligible contribution,so it does not need to be taken into account.Let us express ζ ( s, n ) as a globally convergent Newtonseries, i.e. utilize the Hasse representation [37]: ζ ( s, n ) = 1 s − ∞ (cid:88) i =0 i + 1 i (cid:88) k =0 ( − k (cid:18) ik (cid:19) ( n + k ) − s , (19)valid for s (cid:54) = 1, n >
0. The first term of the outer sum,i.e. for i = 0 and at s = − H , is (cid:88) k =0 ( − k (cid:18) k (cid:19) ( n + k ) H = n H . (20)The second term, i.e. for i = 1, is12 (cid:88) k =0 ( − k (cid:18) k (cid:19) ( n + k ) H = 12 (cid:0) n H − ( n + 1) H (cid:1) . (21)And similarly for higher i . Therefore, for i = 0 one ob-tains an approximation ζ ( − H, n ) ≈ n H − H − , (22)and for i = 1: ζ ( − H, n ) ≈ − H − (cid:18) n H −
12 ( n + 1) H (cid:19) , (23)but since n (cid:29)
1, ( n + 1) ≈ n , so one also obtains ζ ( − H, n ) ≈ n H − H − . (24) For higher i , although given that i (cid:28) n , one has ( n + k ) ≈ n , hence yielding the same approximation. Indeed, i (cid:80) k =0 ( − k (cid:0) ik (cid:1) is the Kronecker δ , δ i , hence only the i = 0term survives. For s = − H − ζ ( − H − , n ) ≈ n H − H − . (25)Therefore, with these approximations one can write: ζ ( − H − , n ) − nζ ( − H, n ) ≈ n H − H − − nn H − H − n H H + 1)(2 H + 1) , (26)so that var (cid:0) B Hn (cid:1) ≈ n H H + 1)(2 H + 1) n ( n − , (27)and since ( n − ≈ n , one finally obtainsvar (cid:0) B Hn (cid:1) ≈ n H H + 1)(2 H + 1) (28)and A fBm ( H, n ) ≈ ( H + 1)(2 H + 1) n − H var (cid:0) G Hn − (cid:1) , (29)which also yields A fBm ( H = 0 , n ) = n/ ( n − (cid:0) G Hn − (cid:1) = 1,the simplest approximation is then obtained as A fBm ( H, n ) ≈ ( H + 1)(2 H + 1) n − H . (30)These various approximations of A fBm ( H, n ) are shownin Fig. 5. For an unreasonably small n = 8 one sees dis-crepancies between the curves for low and moderate H ,although they are rather small. The finite-size effects,introduced by var (cid:0) G Hn (cid:1) from Eq. (12), are significant athigher values of H . For a moderate n = 100, the curvesare indistinguishable for most of the range of H , with theregion of significant influence of var (cid:0) G Hn (cid:1) moved system-atically to higher H , and for even longer time series (e.g. n = 2 ) the consistency of all curves is only strength-ened. Therefore, the approximation from Eq. (30) is adecent one for time series with any reasonable length. C. Variance of DfGn and Abbe value of fGn
To obtain an expression for var (cid:0) Y Hn (cid:1) the samemethodology from Sect. III A is undertaken, i.e. Eq. (11)with G changed to Y is employed. Again, developing thesum for a few values of n , and utilizng the variance andcovariance from Eq. (4e) and Eq. (4f), one can observethe pattern depicted in Table I. Hence one can write exactexact with var ( G )= ( G ) approx. with var ( G )= H ( a ) n = exactexact with var ( G )= ( G ) approx. with var ( G )= - - - - - - H ( b ) n = exactexact with var ( G )= ( G ) approx. with var ( G )= - - - - - - H ( c ) n = FIG. 5. A ( H, n ) for (a) an extremely short ( n = 8) fBm time series, for (b) a moderate ( n = 100), and (c) a long one ( n = 2 ).The ”exact” line (red) corresponds to Eq. (16), ”exact with var( G ) = 1” to the same Eq. (16) but with var (cid:0) G Hn (cid:1) = 1, i.e. set toits asymptotic value (blue), the ”approx. with var( G )” denotes Eq. (29) with the expression for var (cid:0) G Hn (cid:1) included (green), and”approx. with var( G ) = 1” is the most simplified form from Eq. (30) (cyan). The lines overlap so strongly that for visualizationpurposes they are depicted with different thickness. var (cid:0) Y Hn (cid:1) = 2(2 n − − n H + ( n − H − n H + ( n + 1) H n ( n − . (31)The exact expression for the Abbe value, taking into ac- count Eq. (12), is thence A fGn ( H, n ) = n H + ( n − H + n ( n − (cid:0) − H (cid:1) − n − H + 2 − H n −
2) ( n − n H − ) (32) exactapprox. H v a r ( Y n H ) ( a ) log n = exactapprox. H v a r ( Y n H ) ( b ) log n = FIG. 6. var (cid:0) Y Hn (cid:1) for (a) a very short ( n = 32) time series,and (b) for a long one, n = 2 . The ”exact” line (black) isfor Eq. (31), and the ”approx.” line (red) denotes Eq. (33). In the limit n → ∞ , i.e. for any reasonable n (cid:29) (cid:0) Y Hn (cid:1) ≈ − H , (33)i.e. it asymptotically approaches Eq. (4e). Plots ofEq. (31) and (33) are shown in Fig. 6. For very short timeseries there is a certain deviation between the expres-sions, but for higher n the difference is invisible. There-fore, the asymptotic Eq. (33) is adequate in any practicalscenario. The Abbe value is thence A fGn ( H, n ) ≈ − H − var ( G Hn ) . (34)As discussed in Sect. III A, the variance of G Hn can insome instances be approximated by unity. Eq. (34) sim-plifies then to just A fGn ( H, n ) ≈ − H − , (35)independent on the length n of the time series.The asymptotic Eq. (35) ranges from 0 to 3 / H decreases from 1 to 0. The expression from Eq. (34)reaches its maximum of nn +1 at H = 0. The asymptoticminimum, as H →
1, is n − n ln 4ln n = n − n log n
4. For n =2 , these values are 1.49991 and 0.14285, respectively,in perfect agreement with Fig. 3 in [14]. IV. REPRESENTATION OF FBM AND FGN INTHE
A − T
PLANE
The
A−T plane is displayed in Fig. 7. The black pointscome from simulations [14]. In case of fBm [Fig. 7 (a) and
TABLE I. Expressions for var (cid:0) Y Hn (cid:1) for first few n , and theresulting general formula. n n ( n − (cid:0) Y Hn (cid:1) · − · H +1 H − · H +3 H · − · H +2 H − · H +4 H · − · H +3 H − · H +5 H · − · H +4 H − · H +6 H · − · H +5 H − · H +7 H · − · H +6 H − · H +8 H · − · H +7 H − · H +9 H · − · H +8 H − · H +10 H
10 2 · − · H +9 H − · H +11 H n − − n · H +( n − H − n H +( n + 1) H (b)], the red line employs the exact formula for A fBm fromEq. (16), while the cyan line depicts the approximationfrom Eq. (30). Eq. (5) describes T fBm .In case of fGn [Fig. 7 (c) and (d)], the red line corre-sponds to the exact Eq. (32) for A fGn , while the cyan lineto the asymptotic form from Eq. (35). Eq. (6) describes T fGn . Panels (b) and (d) of Fig. 7 employ a logarithmichorizontal axis to fully display the dependence T ( A ) atsmall values of A (i.e. high values of H ). The agreementbetween numerical simulations and the analytic descrip-tion is very good; also the approximations for A fGn ( H, n )work well in the
A − T plane. In particular, the approx-imation from Eq. (34) for fGn is as good as the exactEq. (32).
V. ARMA PROCESSES
The methodology from Sect. II and III is applicablealso to ARMA processes. Generalizing Eq. (8), similarlyas was done in Eq. (9) and (10), the Abbe value of aprocess X is A = 12 var ( dX )var ( X ) , (36)where dX denotes the increments of X . The fraction ofturning points is given by Eq. (1)–(3). It will be con-venient to express ρ ( d ) in terms of the autocorrelationfunction ρ d = E [ X ( t ) X ( t + d )]: ρ ( d ) = 2 ρ d − − ρ d − ρ d ) , (37)so that Eq. (1) becomes p ( d ) = 1 π arcsin (cid:18) (cid:114) − ρ d − ρ d (cid:19) . (38)This formula can be directly applied also to fGn andDfGn, but not to fBm which is nonstationary. For ARMA processes, ρ d and var ( X ) are easily ob-tainable [2, 3]. The variance of the differentiated process,var ( dX ), is calculated asvar ( dX ( t )) = E (cid:2) dX ( t ) (cid:3) = E (cid:104) ( X ( t + 1) − X ( t )) (cid:105) = E (cid:2) X ( t + 1) (cid:3) + E (cid:2) X ( t ) (cid:3) − E [ X ( t + 1) X ( t )] . (39) A. AR(1)
Consider a weakly stationary process X t = a X t − + ε t , − < a <
1, where ε t is a white-noise error term. Then ρ d = a d (40)for all d , thus T = 1 − π arcsin (cid:18) √ a (cid:19) , (41)reaching its minimum of 1 / a →
1, and its max-imum of 1 when a → − X ) = 11 − a (42)and var ( dX ) = 21 + a , (43)and thus A = 1 − a , (44)reaching its minimum of 0 when a →
1, and maximumof 2 when a → −
1. One can then write T explicitly asa function of A : T = 1 − π arcsin (cid:18) √ − A (cid:19) , (45)which is displayed in Fig. 8. The Ornstein-Uhlenbeck(OU) process, a continuous analog of AR(1), yields thesame Eq. (45), but restricted to A ∈ [0 ,
1] (see Ap-pendix ).
B. AR(2)
Consider a weakly stationary process X t = a X t − ++ a X t − + ε t , − < a < ∧ a < a ∧ a < − a .Then ρ d is given by a recurrent relation (the Yule-Walkerequations): ρ d = d = 0 a − a d = 1 a ρ d − + a ρ d − d ≥ ( a ) - - - - ( b ) ( c ) - - - - - - ( d ) FIG. 7. The
A − T plane for n = 2 : (a)–(b) fBm, (c)–(d) fGn. The red lines utilize the exact expressions for A [i.e. Eq. (16)and (32)], and cyan ones use the approximations [Eq. (30) for fBm and the asymptotic Eq. (35) for fGn]. Eq. (5) and (6)describe T fBm and T fGn , respectively. The horizontal and vertical dashed lines mark T = 2 / A = 1, respectively. Panels(b) and (d) display the same as (a) and (c), but with a logarithmic horizontal axis. The discrepancies in (b) at very low A aredue to the deviation of var (cid:0) G Hn (cid:1) from unity when H tends to 1. thus T = 1 − π arcsin (cid:18) √ a − a (cid:19) , (47)reaching its minimum of 0 when a → , a → −
1, andits maximum of 1 along the line a = a + 1.One then obtainsvar ( X ) = 1 − a ( a + 1)( a − a − a + a −
1) (48)and var ( dX ) = 21 + a + a a − a , (49)thus A = 1 + a a − , (50)reaching its minimum of 0 along the line a = 1 − a , andmaximum of 2 along the line a = 1 + a . One cannotwrite T explicitly as a function of A , as ( A , T ) is a two-dimensional region, depicted in Fig. 8. However, it ispossible to give a simple formula for the boundaries ofthis region. Note that for a given a , T is minimal when a → −
1. Hence by setting a = −
1, one can then solveEq. (50) for a , i.e. write a = 2(1 − A ), and insert this into Eq. (47) to obtain the lower boundary as T ( A ) lower boundary = 2 π arcsin (cid:32)(cid:114) A (cid:33) . (51)The upper boundary, T = 1, is attained when a → A = 0,obtained by sweeping a from − a =0. One can then observe, by repeating the above reason-ing for an arbitrary a , that the region of availability, S AR(2) , is formed as a continuum of curves parametrizedby a : S AR(2) = (cid:8) T ( A ) (cid:12)(cid:12) a (cid:9) = (cid:40) − π arcsin (cid:18) (cid:112) ( A − a − (cid:19) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) − < a < (cid:41) , (52)reducing to Eq. (51) when a → − C. MA(1)
Consider X t = b ε t − + ε t , weakly stationary for all b ∈ R . The autocorrelation function is ρ d = d = 0 b b d = 10 otherwise (53)thus T = 1 − π arcsin (cid:32) (cid:115) b − b + b (cid:33) , (54)reaching its minimum of 1 / b = 1, and its maximumof (2 /π ) arcsec √ ≈ .
73 at b = − X ) = 1 + b (55)and var ( dX ) = 2 (1 + b ( b − , (56)thus A = 1 − b + b b , (57)reaching its minimum of 1 / b = 1, and maximumof 3 / b = −
1. One can then write T explicitly as afunction of A : T = 1 − π arcsin (cid:18) √A (cid:19) , (58)which is displayed in Fig. 8. D. MA(2)
Consider X t = b ε t − + b ε t − + ε t , weakly stationaryfor all b , b ∈ R , for which ρ d = d = 0 b (1+ b )1+ b + b d = 1 b b + b d = 20 otherwise (59)thus T = 1 − π arcsin (cid:32) (cid:115) b + b ( b − b + b − b ( b + 1) (cid:33) , (60)reaching its minimum of 2 / b = (cid:0) √ (cid:1) , b = 1,and its maximum of 4 / b = (cid:0) − √ (cid:1) , b = 1.One then obtainsvar ( X ) = 1 + b + b (61) and var ( dX ) = 2 (cid:0) b + b − b ( b + 1) (cid:1) , (62)thus A = 1 − b ( b + 1)1 + b + b , (63)reaching its minimum of 1 − / √ b = √ , b = 1,and maximum of 1 + 1 / √ b = −√ , b = 1. Onecannot write T explicitly as a function of A , as ( A , T ) isa two-dimensional region, depicted in Fig. 8. E. ARMA(1,1)
Consider X t = a X t − + b ε t − + ε t , weakly stationaryfor − < a < b ∈ R , for which ρ d = (cid:40) d = 0 a d − ( a + b )( a b +1)1+2 a b + b otherwise (64)thus T = 1 − π arcsin (cid:32) (cid:115) (1 + a )(1 + a b + b )1 + b ( a + b − (cid:33) , (65)reaching its minimum of 1 / a → , b = 1, andits maximum of 1 when a → −
1, along the line b ∈ R .One then obtainsvar ( X ) = 1 + 2 a b + b − a (66)and var ( dX ) = 2 1 + b ( a + b − a + 1 , (67)thus A = 1 − a + b ( a − a b + b , (68)reaching its minimum of 0 when a →
1, along the line b ∈ R , and maximum of 2 when a → −
1, along the line b ∈ R . One cannot write T explicitly as a function of A ,as ( A , T ) is a two-dimensional region, depicted in Fig. 8.Note that Eq. (65) and (68) are invariant on changing b to 1 /b , so that in the context of geometrical depictionof the region of availability only the range − ≤ b ≤ | b | (cid:29) | b | (cid:28)
1. However, b = ± • when b = −
1, Eq. (68) yields a = 3 − A , whichfulfills the stationarity condition, a ∈ ( − , A ∈ (1 , T ( b = −
1) = 1 − π arcsin (cid:32) (cid:114) − A − A (cid:33) ; (69) • likewise, when b = 1, one obtains a = 1 − A ,which leads to A ∈ (0 , T ( b = 1) = 1 − π arcsin (cid:18) √ − A (cid:19) . (70)In particular, b = 0 reduces an ARMA(1,1) processto an AR(1) one, reproducing respective formulae fromSect. V A.Similarly as was done in Sect. V B for AR(2) processes,the region of availabity for ARMA(1,1) can be describedas a continuum of curves, S ARMA(1 , = (cid:8) T ( A ) (cid:12)(cid:12) b (cid:9) ,parametrized by b : one needs to solve Eq. (68) for a andinsert the solution into Eq. (65). The resulting formula,easy to derive but of a quite complicated and noninfor-mative form, is not displayed herein. VI. APPLICATIONSA. Bacterial cytoplasm
The two-dimensional motion [ x ( t ) , y ( t )] of individualmRNA molecules inside live Escherichia coli bacteriawere tracked in [38]. It was found that they follow anoma-lous diffusion, with
H < .
5, confirmed by other methodsas well [39]. Herein, the time series x and y are treatedseparately. Results for the 27 tracks are displayed inFig. 9 for the x -axis. Similar outcomes were obtained forthe y -axis. The locations in the A−T plane are in agree-ment with an fBm description, and the values extractedusing Eq. (5) and (30), yielding 0 . (cid:46) H (cid:46) .
5, are con-sistent with each other, and confirm that the observedprocess is indeed subdiffusive.
B. Zeros of the Riemann ζ The first 10 + 1 nontrivial zeros, + i γ n , of the Rie-mann ζ function [40] were retrieved [41]. Normalizedspacings between consecutive zeros were computed as [42] δ n = γ n +1 − γ n π ln (cid:16) γ n π (cid:17) . (71)The location of this sequence in the A − T plane is(1 . , . H = 0 .
19. In compar-ison, the discrete wavelet transform (DWT) method [14]returns H = 0 .
06, hence also strongly implying
H < . δ n is not normal, butrather follows the distribution of the Gaussian unitaryensemble (GUE) according to the GUE hypothesis [43],so this sequence is not a Gaussian process, strictly speak-ing. Note, however, that a location in the A − T planecan be computed for any type of data; in this case, dueto A > T > /
3, one gets a clear information thatthe series is—in a sense—(much) more noisy than regularwhite noise. ( a ) AR ( ) ARMA ( ) MA ( ) MA ( ) AR ( ) ( b ) ( c ) ( d ) ( e ) ( f ) FIG. 8. (a) Regions of the
A − T plane available for ARMAprocesses: AR(1) (red line), MA(1) (blue line), MA(2) (darkergray region), AR(2) (light blue in the background), andARMA(1,1) (lighter gray region). The black dotted line in therange
A ∈ (0 ,
1) is Eq. (70) for ARMA(1,1) with b = 1, andthe black dashed line in A ∈ (1 ,
2) is Eq. (69) for ARMA(1,1)with b = −
1. These two lines mark the lower boundary ofthe region of availability for the ARMA(1,1) process. Theyellow dot at (1 , /
3) denotes white noise. In the other pan-els, these regions are shown together with locations of 10 simulated processes: (b) AR(1), (c) MA(1), (d) AR(2), (e)MA(2), and (f) ARMA(1,1). In the simulations, the AR co-efficients were drawn uniformly from the respective regionsfulfilling the stationarity conditions, and the MA coefficientswere drawn uniformly from ( − , C. Sunspot numbers
The curently available from the World Data CenterSunspot Index and Long-term Solar Observations [44]sample of 3250 monthly SNN are described by H = 0 . ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ( a ) ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● ●● H from H f r o m ( b ) x - axis motion FIG. 9. (a) Locations in the
A − T plane of the 27 tracksof x -axis motion of mRNA molecules inside E. coli . The redlines are the fBm lines, corresponding to the length of timeseries 140 ≤ n ≤ n .(b) Estimated H values, obtained from the formulae for A and T , i.e. Eq. (30) and (5), respectively. In both panels,the size of the point is proportional to n . The gray crossessymbolize the (unweighted) means and standard deviations ofthe displayed locations. ζ S S N ChirikovWIG
FIG. 10. Locations of some real-world and generated timeseries. The black lines correspond to fBm and fGn with n =2 . See text for details. tion in the A − T plane is (0 . , . H ≈ . − .
27. This is in agree-ment with some works [46] that also compute
H < . D. Chaos
The Chirikov standard map for a chaotic state (in anunbounded setting) was examined in [14]. Its location inthe
A − T plane is (0 . , . H ≈ .
5, which is in perfect agreement with the estimate withDWT method, yielding H = 0 .
48. This means that, inthe context of long-term memory, it is uncorrelated, andacts like Brownian motion. ( a ) ( b ) r ( d ) rx ( c ) - - m L E ( e ) FIG. 11. Relations between A , T , r , and mLE for the logisticmap. To further illustrate chaotic behavior in the
A − T plane, the logistic map x i +1 = rx i (1 − x i ) is consid-ered. For r ∈ [3 . , . r = 0 . n = 10 were generated and their ( A , T )locations, as well as the maximal Lyapunov exponents(mLEs), were computed. The dependencies of A and T on r , the bifurcation diagram, the A − T plane, andthe relation between A and mLE, are shown in Fig. 11.When chaos is most developed ( r = 4), the trajectoriesapproach the point ( A , T ) = (1 , / r , a gradual decrease in T oc-curs, with wells in periodic windows [Fig. 11(b)]. Notethat T = 1 up to r (cid:46) .
7, because apparently the orbits,even when chaotic, are strictly alternating. A is moresensitive to changes in dynamics, as A = 2 for period-2orbits (see the beginning of Sect. III) before the bifur-cation at r = 1 + √ ≈ .
45, and then systematicallydecreases in the period-4 window before the next bifurca-tion at r ≈ .
54 [Fig. 11(a)]. There are also shallow wellsat periodic windows within the chaotic zone. The path inthe
A−T plane is jagged, and various changes in dynam-ics are manifested, corresponding to changes in mLEs[Fig. 11(e)] and in the bifurcation diagram [Fig. 11(c)].A tight, positive correlation between mLE and H forthe Chirikov map was discovered [13] (see also [47]).Differences between (quasi)periodic and chaotic systemswere observed in the coarse-grained sequences in the A − T plane in case of a sinusoidally driven thermo-stat [18]. These connections between chaos, H , and the A − T plane are nontrivial, and require further research.Describing a universal (if existent) behavior of chaotic1 FIG. 12. Time evolution of the WIG. The green dot marksthe start point, and the red dot denotes the end point. Notethe scale of the horizontal axis. See text for details. systems in the
A − T plane in an interesting perspec-tive. The possibility of differentiating between chaoticand stochastic time series is a hopefully attainable appli-cation.
E. Markets
The efficient market hypothesis [48] is a key concept infinance. If a market exhibits H (cid:54) = 0 .
5, then it might allowarbitrage. A classification of developed, emerging andfrontier stock markets in the
A − T plane was performedin [17], and showed that these three categories of marketsoccupy distinct regions of the ( A , T ) space. Hence, the A − T plane is a useful tool for classifying data withdifferent underlying dynamics.The time evolution of H of many stock markets showsit is oscillating around H = 0 . A and T evolve, and hence trace theirvalues in the A−T plane at different times. As an exam-ple consider the Warsaw Stock Exchange Index (Warsza-wski Indeks Gie(cid:32)ldowy, WIG) [49]. Its location in the
A − T plane, right on the fBm line, is depicted in Fig. 10,and the estimates of H , obtained by solving Eq. (5) and(30) with n = 2166, yield H ≈ . − .
59. The DWTmethod returns H = 0 .
48. The time evolution is depictedin Fig. 12. It is obtained by partitioning the whole timeseries into overlapping segments of size n/
2, advancingeach segment by one point. Throughout its history, theWIG remains confined in a small region of the
A − T plane.
F. Active galactic nuclei
A core focus in astronomy is the investigation of ap-parent variability of various celestial objects such as as-teroids, stars, or galaxies. Within the latter, of partic-ular interest are active galactic nuclei (AGNs), furtherdivided into several types [50]. Among them, blazarsare peculiar AGNs pointing their relativistic jets towardthe Earth. Blazars are commonly divided further intotwo subgroups, i.e. flat spectrum radio quasars (FSRQs) and BL Lacertae (BL Lac) objects, based on character-istics visible in their optical spectra. In a recent study[51] it was found that faint FSRQ and BL Lac candi-dates located behind the Magellanic Clouds are clearlyseparated in the
A − T plane, with means of
A ≈ . A ≈ .
7, respectively. This differentiation, basedsolely on the temporal data in the form of light curves,is another proof that employing the
A − T plane as aclassification tool is a promising approach.
G. Other
The worked-out examples from Sect. VI A–VI F do notexhaust the possible applications of the
A − T plane inregard to constraining the value of H , or classifying timeseries. Some other interesting instances include, but arenot restricted to: ferro- and paramagnetic states of theHeisenberg model that exhibit H ∼ H ∼ .
5, re-spectively [52], and should be easily distinguishable inthe
A − T plane; a photonic integrated circuit yields0 . (cid:46) H (cid:46) . H > .
5, suggesting the accretion is driven by mag-netic fields [54]; football matches can follow the rules offBm with H ∼ . Nitzschia sp. diatoms [57]; Solar wind pro-ton density fluctuations are characterized by H ∼ . H > . A − T plane as well [17].
VII. SUMMARY AND OPEN QUESTIONS
Exact analytical descriptions for the locations of fBmand fGn in the
A − T plane were derived. Working ap-proximations were also obtained in the following forms: (cid:40) A fBm ( H, n ) = ( H + 1)(2 H + 1) n − H T fBm ( H ) = 1 − π arcsin (cid:0) H − (cid:1) (72)and (cid:40) A fGn ( H ) = 2 − H − T fGn ( H ) = 1 − π arcsin (cid:16) (cid:113) H − H +1 − H − (cid:17) (73)and were demonstrated to be adequate for time serieswith any reasonable length n . This allows to classify timeseries of any length, respective to fBm and fGn, withoutrelying on time-consuming numerical simulations. Theseanalyses add to the theoretical results regarding fBmand fGn [60, 61]. The same methodology was applied2to ARMA( p, q ) processes. Analytic descriptions of theavailable regions of the A − T plane were derived andillustrated for p + q ≤ A is required, as it has been rarelyutilized, with some recent, nonextensive examples in as-tronomy [29, 51, 62–64] (but see also [65]). The interre-lations between ordinal patterns, persistence, and chaos[8, 66] are linked even tighter with the bidimensionalscheme of the A − T plane. The presented methodol-ogy is valid for any Gaussian process, but it should beemphasized that the locations ( A , T ) can be computedfor arbitrary time series, serving e.g. as classification orclustering methods for empirical data.Naturally, a question about A − T representations ofother stochastic processes arises. Examples include thefollowing: • Representation of colored noise, i.e. power spectraldensities (PSDs) of the form 1 /f β . fBm can be as-sociated with β ∈ (1 , β ∈ ( − , A − T plane, for any β , is an interesting and challengingproblem. • In some fields, e.g. in astronomy, the observedsignals often yield PSDs of the form 1 /f β + C ,where C is the so called Poisson noise level, comingfrom the statistical noise due to uncertainties in themeasurements; above a certain frequency, the PSDtransitions from a power law to white noise. A rep-resentation of such processes in the A − T plane iscrucial in classifying light curves of several sources,e.g. AGN. • Other continuous-time models, e.g. continuousARMA [67], have been developed. The simplestin this family is the OU process, which is a contin-uous analog of the AR(1) process, with the same
A − T representation. Introducing long-term mem-ory leads to continuous autoregressive fractionally integrated moving average models [68].
ACKNOWLEDGMENTS
The author thanks Ido Golding for providing the dataof the mRNA motion inside the bacteria. Support by thePolish National Science Center through OPUS Grant No.2017/25/B/ST9/01208 is acknowledged.
Appendix: Ornstein-Uhlenbeck process
The OU process with mean µ is given by the stochasticdifferential equationd x t = θ ( µ − x t ) d t + σ d ε t , (A.1)where θ > σ >
0. The correlation function for lag d is ρ d = exp ( − θd ) , (A.2)which via Eq. (38) and (3) leads to T OU ( θ ) = 1 − π arcsin (cid:18) (cid:112) − θ ) (cid:19) , (A.3)reaching its minimum of 1 / θ = 0, and its maximumof 2 / θ → ∞ .The covariance function E [ x t x s ] = σ θ exp ( − θ | t − s | ) (A.4)gives via Eq. (36) and (39) and on simplification A OU ( θ ) = 1 + sinh θ − cosh θ, (A.5)reaching its minimum of 0 for θ = 0, and its maximumof 1 when θ → ∞ . Eq. (A.5) can be solved for θ andinserted into Eq. (A.3), which yields Eq. (45), valid for A ∈ [0 , [1] J. Beran, Statistics for Long Memory Processes (Chap-man and Hall, New York, 1994).[2] P. J. Brockwell and R. A. Davis,
Time Series: Theoryand Methods, 2nd ed. (Springer-Verlag New York, 1996).[3] P. J. Brockwell and R. A. Davis,
Introduction to Time Se-ries and Forecasting, 2nd ed. (Springer-Verlag New York,2002).[4] J. Beran, Y. Feng, S. Ghosh, and R. Kulik,
Long MemoryProcesses: Probabilistic Properties and Statistical Meth-ods (Springer-Verlag, Berlin Heidelberg, 2013).[5] R. Hegger, H. Kantz, and T. Schreiber, Chaos , 413(1999).[6] C. Bandt and B. Pompe, Physical Review Letters ,174102 (2002).[7] J. E. Maggs and G. J. Morales, Plasma Physics and Con- trolled Fusion , 085015 (2013).[8] H. V. Ribeiro, M. Jauregui, L. Zunino, and E. K. Lenzi,Physical Review E , 062106 (2017).[9] I. Simonsen, A. Hansen, and O. M. Nes, Physical ReviewE , 2779 (1998), cond-mat/9707153.[10] A. Carbone, G. Castelli, and H. E. Stanley, PhysicaA: Statistical Mechanics and its Applications , 267(2004).[11] S. Arianos and A. Carbone, Physica A: Statistical Me-chanics and its Applications , 9 (2007).[12] L. Lacasa and R. Toral, Physical Review E , 036120(2010).[13] M. Tarnopolski, Physica A: Statistical Mechanics and itsApplications , 834 (2018).[14] M. Tarnopolski, Physica A: Statistical Mechanics and its Applications , 662 (2016).[15] C.-K. Peng, S. V. Buldyrev, S. Havlin, M. Simons, H. E.Stanley, and A. L. Goldberger, Physical Review E ,1685 (1994).[16] C.-K. Peng, S. Havlin, H. E. Stanley, and A. L. Gold-berger, Chaos , 82 (1995).[17] L. Zunino, F. Olivares, A. F. Bariviera, and O. A. Rosso,Physics Letters A , 1021 (2017).[18] Y. Zhao and G. J. Morales, Physical Review E , 022213(2018).[19] L. Zunino, F. Olivares, F. Scholkmann, and O. A. Rosso,Physics Letters A , 1883 (2017).[20] F. Traversaro, F. O. Redelico, M. R. Risk, A. C. Frery,and O. A. Rosso, Chaos , 075502 (2018).[21] C. Bandt and F. Shiha, Journal of Time Series Analysis , 646 (2007).[22] In the derivation of T fGn and T DfGn one can actually em-ploy only the relations for fBm; the calculations wouldthen become quite burdensome, though.[23] M. Kendall and A. Stuart,
The advanced theory of statis-tics (London: Griffin, 3rd ed., 1973).[24] The first and last points cannot form a turning point,hence the subtraction of 2 in µ T .[25] M. Sinn and K. Keller, arXiv e-prints , arXiv:0801.1598(2008), arXiv:0801.1598 [math.PR].[26] J. von Neumann, The Annals of Mathematical Statistics , 153 (1941).[27] J. von Neumann, The Annals of Mathematical Statistics , 367 (1941).[28] M. G. Kendall, Biometrika , 369 (1971).[29] N. Mowlavi, Astronomy & Astrophysics , A78 (2014).[30] J. D. Williams, The Annals of Mathematical Statistics , 239 (1941).[31] C. Bingham and L. S. Nelson, Technometrics , 285(1981).[32] R. Bartels, Journal of the American Statistical Associa-tion , 40 (1982).[33] A. Mateus and F. Caeiro, “Exact and approximate prob-abilities for the null distribution of bartels randomnesstest,” in Recent Studies on Risk Analysis and Statisti-cal Modeling , edited by T. A. Oliveira, C. P. Kitsos,A. Oliveira, and L. Grilo (Springer International Pub-lishing, Cham, 2018) pp. 227–240.[34] D. Deligni`eres, Mathematical Problems in Engineering , 485623 (2015).[35] See also .[36] T. M. Apostol,
Introduction to Analytic Number Theory (Springer Science+Business Media New York, 1976).[37] H. Hasse, Mathematische Zeitschrift , 458 (1930).[38] I. Golding and E. C. Cox, Physical Review Letters ,098102 (2006).[39] M. Magdziarz, A. Weron, K. Burnecki, and J. Klafter,Physical Review Letters , 180602 (2009).[40] D. J. Platt, Mathematics of Computation , 1521(2015).[41] .[42] A. M. Odlyzko, Mathematics of Computation , 273(1987).[43] H. L. Montgomery, in Analytic number theory , Proc.Sympos. Pure Math., XXIV, Providence, R.I.: AmericanMathematical Society, Vol. 24, edited by H. G. Diamond(1973) pp. 181–193.[44] SILSO World Data Center, Royal Observatory of Bel-gium, Brussels, . [45] Similar value can be obtained with data in [14] when thetwo outliers on the left of Fig. 6b therein are discarded.[46] M. S. Movahed, G. R. Jafari, F. Ghasemi, S. Rahvar,and M. R. R. Tabar, Journal of Statistical Mechanics:Theory and Experiment , P02003 (2006).[47] W.-H. Steeb and E. C. Andrieu, Zeitschrift Natur-forschung Teil A , 252 (2005).[48] R. N. Mantegna and H. E. Stanley, An Introduction toEconophysics: Correlations and Complexity in Finance (Cambridge University Press, New York, 2000).[49] , accessed on November 10,2019.[50] C. M. Urry and P. Padovani, Publications of the Astro-nomical Society of the Pacific , 803 (1995).[51] N. ˙Zywucka, M. Tarnopolski, M. B¨ottcher, (cid:32)L. Stawarz,and V. Marchenko, arXiv e-prints , arXiv:1912.03530(2019), arXiv:1912.03530 [astro-ph.GA].[52] Z. Zhang, O. G. Mouritsen, and M. J. Zuckermann,Physical Review E , R2327 (1993).[53] K. E. Chlouverakis, A. Argyris, A. Bogris, andD. Syvridis, Physical Review E , 066215 (2008).[54] G. Anzolin, F. Tamburini, D. de Martino, and A. Bian-chini, Astronomy and Astrophysics , A69 (2010),arXiv:1005.5319 [astro-ph.SR].[55] A. Kijima, K. Yokoyama, H. Shima, and Y. Ya-mamoto, European Physical Journal B , 41 (2014),arXiv:1402.0912 [physics.pop-ph].[56] N. Makarava, S. Menz, M. Theves, W. Huisinga, C. Beta,and M. Holschneider, Physical Review E , 042703(2014).[57] J. S. Murgu´ıa, H. C. Rosu, A. Jimenez, B. Guti´errez-Medina, and J. V. Garc´ıa-Meza, Physica A: Statis-tical Mechanics and its Applications , 176 (2015),arXiv:1410.3135 [q-bio.QM].[58] F. Carbone, L. Sorriso-Valvo, T. Alberti, F. Lepreti,C. H. K. Chen, Z. Nˇemeˇcek, and J. ˇSafr´ankov´a, As-trophysical Journal , 27 (2018), arXiv:1804.02169[physics.space-ph].[59] C. Witton, S. V. Sergeyev, E. G. Turitsyna, P. L. Fur-long, S. Seri, M. Brookes, and S. K. Turitsyn, Journalof Neural Engineering , 056019 (2019).[60] L. Zunino, D. G. P´erez, M. T. Mart´ın, M. Garavaglia,A. Plastino, and O. A. Rosso, Physics Letters A ,4768 (2008).[61] T. Sadhu, M. Delorme, and K. J. Wiese, Phys. Rev.Lett. , 040603 (2018).[62] M.-S. Shin, M. Sekora, and Y.-I. Byun, Monthly Noticesof the Royal Astronomical Society , 1897 (2009).[63] K. V. Sokolovsky, P. Gavras, A. Karampelas, S. V.Antipin, I. Bellas-Velidis, P. Benni, A. Z. Bonanos,A. Y. Burdanov, S. Derlopa, D. Hatzidimitriou, A. D.Khokhryakova, D. M. Kolesnikova, S. A. Korotkiy, E. G.Lapukhin, M. I. Moretti, A. A. Popov, E. Pouliasis,N. N. Samus, Z. Spetsieri, S. A. Veselkov, K. V. Volkov,M. Yang, and A. M. Zubareva, Monthly Notices of theRoyal Astronomical Society , 274 (2017).[64] M. F. P´erez-Ortiz, A. Garc´ıa-Varela, A. J. Quiroz, B. E.Sabogal, and J. Hern´andez, Astronomy & Astrophysics , A123 (2017).[65] J. Lafler and T. D. Kinman, Astrophysical Journal Sup-plement , 216 (1965).[66] M. Zanin, L. Zunino, O. A. Rosso, and D. Papo, Entropy , 1553 (2012).[67] P. J. Brockwell, Annals of the Institute of Statistical Mathematics , 647 (2014).[68] H. Tsai, Bernoulli15