Microsoft Word - articulo Alexander Caneva corregido331marzo.doc EARTH SCIENCES RESEARCH JOURNAL Earth Sci. Res. J. Vol. 8, No. 1 (Dec. 2004) : 3 - 9 Manuscript received May 2004. Paper accepted October 2004. 3 USING THE FRACTAL DIMENSION OF EARTHQUAKE DISTRIBUTIONS AND THE SLOPE OF THE RECURRENCE CURVE TO FORECAST EARTHQUAKES IN COLOMBIA Alexander Caneva 1 and Vladimir Smirnov 2 E-mail address: alexander.caneva@uan.edu.co 1 Research Center, University Antonio Nariño, Bogotá, Colombia – Geophysics – 2 Earth Physics Division, Faculty of Physics, Lomonosov Moscow State University, Moscow, Russia ABSTRACT An estimation and analysis of temporal and spatial variations (on the surface, and depending on depth) of the parameters of the seismic regime (fractal dimension of the earthquakes distribution d and slope of the seismic recurrence curve b) in Colombia are carried out, considering the fractal properties of seismicity. The variations of the difference (b – d/a), in time are analyzed (where a is the exponent of the power law E ~ la which establishes the relation between the energy E of the earthquake and the size l of its focus). This difference describes the deviation of the geophysical medium from a “stable” state for intervals of time in which strong earthquakes occur. The possibilities to use the variations of the parameters of seismic regime to forecast earthquakes are discussed. Key Words: fractal dimension, earthquake prediction, seismic recurrence, seismic regime. RESUMEN Se lleva a cabo la estimación y el análisis de las variaciones en el tiempo y en el espacio (en superficie y en función de la profundidad) de los parámetros del régimen sísmico (dimensión fractal de la distribución de sismos d y la pendiente del gráfico de recurrencia sísmica b) de Colombia considerando las propiedades fractales de la sismicidad. Se analizan las variaciones de la diferencia (b – d/a) en tiempo (donde a es el exponente de la expresión E ~ la , que vincula la energía E del sismo con el tamaño l de su foco). Esta diferencia describe una desviación que experimenta el medio geofísico con respecto a un estado estable en los intervalos de tiempo en que ocurren terremotos fuertes. Se analizan las posibilidades de utilizar las variaciones de los parámetros de régimen sísmico para pronosticar sismos fuertes. Palabras clave: dimensión fractal, predicción de terremotos, recurrencia sísmica, régimen sísmico.  2004 ESRJ - Unibiblos. Alexander Caneva and Vladimir Smirnov 4 INTRODUCTION The term seismic regime is understood as an assemblage of earthquakes distributed in space and time. The quantitative characteristics of the regime are the statistical parameters of this assemblage (Smirnov, 1995). It is known that earthquakes are distributed according to their energies, in space, and time in an irregular way (Pisarenko, 1989). The distribution of earthquakes according to their energies is given by the Gutenberg-Richter law, which represents the distribution of the number N of earthquakes with respect to their energies E as a power law: N ~ Eb (Aki & Richards, 1980). According to this law the number of small events in a given region of space is greater than the number of strong events. On the other hand, the distribution of defects N, according to their sizes l in a certain region of space also follows a power law: N ~ l-d (Aki, 1981; Fukao and Furumoto, 1985). The heterogeneity of the assemblage of earthquakes corresponds to their self-similar, hierarchic (fractal) structure (Keilis-Borok et al. 1989). The fractal character of the spatial and temporal heterogeneity of seismicity corresponds to some general properties of the seismic regime (Rykunov et al. 1987), and these fractal properties together with the law of recurrence are the basis for understanding the self-similarity of the seismic process (Sadovsky & Pisarenko, 1991). Considering that a fractal object is heterogeneous on all scales (Crownover, 1995), it is indispensable to reformulate the seismic recurrence law paying attention to the fractal properties of earthquakes assemblage. In the same way, we have to connect the spatial, temporal and energetic statistical properties of the seismicity one to each other. These properties are the expression of a physical process: the evolution of the lithosphere (Smirnov, 1995). The quantitative estimation of the parameters defining the seismic regime allows clarifying their features. Let us analyze these aspects in more detail. THE SEISMIC REGIME PARAMETERS Slope of the curve of recurrence b The characteristic of self-similarity of seismic process in the energetic aspect is the slope of graphic of recurrence. This graphic shows the relation, existing between the numbers N of earthquakes in a region having different values of energy E: ( )1 2 1 2 1 +−       = b E E N N (1) where Ni is the number of events with energy between (Ei, Ei+dE), b is the slope of the graphic of recurrence. From equation (1): ( ) AKKbN loglog 0 +−−= , (2) where EK log= is called the earthquakes energetic class (Pisarenko, 1989), and N is the number of events in the range K+dK. In the classic formulation (Aki & Richards, 1980) the law of seismic recurrence expresses the relation between the number of earthquakes N and the energy of those earthquakes E. This law characterizes the probability of occurrence of earthquakes with certain energy. In this case it is assumed that the earthquakes are distributed uniformly in space and time, so the classic law of seismic recurrence reflects just the energetic properties of seismicity without considering the earthquakes distribution in space and time. Fractal dimension d of the earthquakes spatial distribution The fractal dimension d can be understood as a measure of geometric self-similarity. It determinates the relationship between the numbers of non-empty cells Ni (cells which contain at least one element of the assemblage) of different size li inside the region of space (Fukao & Furumoto, 1985): d l l N N −       = 2 1 2 1 . (3) The number of non-empty cells of size l satisfies: ( ) dllN −= (4) whereas the density of events in a Euclidean space of dimension r is defined as: rl n=µ , (5) being n the average of events in cells of size l. If there are altogether m events distributed among N cells in a certain region of space then: N m n = (6) and from equation (4) rdl −≈µ . When the events are distributed uniformly rd = and µ = constant. If the assemblage of events has a fractal structure then d < r and al −≈µ , where 0>−=α dr , indicating that with the reduction of the cell’s size the density of events increases. So, the density estimation over events according to equation (5) loses its sense for a Using the fractal dimension of earthquake distribution and the slope of the recurrence curve to forecast earthquakes in Colombia fractal assemblage since its value depends on the cell size l. GENERALIZED LAW OF SEISMIC RECURRENCE A previous work (Smirnov, 1995) tried to reconcile the classical definition about seismic regime parameters (in terms of b and d – values) with the fractal features of seismicity. In an attempt to formalize certain statistical characteristics of seismic recurrence, (Smirnov, 1995) Smirnov introduced the so-called generalized law of seismic recurrence. The analytical expressions for the seismic regime parameters were derived from the following aspects: After an earthquake had occurred in a certain region of space, there is an interval time and a certain rank of energy values such that no further earthquake with the same energy can occurs during this period. Assume that the size of the region is given by: αλ= 0lR , (7) where l0 is a linear size of the focus of the earthquake (characteristic size of the destroyed region of the physical medium defining the scale of the event) (Rice, 1982); λ and α are constants. Physically, R is the size of the region that had distended as a result of the earthquake. Assume that the time of prohibition of new earthquakes in region R is: βθ=τ 0l , (8) where θ and β are constants. The seismogenic region (the assemblage of seismic focuses, for epicenters as for hypocenters) has a fractal spatial structure (Aki, 1981). If we divide a region of size L into smaller regions (for example spheres) of ∆ size, the number of such sub-regions that contain earthquakes becomes: ( ) d L nn       ∆ =∆= , (9) where the constant d is the fractal dimension of the seismogenic region (Crownover, 1995). Note that the sismogenic region can be related to a system of faults, which also has fractal properties (Ulomov, 1993). The energy of the elastic earthquake waves depends on the size of the focus according to a power law: alE 0ε= . (10) This formula is based on an empirical relation between the energy of the earthquake E and the linear size of focus l0 (Kasahara, 1981). From these assumption it follows: In a region of space of size L, during a time τ, a number N of earthquakes can occur, each one with focus of size l0, such that the regions where these stresses are released as a result of earthquakes, do not intersect each other. The number of these non-intersecting regions, or the number of distinct earthquakes N, is defined by the fractal geometry of the seismicity of this region. From the anterior postulates we may formulate a new law of seismic recurrence, which takes into consideration the fractal properties of seismicity. Covering the studied region with spheres of size ( )0lR=∆ (that do not intersect each other), and taking the equation (9) ( )dR LN = . If during the time interval τ, on the average, N earthquakes with linear dimensions l0 take place, then, during a time interval T, on the average, T/τ times more events will occur. Therefore in the region of size L, during a time interval T, ( )dRL TN τ = earthquakes will be observed with linear l0 focal size. From equation (7) and equation (8), and considering that ( )aEl 1 0 ε= , the earthquakes number becomes ( ) ( ) a d a dd E TL N β+α−β+α ε      θ       λ = . From this expression we obtain the generalized law of seismic recurrence that considers the fractal properties of the seismicity: BTLdbKN +++−= logloglog (11) Where EK log= , ( ) a db β+α= and θ−λ−ε= logloglog dbB From equation (11), analytical expressions are obtained for the parameters of the seismic regime: TL E N b , lg lg       ∂ ∂−= (12) and TE L N d , lg lg       ∂ ∂= (13) where the subscripts L (size of the region), T (the interval of time) and E (the energy of the earthquakes) indicate the constancy of parameters. A similar to the equation (11) expression was analyzed and proposed (Chelidze, 1993) based on results presented in papers (Mazhkenov, 1989; Smalley et al. 1987). The fractal seismicity properties, besides the power law of seismic recurrence, assure the self-similarity of seismic processes. From the generalized law of seismic recurrence: Alexander Caneva and Vladimir Smirnov 6 N ∼ db LE ⋅− ∼ dab ll 00 ⋅− ∼ ( )abdl − 0 . (14) The worldwide average values of the parameters of the generalized law of seismic recurrence are, respectively: 5.0≈b and 5.1≈d (Aki, 1981; Hirata, 1989), whereas the value commonly accepted of a is near 3 (Kasahara, 1981), therefore the value of the exponent in formula (14) is 0≈− abd , then we have: abd ≈ . (15) This empirical expression (first proposed by Keilis-Borok et al. 1989) should not be understood like a correlation between the b and d values. This expression describes the relation between the average values (taking over sufficiently long intervals of time) of these parameters. Several references where the seismic activity after a main shock is analyzed (e.g. Hirata, 1989) point out immediately after the main shock the deviation with respect to the theoretical relation d ≈ ab reaches the maximum value. This means, d ≠ ab when the geophysical medium is far from its stationary (background, common, normal) level of seismicity. With the attenuation of aftershocks, the difference ( )abd − returns again to its average value near zero. There is also published information about variations of the same type in case of induced seismicity in dams (Smirnov et al. 1994). In this case, the distance of the magnitude |d - ab| is maximal during the period of dam filling. It seems that the relation d ≈ ab depends on a certain “stable” state of the geophysical medium. As a result of a strong earthquake on a smaller scale in comparison with the characteristic scale of the global rupture, the “stable” dynamic state loses its original balance, which leads to a redistribution of tensions in that region forming, after a certain time, a new “stable” scheme. This way, it seems that the deviation from the relation d ≈ ab is a measure of the (seismic) stability of the geophysical medium. ESTIMATING THE SEISMIC REGIME PARAMETERS Estimating the b – value The maximum likelihood estimation (Hudson, 1964) of the parameters in the law of seismic recurrence was described in (Sadovsky & Pisarenko, 1991) expressing this law mathematically as follows: AKbN jj loglog +∆−= , Where j = r, r+1,...,r+n, Nj, is the average of the number of earthquakes in the interval ( ) ( )       +∆+−∆+ 2 1 02 1 0 10,10 jKKjKK , normalized at the interval T in study; ( )EK log= , which corresponds to the parameter of seismic activity A; ∆K is the variation of the energetic class for the construction of the histogram. The maximum likelihood estimates of b and A is found from the system of equations: ( ) ( )bTAMK bTAM ϕ′+∆ =ϕ− 1 0 10lg 0 (16) where: ( ) ( ) ∑ ∑ ∑ ∑ + = ∆− = ∆− = + = ∆−ϕ′ =ϕ = = nr rj kb n rj kb n rj j nr rj j jKb b jmM mM j 1010log 10 1 0 (17) and {mj} is the histogram of distribution of earthquakes with respect to energy (number of earthquakes in the j - cell). This way, the empirical data enter into the equations (12) in form of the parameters: M0 and M1. Estimating the d – value To estimate the fractal dimension of an assemblage of earthquakes we used the correlation dimension: ( ) l lC d l log log lim 0→ −= , with l the linear size of the cells into which the Euclidian space is divided. C(l) is the integral of correlation: ( ) ( )lrrNlC ji ≤−= , where N is the number of pairs of events separated by a distance not greater than l (Feder, 1988). For empirical data, we cannot take the limit for l → 0. To solve this problem we use a scaling approach, which consists of selecting a rank of l values where the relation between ( )Clog and ( )llog is close to linear. The region of scaling is limited below by the accuracy of the data, and above by the size of the studied region or by real variations of the assemblage structure (nonlinearities of the curve ( ) ( )( )lfC loglog = . Such nonlinearities indicate nontrivial changes with size of the fractal dimension. Estimating fractal dimension of the earthquakes assemblage hypocenters was carried out by constructing the correlation integral histogram; then - Using the fractal dimension of earthquake distribution and the slope of the recurrence curve to forecast earthquakes in Colombia based on it- the function: ( ) ( )( )lfC loglog = was constructed; the region of scaling was identified; and the value of d was calculated by regression. Table 1. Colombian seismic regime parameters values Region b σb d σd Colombia 0.59 0.02 1.55 0.03 Colombia (shallow) 0.56 0.02 1.53 0.03 Colombia (deep) 0.57 0.03 1.61 0.03 RESULTS The catalogue of earthquakes in Colombia was provided by the: Institute of Geoscientific, Mining- Environmental and Nuclear Investigation and Information of Colombia, Ingeominas (Responsible for the National Seismological Network of Colombia) (Ingeominas, 1999). The catalogue contained 17221 records since 1993 to 1999. After eliminating double (repeated) records, records with erroneous format, and aftershocks the number of records diminished to 10898. The estimation of the parameters of seismic regime was carried out from this “clean” catalogue. Estimating the b and d– values for Colombia The estimation of b and d values was carried out using specially designed software. Table 1 shows the results for the entire Colombian territory, including a differentiation for shallow and deep seismicity. Figure 1. Variations of b-value � , d-value Ο and in dashed line (b-d/a) for earthquakes as function of time for the entire territory of Colombia. The symbol indicates strongest earthquakes. When these curves are analyzed it is difficult to identify some regularity except the counter-phase variation of parameters in the intervals 1996-1997 and 1998-1999. The b-value is related to the existing proportion between the number of weak and strong events An increase in b–value corresponds to a rise in the number of weak events and a decrease in the number of strong events, while a reduction of b–value indicates a rise in the number of strong events and a decrease in weak events. The d–value characterizes the level of spatial events homogeneity. The less the value of the fractal dimension, the more located are the events; that is a decreasing of d corresponds to a grouping of the events. Figure 2. Figure 1. Variations of b-value � , d-value Ο and in dashed line (b-d/a) for earthquakes as function of time for the entire territory of Colombia. The symbol indicates strongest earthquakes. Figure 1 depicts the variations of the curve b - d/a, which is equivalent to the expression (ba - d) that corresponds to a stable state of the medium. In intervals of time in which strong earthquakes appear, no peculiar |b - d/a| variations are observed, possibly due to the heterogeneity of events when considered over the whole Colombian territory (as the type of seismicity, whether of subduction, volcanic, types are different for distinct regions in the country). To reduce the degree of heterogeneity in depth we separately considered the events with H < 100 km (Figure 2). In this case the counter-phase variation of the b and d – values is observed for almost all the time intervals. But this behavior, as well as a greater deviation of |b - d/a| from zero, can be observed starting at 1996. After reaching the maximum at the end of 1997, the deviation tends again to zero. Exactly in this interval of time there occurred the strongest events with mb = 5.8 and H = 0 km in the region of Urabá in northwest of Colombia (this region is characterized by shallow seismicity). Alexander Caneva and Vladimir Smirnov 8 For deep seismicity no particular behavior was observed. Although in this case the coincidence of the presence of a strong earthquake and the deviation of the difference |b - d/a| from zero is evident. This unique example is not enough for considering this deviation as a prognostic element. It is important to note that there are several strong earthquakes at the starting time of the study and one more at the end (Figure 2). However, as for statistical analysis it is not appropriate to evaluate the initial and final parts of a catalogue, we ignored those events and considered as relevant only one event, which took place in the middle of the analyzed period. This way, the relevance for earthquake prognosis of the value of the |b - d/a| difference could not be checked in the case of Colombia. In a paper (Smirnov, 1995), the deviation from zero of the |b - d/a| difference is described, in relation to catastrophic earthquakes with Ms ~ 7 – 8, whereas the events considered as strong in this work are only Ms ~ 5 - 6. Figure 3 shows the variations of the b–value in space (on the surface). In the region near 75º Western longitude and 5º North latitude, an “island” of diminished values of b can be observed. Figure 3. Distribution of b-value on the surface for the entire territory of Colombia This is a region of low seismicity located in the valley of Magdalena River. When studying the b-value behavior for deep and for shallow seismicity, respectively, we must take into account that most earthquakes with H > 100 km take place in the seismically active zone called “nest of Bucaramanga” (near 73º Western longitude and 5º North latitude) while the seismicity in the rest of the country is mostly shallow. This way the distribution of the values of b on the surface for shallow seismicity does not change in relation to Figure 3. Figure 4 shows the variations of the average values of b and d as a function of the depth for the entire Colombian territory. It is interesting to note the increase of b and the decrease of d starting at around of 100 km depths, related to the seismicity in the “nest of Bucaramanga”. The decrease of d– value indicates an increase of the spatial earthquakes concentration at this depth, matching with the characteristics of a nest (region) of deep seismic activity. These counter-phase variations with the depth, starting at 100 km, are possibly related to specific conditions of pressure and temperature at those depths (particularly in the “nest of Bucaramanga” that, to a great extent, defines the deep seismicity in Colombia (Caneva, 2000)). Figure 4. Variations of b-value Ο and d-values � in function of depth for the entire territory of Colombia. Using the fractal dimension of earthquake distribution and the slope of the recurrence curve to forecast earthquakes in Colombia CONCLUSIONS Based on the catalogue of earthquakes of Colombia we estimated the parameters of its seismic regime. The average values of the recurrence slope graph b were obtained in general for Colombia (b ≈ 0.6) and separately for the shallow and for deep seismicity. The analysis of the b – values variations in time and in space (both on surface and with depth) has also been accomplished. The average fractal dimension value d of epicenters for Colombian earthquakes was found (d ≈ 1.6). Counter-phase variations of the b and d – values in time and with depth were observed. The variations of the (b – d/a) difference were analyzed. This difference corresponds to a deviation from a stable state of the geophysical medium, in such intervals of time in which strong earthquakes occur. The absence of significant correlation does not allow, yet, using the variations of the (b – d/a) difference to forecast earthquakes in Colombia. REFERENCES Aki, K. (1981). A Probabilistic Synthesis of Precursory Phenomena. Earthquake prediction. Amer. Geophys Union. Washington. p. 556 -574. Aki, K. and Richards, P. (1980). Quantitative Seismology. Theory and Methods. W. H. Freeman. Caneva, A. (2000) Fractal Properties of Seismicity of Colombia. Ph.D. Thesis. MSU, Moscow. In Russian. Chelidze, T. (1993). Fractal Damage Mechanics of Geomaterials. Terra Nova 5, pp. 421 – 437. Crownower R. (1995). Introduction to Fractals and Chaos. Jones and Barlett Publishers.. Feder, J. (1989). Fractals. Plenum Press, New York. Fukao,Y., Furumoto, M. (1985). Hierarchy in Earthquake Size Distribution. Phys. Earth Planet. Inter. V. 37, pp.149 – 168. Hirata,T. (1989). A Correlation Between the b-value and Fractal Dimension of Earthquakes. J. of geophys. res. 94. B6. pp. 7507 - 7514. Hudson, D. J. (1964). Statistics. Lectures on Elementary Statistics and Probability. Geneva. 2. Ingeominas. Red Sismológica Nacional de Colombia. Ingeominas.http://www.ingeomin.gov.co. Bogota, Colombia. Kasahara, K. (1981). Earthquake Mechanics. Cambridge University Press. 264 p. Keilis-Borok, V. I., Kosobokov V. G. & Mazhkenov V. G. (1989). Similarity in Spatial Distribution of Seismicity. Vychislitelnaya seismologia. V. 22. pp. 28 – 40. In Russian. Mazhkenov, S. A. (1989). Diagnostics of Times of Increased Probability of Earthquakes in Northern Tjan-Shan. Ph.D. Thesis. Moscow. In Russian. Pisarenko, V. F. (1989). The Law of Recurrence of Earthquakes. .diskretnye svoistva geofizicheskoi sredy (pod red. Sadovskogo M. A.). М.: Nauka. In Russian. Rice, J. (1982). Mechanics of Earthquake Focus. M. Mir. Rykunov, L. N., Smirnov V. B., Starovoyt Y. O. & Chubarova, O. S. (1987). Self-similarity of Seismic Emition on Time. Dokl. AN SSSR. V. 297. No 6. pp. 1337 – 1341. In Russian. Sadovskii, M. A. & Pisarenko V. F. (1991). Some Peculiarities of Seismic Process in Block Media. М.: Nauka. 96 p. In Russian. Smalley, R. F., Chatelain H. -L., Turcotte D. L. & Prevot R. (1987) A Fractal Approach to the Clustering of Earthquakes. Bull. Seismic. Soc. Am. 77. pp. 1368 – 1382. Smirnov, B. V. (1995). Recurrence of Earthquakes and Parameters of Seismic Regime. Vulkanologia and seismologia. No 3. pp. 59 – 70. In Russian. Smirnov, V. B., Cherepantsev A. S. & Mirzoev V. K. (1994). Fractal Properties of Induced Seismicity in the Region of the Dam Nurekskoe. Induced Seismicity. М.: Nauka. p. 138 – 147. In Russian Ulomov, V. I. (1993). Global Order of the Seismic- geo-dynamic Structures and Some Peculiarities of the Seismic Zoning and Forecasting of Earthquakes. Seismicity and Seismic Zoning of Northern Euro-Asia. Pod red. Ulomova V. I. М.: IFZ RAN.pp. 24 – 44. In Russian.