Modeling statistics of the natural aggregation structures and processes with the solution of generalized logistic equation
11 Modeling kinetics of the natural aggregation structures and processes with the solution of generalized logistic equation
Lev A Maslov , Vladimir I Chebotarev University of Northern Colorado, Aims College, Greeley, CO, USA, [email protected] Computer Center RAS, Khabarovsk, Russia, [email protected]
Abstract
The generalized logistic equation is derived to model kinetics and statistics of natural processes such as earthquakes, forest fires, floods, landslides, and many others. This equation has the form AANANsdAAdN q )())(1()( , q>0 and A>0 is the size of an element of a structure, and . The function ( ) N A can be understood as an approximation to the number of elements of size A. The general solution of this equation for q =1 is a product of an increasing bounded function and power-law function with stretched exponential cut-off; the power-law distribution is asymptotically nested in the stretched exponential distribution. The relation with Tsallis non-extensive statistics is demonstrated by solving the generalized logistic equation for q >0. In the case 0< q <1 the equation models super-additive, and the case q >1 it models sub-additive structures. The Gutenberg-Richter (G-R) formula results from interpretation of empirical data as a straight line in the area of stretched exponent with small α . The solution is applied for modeling distribution of foreshocks and aftershocks in the regions of Napa Valley 2014, and Sumatra 2004 earthquakes fitting the observed data well, both qualitatively and quantitatively. Keywords: generalized logistic equation, power-law distribution, stretched exponential distribution, non-extensive statistics, hierarchical aggregation processes, earthquakes
Introduction
A probability distribution function with a power-law dependence in terms of a process or object sizes is the very common statistical feature of many natural processes and systems (Clauset, et al., 2009; Sornette, 2009). To account quantitatively for fat tail distributions, J. Laherr`ere and D. Sornette (1998), proposed the combination of power-law function with stretched exponential cut-off. Statistical studies (Clauset, et al., 2009) show that in some cases the power law is ruled out, but in the other cases, conjectures are consistent with a power law-distribution. However, once the evidence for a power law distribution is demonstrated, the uncertainty – what physical process is behind this distribution – remains. A proper solution of kinetic equation modeling the process should be complementary to the results of statistical research. Examples of this approach are fracture kinetics equations developed to explain clustering of earthquake sequences (Newman, Knopoff, 1983) and fracture kinetics equations developed by Z. Czechowski (1991). The authors of (Ovchinnikov, Eliseev, 1981) use the logistic equation to describe the process of corrosion fracture registering signals of acoustic emission. J. D. Achter and C. T. Webb (2006) and J. Goldenberg with coauthors (2000) use logistic equation to describe aggregation processes in biology and in marketing. F. Vallianatos and coauthors (Vallianatos, Michas, Papadakis, 2015) presented kinetic equation to model evolution of structures in the frame of non-extensive physics statistics. However, the main lack of equations mentioned above is that they do not have solutions with the fractional power-law exponent. Publications (Gabrielov, et al., 1999; Yakovlev, et.al., 2005) demonstrate the other approach to modeling the “power-law” structures and processes. The authors developed the formal model of hierarchical aggregation in the inverse cascade process. In this model, the aggregation of elements, similar to that of aggregation of clusters in percolation, produces a fractal hierarchical structure. These authors showed that the distribution of elements in their model approaches asymptotically the power-law distribution with fractional exponent α : , 0, i i N A i >>
1. In the current work we: a) combine two approaches to derive kinetic equation to model hierarchical aggregation structures and processes, b) analyze solution of this equation and to show its relation with power-law function with stretched exponential cut-off, and c) demonstrate relation of general solution with Tsallis non-extensive statistics.
Kinetic equation for q = 1 Let
N(A) is a differentiable function of a real variable A ≥ 0. This variable represents the “size” of an element in a structure (it can be an area, mass, magnitude of earthquake, etc.).
The function ( )
N A can be understood as an approximation to the number of all elements with the size less than A. Considering A as a “natural variable” similar to that used in thermodynamics, we suggest that the rate of growth of N with A is proportional to product ( ) N A A : AANsdAAdN gain )()( , q =1, s > , (1) where A is the size of an element of a structure defined by its characteristic diameter and 0≤ α ≤1. A represents the “fractal correction” to the area of an element of a structure. The rate of loss of elements due to their coalescence is AANNLdAAdN loss )()()( (2) Function
L(N) can be derived analytically for some relatively simple models of aggregation. In this work, we suggest that
NkNL )( , k is a positive constant. Thus, the balance equation is: ( ) [ ( )] ( ) dN A s k N A N A AdA (3) Normalizing N by s/k , s > 0 , k > 0 , we arrive at the equation (1 ) dN s N N AdA , (4) with skNN / , N , and N for A , N , i.e. )( AN is the cumulative distribution function of the random size of an element of a structure. Equation (4) is presented first in (Maslov, Anokhin, 2012). If 0≤ α <1, then equation (4) can be written in the form of a classical logistic equation: NNsdxNd ~)~1(~ (5) with )()(~ )1/(1 xNxN , Ax and )1/( ss . If α =1 equation (4) can be written in the form (5) with AxeNxN x ln),()(~ , and ss . In the rest part of the work we omit the bar sign over N in (4). The boundary condition for (4) is )(lim)0( NANN A . (6) The solution to the problem (4), (6) for is: )0( 0,1exp1)( NN AAsCAN (7) with the derivative s As A dN eC s AdA Ce , A >0. (8) Function (8) is not a density function because NNNAdNdAdAAdN , (9) but for sufficiently large C , dAdN can be considered approximately as a probability density function on half-axis A . Note that the function
0, 0,( ) ( ), 0
AF A N A A is the cumulative distribution function with the jump N at A . For α = 0 in (7) we have )0( ,0,)1()( NN AeCAN sA (10) with NNC . Note that )1( sA eC is a distribution function of a logistic law if to consider it on the axis ),( . For α = 1 the solution of (4) is bNbN AACAN b s (11) which is the log-logistic distribution function with b bs NNbC . The complementary to N(A) function in (7) :1 * NN .1)0( ,0,1)(
0* 11* NN AeCeCAN
AsAs (12) Analysis and interpretation of the solution
The function (8) can be written in the form: )()()(
AEAPALdAdN , (13) where ( ) (1 ) s A C sL A Ce , ( ) P A A and ( ) s A E A e (14) L(A) is the increasing bounded function varying from
Cs/(1+C) for A = 0 until Cs for A . The product ( ) ( ) s A P A E A A e (15) is the power function with stretched exponential cut-off. These two functions define two stages of a structure building. In the first stage, the elements (e.g. small fractures in the Earth’s crust, microcracks in rocks and metals, etc.) grow in number and size independently because of their relatively low density in a substance. As the same mechanism is applied to each one of them, this process is scale invariant and is described by the power function ( ) P A A (Laherr`ere, Sornette 1998). The second stage is interpreted conceptually as interaction of elements, their coalescence and aggregation. This process is inverse to the multiplicative process described by stretched exponential function (Sornette, 1998; Naumis, Cocho, 2007). As these two processes are formally undistinguishable, same formula is applied to the second stage. Thus, one solution of generalized logistic equation models the process from birth and growth of elements to their coalescence and formation the whole structure. Figure 1 shows graphs of the function AsAs eCeCAN , C = 1000, s = 2.3, α = 0.1 (16) altogether with data from regions of Sumatra 2004 and Napa Valley 2014 earthquakes. Figure 1 a. Theoretical function )(log *10 AN and real data )(log *10 mN from Sumatra 2004 earthquake region; b. Theoretical function )(log *10 AN and real data )(log *10 mN from Napa Valley 2014 earthquake region; N * (m) number of earthquakes with magnitudes M≥m . Different intervals of earthquake magnitudes, time and depth intervals are tested; c. The “roll-off” part of )(log *10 AN for A ≤1. The parameters C, s, and α in (7) for theoretical curves in Figures 1a, and b are found by solving the system of three equations (7) with points ( A i , N i ) , i =1,2,3, taken arbitrary from the observed data. The “roll-off” part of theoretical graphs and of the real data plots displayed on Figure 1 is not due to lack of data observed as it is stated in a number of publications, e.g. (Christensen, et al., 2002). It is the innate property of hierarchical aggregation processes. This part of a curve reflects the independent growth of elements in number and size, and because of their relatively low density in a substance it is described by the power function in (8) and (15). Earthquakes is one of examples of processes of that kind. Interpretation of the Gutenberg-Richter formula
In this section, we substitute a variable A in formula (12) by traditional in seismology letter m , magnitude of earthquake. The empirical G-R formula is: *10 log ( ) N m a b m (17) N * is the number of all earthquakes of magnitudes M ≥ m ; a and b are constants. Using an expansion of (12) in powers of exp( / (1 ) ) C s m for large m , we obtain )( ms eCmN , (18) and * 110 10 10 log ( ) log log1 sN m C m e . (19) For bsCa ,0,log we arrive at formula (17). Expression (19) is the analogy of Gutenberg-Richter empirical formula with b-value log1 sb e . (20) In (19) α defines the concavity of the part of the )(log *10 mN curve dominated by the stretched exponent E(A) . For α ≈0 this part of )(log *10 mN is close to a straight line and interpreted in empirical data as a power dependence between magnitude m and total number of earthquakes with magnitudes mM . Y. Malevergne and coauthors (Malevergne, Pisarenko, Sornette, 2005; Malevergne, Sornette, 2006) in the course of statistical analysis of empirical data showed that the power law family is asymptotically nested in the stretched exponential family. That is the exponent b of the power law is obtained from the parameters of the stretched exponential fit ))/(exp( dm when . In the present work this result follows directly from solution (12) for large m and . Relation with Tsallis non-extensive statistics
For the unconstrained growth of a number of elements: ANsdAdN qgain , q > 0. (21) In Tsallis non-extensive statistical physics parameter q is interpreted as an “entropic index”, corresponding if q =1 to Boltzmann-Gibbs distribution (Vallianatos, et al., 2015). For 0< q<1 the general solution of (21) is 0 q AsqCAN , (0, ) A , (22) and C is an arbitrary positive constant. If q >1, then the equation (22) is valid for sqCA . If q =1 and , the general solution of equation (21) is )( As eCAN . (23) For q =1 and α =1 the solution of (21) is a power-law function ( ) s N A C A . (24) If 0< q<1 , and α =1, the general solution of equation (21) is ( ) (1 ) ln q N A C s q A , (25) )1( qs C eA , and if q >1, )1( qs C eA . In case of a constrained growth of a structure due to loss of elements in the coalescence and aggregation, the main equation has the form AANANsdAAdN q )())(1()( , q > 0. (26) The general solution of (26) for q , 0≤ α< )(1)(1 )(ln1 CANnANANsA qn n (27) provided the base to exponent 1/(1- α ) is not negative. Solution of the same equation for 0< q <1, α ≠1 is 1 )(1 n qn Cqn ANsA (28) provided the base to exponent 1/(1- α ) is not negative. Figure 2 shows graphs of functions N(A) from (27) and (28) for C =30, C = , α= , s= , q , q =1, and q =2. These graphs are plotted by reflecting graphs A=N(A) about the line
N(A)=A . Figure 2 Solution of (28) for q =1/2 and solutions of (27) for q =1, q =2. These graphs demonstrate distinction of aggregation processes in structures with different properties. If to take (8) as a reference case corresponding to Boltzmann-Gibbs extensive statistics, q =1, then solution of (28) represents the “super-additive” aggregation with the entropy of a new system greater than entropies of components. Correspondingly, the solution of (27) 2 represents the “sub-additive” aggregation with the entropy of a new system less than entropies of components before the aggregation (Vallianatos, et al., 2015). Discussion and Conclusion
The generalized logistic equation to model kinetics of natural processes such as earthquakes, forest fires, floods, landslides, and many others belonging to the class of hierarchical aggregation processes is derived. The solution of this equation is a product of an increasing bounded function and power-law function with stretched exponential cut-off. A. Clauset and coauthors (2009), using maximum likelihood fitting methods, tested twenty-four real world data sets including earthquake data. They demonstrated that the power-law with stretched exponential cut-off is “clearly favored over the pure power-law”. Solution of generalized logistic equation confirms results of statistical analysis made by A. Clauset and coauthors (2009) and gives physics-based qualitative and quantitative approaches to the modeling of the hierarchical aggregation structures and processes. From solution (12) of this equation if follows that for large m and it asymptotically approaches the power function. This confirms result made in the course of statistical analysis of empirical data (Malevergne, Pisarenko, Sornette, 2005; Malevergne, Sornette, 2006). The exponent α defines the concavity of the part of the )(log *10 mN real data curve dominated by the stretched exponent. See example on Figure 1b. Gutenberg-Richter formula results from interpretation of empirical data as a straight line in the area of stretched exponent with small α . The solution of generalized logistic equation is applied for modeling distribution of foreshocks and aftershocks in the regions of Napa Valley 2014, and Sumatra 2004 earthquakes fitting the observed data well qualitatively and quantitatively. The relation of the theory developed with the Tsallis non-extensive statistics is demonstrated by solving the 3 generalized logistic equation for q >0. It is shown that the case 0< q <1 models “super-additive”, and the case q >1 models “sub-additive” structures and processes. In the case q =1 the probability BG p of distribution of a number of elements in a structure corresponds to Boltzmann-Gibbs distribution. For 0< q <1 and 0< N <1 the ratio q N of probability of distribution of a number of elements to BG p is greater than 1, and the rate of change of number of elements with size is greater than the rate of “growth” of “Boltzmann-Gibbs structure”, see for illustration graph for q in the Figure 2. Also, for q>1 the ratio q N is less than 1 which results in a slower rate of change of number of elements with size, compared to that rate of in a “Boltzmann-Gibbs structure”, graph for q =2 in the Figure 2. Data analysis
References