versión impresa ISSN 0016-7169
Geofís. Intl v.48 n.4 México oct./dic. 2009
Imaging the earthquake cycle with the MomentRatio Method: An exploration
C. Zhang1,2*, Y. Shi2, L. Ma3, J. Huang3 and C. Lomnitz4
1 China Earthquake Network Center, Beijing, 100045, China.
2 Laboratory of Computational Geodynamics, Graduate University of the Chinese Academy of Sciences, Beijing, 100049, China.
3 Institute of Earthquake Science, China Earthquake Administration, Beijing, 100036, China.
4 Instituto de Geofísica, Universidad Nacional Autónoma de México, Del. Coyoacán, 04510 Mexico City, Mexico. *Corresponding author: firstname.lastname@example.org.
Received: May 23, 2009
Accepted: July 29, 2009
El Método de Momentos Residuales (MRM) puede servir para monitorear en forma gráfica el ciclo del momento sísmico en un borde de placa. Es una técnica exploratoria que permite observar las deficiencias del valor acumulado de los momentos sísmicos para evaluar el ciclo de degradación y recuperación de la energía libre en una región sísmica activa. Se adopta una vida media de 20 años en el proceso de cicatrización de las rupturas sísmicas en el borde de la Placa del Pacífico. Se estima el momento residual para las regiones de Japón y México a intervalos de veinte años con el Catálogo de Pacheco y Sykes de sismos someros de magnitud M>7 en el siglo 20. Los picos de MRM tienden a mostrar alguna correlación con sismos de magnitud mayor que 7 que ocurrieron en el intervalo siguiente de 20 años.
Palabras clave: Momento sísmico, sismicidad, zonas de subducción, predicción de sismos.
The seismic moment cycle at plate boundaries may be imaged by the MomentRatio Method (MRM). In this exploratory technique the fluctuation of cumulative moment deficiency is used as an indicator of the cycle of degradation and recovery of free energy along an active tectonic region. A 20year halflife is assumed for the healing process of earthquake ruptures along the Pacific plate boundary. Using the PachecoSykes Catalog of large shallow earthquakes (M>7), we estimate the moment ratio for Japan and Mexico at 20year intervals during the 20th century. Momentratio peaks for 1980 correlate with subsequent earthquakes of magnitude M>7 during the following 20year period in Mexico and Japan.
Key words: Cumulative seismic moment, subduction zone seismicity, earthquake prediction.
Earthquakes share with avalanches the property that the strain configuration remains stationary until an event is triggered, after which the stress configuration changes abruptly. Each macroscopic strain event involves a release of strain energy. Postseismic phenomena include aftershocks and other transients of energy and matter, after which the strain configuration freezes again until a new event is produced. This process appears to be a common feature of a class of extended nonlinear systems (Fig. 1).
Earth is a steadystate, nonequilibrium, complex system that performs work by degrading sources of free energy, thereby increasing entropy (Kleidon and Lorenz, 2005). When fixed boundary conditions are assumed the time evolution of regional states depends on the initial conditions of the system (Hasselmann, 2002). A theory is worthless if it cannot predict (Lomnitz and Zhang, 2009). However, earthquake prediction is notoriously uncertain, as no known method can reliably predict when, where, and how large future earthquakes will occur (National Research Council, 2003, p. 8).
The momentratio method is useful for describing and imaging the earthquake process at plate boundaries (Lomnitz, 1985; 1993; 1994). The moment ratio parameter MRM is defined as
where MTM is the average tectonic moment in a segment m of the plate boundary and CSM is the cumulative seismic moment release at occurrence of the nth event in the catalog. The value of MTM is assumed to be stationary in time. Thus the moment ratio MRM fluctuates cyclically about a mean value in time. When the cumulative seismic moment CSM goes through a minimum a peak in MRM arises. In this paper we examine the correlation of MRM peaks with the occurrence of large earthquakes.
Let large earthquakes be represented as a point process in timespace, i.e., a method of randomly allocating events to hyperrectangles in a fourdimensional Euclidean space. A stationary Poisson process is a stochastic process such that
is the probability of finding no events in some interval of length x. Eq. (2) means that the interval between events has an exponential distribution. It may be shown (Daley and VereJones, 2002) that the numbers N of events in disjoint intervals must be independent random variables. Independence may be approximated if the events are rare, i.e., far between in space and in time. Large earthquakes are rare events: therefore they tend to be distributed as a stationary Poisson process.
In the real Earth, however, the numbers of earthquakes per interval are not independent random variables. The point process of large earthquakes may be shown to tend asymptotically to a stationary Poisson process under certain conditions (cf. Daley and VereJones, 2002). Suppose the strain release may be modeled as a Hawkes process or an ordinary ETAS model (VereJones and Ozaki, 1982; Ogata, 1988). Each plate boundary may be divided into a finite number of tectonic regions (see Lomnitz, 1974; Maps 1 to 3), much as a hydrological map may be divided into finite watersheds. If the mainshocks form a compound Poisson process with constant rate and fixed magnitude distribution, a thinning procedure in terms of magnitudes will tend to eliminate the aftershocks leaving only mainshocks. Thus the remaining events tend to form a Poisson process.
The momentratio method
Consider now the seismic gap model in its various versions (Fedotov, 1968; McCann et al., 1979; Nishenko and McCann, 1981; Shimazaki and Nakata, 1980; Schwartz and Coppersmith, 1984; Main, 1999; Mignan et al., 2006). This model is based on the idea that plate boundaries may be represented as a finite number of adjacent fixed segments of limited spatial extent. Actually the locations of seismic moment minima shift in space as well as in time (Lomnitz, 1993, 1994; Jackson and Kagan, 2006). Thus seismic ruptures overlap partially in timespace rather than breaking repeatedly the same segments.
This finding may be represented as a diagram of healing vs time (Fig. 2). It reflects the discovery that the seismic moment recovery rate, or "healing rate," closely matches the rate of moment release in large earthquakes. There is moment conservation along plate boundaries, as independently confirmed by Kagan (2002). This discovery leads to the formulation of the master equation
which is the basis of the MomentRatio Method (MRM) (Lomnitz. 1993). Here L(t) is the length of a seismic gap in the Central Chile plate boundary which is healing after a rupture of initial length L0 has occurred at time t =0, and is the halflife of the gap in years.
Notice in Fig. 2 that the halflife for all seismic gaps at the plate boundary is a constant on the order of =20 years. The numerical constant ln 2=0.69315 is the halflife factor. The location of successive gaps usually shifts in time but moment is conserved. The imaging procedure called Moment Release Method (MRM) attempts to reconstruct a macroscopic model of the earthquake process from the sequence of historical earthquakes. The observations from sketches such as Fig. 2 suggest some rules of construction: (a) two sequential ruptures cannot break the same segment, (b) every rupture heals exponentially in time from the edges inwards, (c) ruptures tend to nucleate at the edge of a previously healed area and propagate outward (Lomnitz, 1993, 1994). For example, in Fig. 2 we may observe that the large scar of the 1906 Valparaiso earthquake, though skewed, is not ruptured by the 1928 Talca earthquake or by the 1985 Valparaiso earthquake. Only scars less than 30 km wide are ruptured.
Thus ruptures take place sequentially along a plate boundary in such a manner that only segments that have previously healed according to Eq. (3) can rupture again.
A "scar" is a segment that has not yet healed. Scars can be outflanked by ruptures provided that the area of the scar has shrunk to less than the width of the plate interface.
Now consider the regional seismic process in time. The cumulative resident moment state CSM in a segment of finite length at time t may be written approximately as
where m0 is the seismic moment of an earthquake within the segment at time t=t1. Let i=1 denote the sequential number of the earliest earthquake in the catalog and assume that smaller earthquakes (M<7) can be neglected. Then
is the updated resident moment state in Equation (1), where i=n refers to the nth earthquake in the catalog. In this paper we use the catalog of world earthquakes of magnitude M>7 by Pacheco and Sykes (1992) which covers the period of 19001999. For purposes of testing the model we use the Harvard Catalogue of magnitudes for M>6.
Let us normalize the cumulative seismic moment over the mean moment for the Pacific plate boundary.
We assume that the mean moment state MTM remains constant in time but differs from region to region. We divide the plate boundary into segments of 100 km length and we adopt the following averaging algorithm:
where m designates the mth segment within the discretized onedimensional map of the plate boundary. Eq. (6) represents a smoothed space average of the updated cumulative seismic moment over a moving window of constant width comprising three adjacent segments of 100 km each.
The residual moment in a given year for any segment of the boundary may now be evaluated as referred to the year 1900 and the mean tectonic moment may be computed over a range of 150 km from any specific epicenter. Finally, Eq. (1) provides the expression for the normalized moment ratio.
A test of MRM
In order to evaluate the performance of the MRM approach we compute values of MRM along the ring of epicenters and plate boundaries surrounding the Pacific Ocean (Fig. 3). The Pacific Plate boundary has a total length of about 40,000 km. We discretize the boundary in units of 50 km and we filter the PachecoSykes Catalog to exclude events deeper than 80 km or intraplate events not belonging to the seismogenic zone of the Pacific Rim. Only events of magnitude 7 and greater are included.
Every earthquake in this sample was processed by distributing its seismic moment over the adjacent segments as follows: 0<Δ<50km, 22%; 50<Δ<100km, 12%; 100<Δ<<150km, 6%, where Δ is the epicentral distance in kilometers from the earthquake.
Computation of MRM values was performed by using MTM values from all preceding earthquakes in the same area and introducing them in Eq (1). We assume zero cumulative moment at t=0 as the initial condition. The sample size increases with time; thus the estimation of residual moments improved significantly over time in this test.
The results of our calculations are summarized in Figs. 4, 5 and 6, representing MRM plots for Japan, Mexico and the Pacific Ring. Fig. 4 shows the MRM plot in Japan, an area of high rate of seismic energy release. Fig. 5 (Mexico) is typical of an area of intermediate energy release, and Fig. 6 summarizes the MRM results for the entire Pacific Ring.
(a) The Japan region.
Fig. 4 shows saddle points in cumulative seismic moment that tend to appear as MRM peaks. Neighboring regions are activated or deactivated, and saddle points shift in time and space. These fluctuations affect the values of MRM. Every large earthquake depresses the MRM value, causing new peaks to arise in the area. Thus the effect of a large moment release is not confined to the immediate neighborhood but extends beyond the seismic gap.
Fig. 4 shows nine major peaks in MRM for the time datum 1980 (red line). These peaks appear to be associated with earlier MRM peaks as they appeared twenty years earlier (light blue line). Thus MRM peaks tend to persist over time periods exceeding 10 years. At least five of these peaks (south Kyushu, north Kyushu, Bungo Strait, Hakodate, and Hokkaido) are associated with maxima of earthquake activity during the following twentyyear period (19802000, dark blue/green). Another peak appears to be associated with the Kobe earthquake of 1985 though this earthquake does not appear in the PachecoSykes Catalogue as its magnitude was slightly below 7. Note that both trust and strikeslip events are included.
Thus six out of nine MRM maxima appear to be associated with relevant earthquake activity within a twentyyear period and within the 300 km location uncertainty caused by discretization. The more relevant peaks of MRM in terms of predictive power seem to appear decades before a large event. In the case of the Kobe earthquake the MRM peak appeared in 1960 (MRM=7), that is, thirtyfive years before the earthquake. If such observations may be confirmed in experiments to be performed in the future this might mean that significant earthquakes can build up moment deficiencies for decades.
Consider the Japanese results in more detail. For the 2000 MRM the three major MRM peaks in the Japanese area (from N to S) are Hokkaido, the Nankai Trough area (about 300 km southwest of Tokyo), and the Bungo Strait area between Honshu and Kyushu. All three peaks are large (between MRM=7 and MRM= 10). There were no major events in 19802000 in these areas. On the other hand, both Hokkaido and Bungo were areas of high average activity in 19802000, and the Nankai Trough area remains an area of high seismic risk.
(b) The Mexico area.
The Kobe earthquake occurred unexpectedly. This was also the case of the 1985 Michoacan, Mexico earthquake (MRM=5, Fig. 5). As in the case of the 1995 Kobe earthquake it was a surprise to seismologists as the "Michoacan gap" was believed to be inactive and aseismic (LeFevre, et al., 1985). Both events caused the worst disasters in the history of their respective regions.
In the Mexican region (Fig. 5), there were four major MRM peaks including (from N to S) Michoacan, Acapulco, Chiapas, and northern Guatemala. The Acapulco MRM peak is in the "Guerrero Gap" which has been singled out as site for a prospective major earthquake. The Michoacan peak almost reached MRM= 10 in 1960 but it decayed to MRM=6 after the 1985 earthquake. It has generated several strong earthquakes since the earthquake and might still be hazardous. The Chiapas peak is below MRM=6 and the MRM situation in Oaxaca seems to fluctuate around similar background values.
From the MRM diagram, it seems likely that the 1985 Michoacan earthquake (M8.1) did not completely remove the hazard of large earthquakes in the area. The Guerrero gap, like the Nankai Trough in Japan, has long been regarded as hazardous. In 2000 it is beginning to look threatening, with a peak value of MRM=11. Central Oaxaca might also be due for a mediumsized shock.
In MRM computations the state of the system is predicted at all times from the initial conditions and the boundary conditions. As the actual system is stochastic any attempt to compare the input vs the output could produce spurious high correlations for any algorithm one might wish to consider. In other words, MRM is a highly modeldependent procedure. The output of the MRM procedure reveals unexpected peaks as it takes into account the variations of seismicity in neighboring segments.
On the other hand, MRM is an initial attempt at modeling moment deficiencies in terms of the influence of seismicity on neighboring regions. Here we explicitly assume that the subduction zone is uniform along strike. In the future we may attempt to introduce inhomogeneous stress transfer in the strike and downdip directions. The input time series is a catalog of large earthquakes (M>7) which is a highly nonlinear sample of the Earth's seismicity.
Extrapolating such shorttime samples, on the order of one to two decades, to seismicity estimates is a hazardous proposition. The observed coincidences between MRM and relevant earthquakes may be promising but should be considered highly tentative and preliminary.
The MRM image of the Pacific Ring (Fig. 6) is significant as it reflects the nearchaotic nature of the data. Major peaks with amplitudes of MRM=12 or more do occurbut they often correspond to regions of known high activity such as the Kurile Islands, where large events are common. As the data are scarce the uncertainties are high. Introducing probabilities explicitly at this stage is not thought to be useful, as the results might be interpreted as predictions (Jackson and Kagan, 2006). Plate boundaries will rupture sooner or later. The information provided by MRM is intended as a tentative model to alert seismologists to possible scenarios of what might happen, and not as a forecast of any future event.
These caveats must be kept in mind. The present exercise does seem to indicate a likely direction of future research. As someone famously said, prediction is uncertainespecially when it concerns the future.
(a) The halflife r of a seismic gap.
In this work we arbitrarily use a mean value =20 years for the Pacific Plate boundary, but it seems likely that the halflife could vary with the seismicity.
(b) The spatial window of seismicmoment diffusion.
We assume that the seismic moment m0 of any earthquake affects a mean rupture length of 300 km and that the moment release is diffused according to a fixed ratio over the central and peripheral portions of the rupture. Actually the rupture length depends on the size of the earthquake. Thus the diffusion pattern of the seismic moment should be made to depend on the magnitude of the event.
We assume that the moments of successive earthquakes are additive by linear superposition, and that the residual moment "heals" exponentially with time. Healing of seismic ruptures may be regarded as an empirical fact, but it might be a result of strength decay of the plate boundary as a result of corrosion (Lomnitz and Zhang, 2009). When a segment of plate boundary ruptures in an earthquake it features a period of heightened seismic activity during the aftershock process but eventually the area of heightened activity evolves into a seismic gap. The location and extent of any gap depends on the evolution of adjacent gaps.
In this work we propose to normalize the residual seismic moment in order to obtain the moment ratio. The underlying assumption is that the probability of a rupture depends, among other factors, on the mean tectonic moment which varies strongly along the plate boundary. Averaging the residual moments over the width of the window of moment diffusion may represent an excessive simplification of the physics of the process.
(e) Discretization scheme.
For computational purposes we discretize the plate boundary as a onedimensional string or chain which extends clockwise from the Macquarie Islands, south of New Zealand, to the tip of Patagonia. The arbitrary unit of discretization is 100 kilometers. Actually, of course, a plate boundary is a threedimensional object.
(f) Stochastic modeling.
The present approach is halfway deterministic. Probabilistic forecasts should involve a system of coupled differential equations in the state variables, known as prognostic equations (Hasselmann, 2002; VereJones, 1994). Numerical integration of such a system is a challenging proposition in many respects, including computer time. A stochastic equation is of the general form
where F is the mean coupling force between the cumulative residual moment and a slowly varying MRM, and is a stochastic forcing function generated by the shortterm variability of seismic weather.
In conclusion, the Earth's lithosphere is a steadystate complex system that is subject to instabilities but which may admit some restricted forms of forecasting.
Procedures of imaging the state of the lithosphere, including MRM, are worth exploring.
Dr. Zhou Yuanze offered valuable ideas and advice. One of us (C.Z.) was supported by a scholarship grant of the Mexican State Department (SRE) and UNAM. We gratefully acknowledge grant 2004CB418405 (973 project) of the National Basic Research Program of China and grant 90814014 of the National Natural Science Foundation of China.
Daley, D. J. and D. VereJones, 2002. An introduction to the Theory of point processes, Vol 1, 2nd. Ed., Springer, New York, 469 pp. [ Links ]
Fedotov, S. A., 1968. On the seismic cycle, feasibility of quantitative seismic zoning and longterm seismic prediction, in Seismic Zoning of the USSR, Nauka, Moscow, 121/150 (in Russian.) [ Links ]
Harris, R. A. and J. R. Arrowsmith, 2006. Introduction to the Special Issue on the 2004 Parkfield earthquake and the Parkfield Prediction Experiment, Bull. Seismol. Soc.Am., 96, S1S10. [ Links ]
Hasselmann, K., 2002. Is climate predictable? In The Science of Disasters, Bunde, A., J. Kropp and H. J. Schellnhuber, eds. (Springer, Heidelberg, pp. 141169). [ Links ]
Jackson, D. D. and Y. Y. Kagan, 2006. The 2004 Parkfield earthquake, the 1985 prediction, and characteristic earthquakes: Lessons for the future, Bull. Seismol. Soc. Am., 96, S397S409. [ Links ]
Kagan, Y. Y., 2002. Seismic moment distribution revisited: II. Moment conservation principle, Geophys. J. Int., 149, 731754. [ Links ]
Kagan,Y. Y. and D. D. Jackson, 1992. Seismic gap hypothesis: ten years after (J). J. Geophys. Res., 96, 21, 41921,431. [ Links ]
Kanamori, H., 1981. The nature of seismicity patterns before large earthquakes, in Earthquake Prediction, Maurice Ewing Series, 4, AM. Geophys. Union, Washington, DC., 119. [ Links ]
LeFevre, L. V. and K. C. McNally, 1985. Stress distribution and subduction of aseismic ridges in the Middle America subduction zone. J. Geophys. Res., 90, B6:44954510; DOI: 10.1029/JB090iB06p04495. [ Links ]
Lomnitz, C., 1993. Momentratio imaging of seismic regions for earthquake prediction. Geophys. Res. Letters, 20, 21712174. [ Links ]
Lomnitz, C., 1994. Fundamentals of Earthquake Prediction (John Wiley & Sons, New York). [ Links ]
Lomnitz, C., 1985. Tectonic feedback and the earthquake cycle. Pageoph, 123, 667682. [ Links ]
Lomnitz, C. and C.J. Zhang, 2009. Parkfield revisited: I. Data retrieval, Lithosphere, 1, 227234, doi: 10.1130/ L14.1. [ Links ]
Lorenz, E. N., 1975. Climate predictability: The physical basis of climate and climate modeling. World Meteorological Organization, Rept. 16, p. 132. [ Links ]
Main, I. (Moderator), 1999. Is the reliable prediction of individual earthquakes a realistic scientific goal? Debate in Nature, www.nature.com/nature/debates/earthquake_frameset.html. [ Links ]
McCann, W. R., S. P. Nishenko, L. R. Sykes and J. Krause, 1979. Seismic gaps and plate tectonics: Seismic potential for major boundaries, Pure Appl. Geophys. 117, 10821147. [ Links ]
Mignan, A., G. King, D. Bowman, R. Lacassin and R. Dmowska, 2006. Seismic activity in the SumatraJava region prior to the December 26, 2004 (Mw=9.09.3) and March 28, 2005 (Mw=8.7) earthquakes, Earth Planet. Sci. Lett., 244, 639654. [ Links ]
Ogata, Y., 1988. Statistical models for earthquake occurrences and residual analysis for point processes. J. Amer. Statist. Assoc., 83, 927. [ Links ]
Pacheco, J. and L. R. Sykes, 1992. Seismic moment catalog of large shallow earthquakes, 1900 to 1989, Bulletin of the Seismological Society of America, 82, 13061349. [ Links ]
PeraltaFabi, R., E. MoralesGamboa, V. RomeroRochin, and J. LomnitzAdler, 1994. Development of a theory for 2D avalanches, in Lectures on Thermodynamics and Statistical Mechanics, M. Lopez de Haro and C. Varea, Eds., World Scientific, Singapore, 241253. [ Links ]
Schwartz, D. P. and K. J. Coppersmith, 1984. Fault behavior and characteristic earthquakes: Examples from Wasatch and San Andreas fault zones, J. Geophys. Res., 89, 56815698. [ Links ]
Shimazaki, K. and T. Nakata, 1980. Timepredictable recurrence model for large earthquakes, Geophys. Res. Lett., 7, 279282. [ Links ]
VereJones, D. and T. Ozaki, 1982. Some examples of statistical inference applied to earthquake data. Ann. Inst. Statist. Math. 34, 189207. [ Links ]
VereJones, D., 1994. Statistical models for earthquake occurrence: clusters, cycles and characteristic earthquakes, in Proc. First USJapan Conf. Frontiers of Statistical Modeling, H. Bozdogan, ed., Kluwer Academic, Dordrecht, 105136. [ Links ]