Open access peer-reviewed chapter

Aftershock Identification Through Genetic Fault-Plane Fitting

By F.A.Nava, V.H.Márquez and J.F.Granados

Submitted: November 29th 2010Reviewed: April 25th 2011Published: August 1st 2011

DOI: 10.5772/23733

Downloaded: 2458

1. Introduction

Since the earlier times of documented seismological observations, it was noticed that an earthquake (usually a large one) was followed by a sequence of many smaller earthquakes, originating in the epicentral region; the first, larger, earthquake is called the mainshock, or main shockor main event, and the following, smaller, earthquakes are called aftershocks. These sequences, and their spatial and temporal distributions, depend on the characteristics of the mainshock and on the physical properties and the state of stress, strain, temperature, etc., of the region of occurrence (Kisslinger, 1996).

The observation that in many aftershock sequences the magnitude of the largest aftershock is about ΔM1.2less than that of the mainshock is known as Båth’s law (Helmstetter & Sornette, 2003; Richter, 1958); it implies an energy ratio E(Maft)/E(Mmain)~0.007between largest aftershock and mainshock. An earthquake following some mainshock but not small enough to be considered an aftershock by the Båth’s law criterion is often considered to be a mainshock in its own right and to constitute, together with the mainshock, a multiple eventor multiplet. A large event followed by an even larger one is demoted from mainshock to foreshock; a probabilistic calculation, based on the observation that aftershocks follow the Gutenberg-Richter relation (Gutenberg & Richter, 1954; Kisslinger & Jones, 1991; Shcherbakov et al., 2005), indicates that for ΔM=1.2the seismicity following a given mainshock has a ~ 6.3% probability of including a “daughter” event larger than the mainshock (Helmstetter et al., 2003; Holliday et al, 2008).

Thus, the simple definition of aftershock as an earthquake occurring after a mainshock and in its epicentral region, although implying some causal relation with the mainshock, is partly semantic and largely circumstantial. Indeed, smallish earthquakes that constitute the background seismicityoccur all the time in a seismic region in the absence of large events, and continue occurring whether or not a large earthquake occurs, so that not all earthquakes occurring in the region after a mainshock are necessarily aftershocks.

If aftershocks are a result of the occurrence of the mainshock, then they should be related in a physical way with its rupture process. The aftershock-producing mechanism is not yet known, but it is conceivably related with adjustments to the post-mainshock stress field (Lay & Wallace, 1995) possibly through viscolelastic processes or through fluid flow (Nur & Booker, 1972); whichever the actual process, aftershocks should be related with the rupture plane of the mainshock. Kisslinger (1996) qualitatively defines three kinds of aftershocks: Class 1 is those occurring on the ruptured area of the fault plane or on a thin band around it. Class 2 is events that occur on the same fault but outside of the co-seismic ruptured area. Class 3 is events occurring elsewhere, on faults other than the one ruptured by the mainshock; these events, whether in the same region or not, will not be considered here as aftershocks, but rather will be classified as triggeredearthquakes.

The number of aftershocks decreases with time after the mainshock according to the modified Omori relation (Utsu, 1961, 1969, 1970) as


where tis time measured from that of the mainshock, and k, c, and pare positive parameters which vary with the lithologic, tectonic, and other conditions of the study area. Commonly, pranges from 0.9 to 1.8, with most instances between 1.1 and 1.4 (Utsu, 1961); values which do not show dependence on the magnitude of the mainshock (Utsu, 1962).

Aftershocks occurring within 24 to 48 hours after a large earthquake locate mostly over the co-seismic rupture area, and provide a good estimation of it (Lay & Wallace, 1995), which indicates that seismicity at the time is mostly class 1; over longer times the aftershock area increases (Felzer & Brodsky, 2006; Helmstetter & Sornette, 2002; Mogi, 1968; Tajima y Kanamori, 1985; Valdés et al., 1982), including at the edges class 2 events.

Among other reasons why aftershock identification is important, we can mention the following few examples. Aftershocks can give important information about the rupture area; also, from estimations of co-seismic slip on the fault plane by inversion of seismic waves, several authors have found that aftershocks are scarce in areas of maximum slip and concentrate around their edges (Dreger et al., 1994; Engdahl et al, 1989; Hauksson et al., 1994; Mendoza & Hartzell, 1988), so that aftershocks give information about the rupture process of the mainshock. Aftershocks can also yield information about the properties of the epicentral region (Knopoff et al., 1982; Kisslinger, 1996; Kisslinger & Hasegawa, 1991; Figueroa, 2009), and about possible triggering mechanisms (Roquemore & Simila, 1994).

Since aftershocks can be large enough to contribute to the damage (particularly after structures have been debilitated by a mainshock), it is important to estimate the hazard associated with aftershock activity (Felzer et al., 2003; Reasenberg & Jones, 1989).

It has been proposed that, since aftershock activity depends among other factors on the stress status of the region, there is research on whether some characteristics of the aftershock activity from intermediate-sized earthquakes can be useful as precursory data for large earthquake hazard estimation (e.g. Jones, 1994; Keilis-Borok et al., 1980; Wyss, 1986).

Finally, for some studies concerned with large earthquakes, aftershocks can be considered as noise, and have to be eliminated from the catalogs (e.g. Gardner. & Knopoff, 1974; Habermann & Wyss, 1984).

In order to use or eliminate aftershocks it is first necessary to identify them. Many methods have been used, ranging from visual inspection (Molchan & Dmitrieva, 1992) to sophisticated numerical techniques. A common method identifies aftershocks as those shocks locating within temporal and spatial windows having lengths which usually depend on the magnitude of the mainshock (Gardner & Knopoff, 1974; Keilis-Borok et al., 1980; Knopoff et al., 1982). More sophisticated methods identify aftershocks as belonging to a spatial cluster, consisting of events within a given distance of at least one other event belonging to the cluster, which includes the mainshock (Davis & Frohlich 1991a,b; Frohlich & Davis, 1985). A variation of the window method considers events from larger to smaller magnitudes with the size of the spatial windows a function of the magnitude, and density of events (Prozorov & Dziewonski, 1982; Prozorov, 1986). Other methods include recognizing some statistical property (e.g. Kagan & Knopoff 1976, 1982; Vere-Jones & Davies 1966) or interpreting the relations between events according to some statistical chain or branching model (Molchan & Dmitrieva, 1992; Reasenberg, 1985).

Our method includes some of the above mentioned techniques used to discard events which cannot be aftershocks, and then proceeds to identify aftershocks based on the physical model of a rupture plane and on recognized statistical relationships. An early unsophisticated application of the rupture plane model, which proved that this principle of aftershock identification was feasible, was part of an unpublished MSc. thesis (Granados, 2000).


2. The method

We work with seismic catalogs containing occurrence time (days), hypocentral x(East), y(North), and z(up), and, optionally, horizontal and vertical location uncertainties uhanduv, all these in kilometers. If location uncertainties are not included in the catalog, optional horizontal and vertical uncertainties are assigned equally to all events.

Any events occurring before the event with the largest magnitudeMmax, the mainshock, are eliminated. All spatial coordinates are then referred to those of the mainshock.

A rough time cutoff, eliminates events occurring after more than an optional cutoff time (default isnyr=4years), because after a few years it is difficult to distinguish aftershock activity from background seismicity.

The extent of the aftershock area depends on the energy, i.e. on the magnitude, of the mainshock (Utsu, & Seki, 1955), as does the rupture area. A first rough spatial discrimination, based on an average of the empirical magnitude Mvs. source-length r(km) relationships of Wells and Coppersmith (1994):


or (Kagan, 2004):


eliminates events farther away from the hypocenter than 1.5 times the rflength estimated by (2).

Next, a spatial clustering analysis, where events separated by no more than a given critical distance r, of the order of hundreds of meters to a few kilometers, depending on the spatial coverage of the catalog, are considered to be related, is used to eliminate events which do not relate to the mainshock or to other possible aftershocks.

The parameters in the modified Omori’s law (1) are not known a priori; they are estimated from the statistics of the aftershocks (Davis & Frohlich, 1991b; Guo & Ogata, 1997; Ogata, 1983), which we do not yet have. However, this relation tells us that for long enough times after the occurrence of the mainshock the number of aftershocks decreases until seismic activity in the epicentral region returns to its background level (Ogata & Shimazaki, 1984). When aftershock occurrence shows gaps comparable to those characterizing the background seismicity, we can consider that the aftershock activity is, if not ended, at least scarce enough to be comparable to the background activity and can no longer be distinguished from it. The critical gap length depends on the region and the magnitude threshold of the observations; we use a default critical gap lengthtgap=10days. The gap is measured as the average of a given number of inter-event times (defaultngap=10). Events occurring after a critical gap are discarded from the possible aftershocks.

Weights, based on relation (1), are optionally assigned to the remaining events, using typical values for c(defaultc2) and p(defaultp1.0) and by settingk=cpN(t=0)=1, so that those events which have a large likelihood of being aftershocks have weights ~1, while later events have smaller weights:


Next, plane fitting is carried out iteratively; at each iteration, a plane that passes through the mainshock hypocenter is fitted to all remaining events, through a genetic scheme described below, and fit outliers (events too far away from the plane) are discarded. Iteration continues until the goodness-of-fit criterion is met (successful fit) or until a preset maximum number of iterations is attained (unsuccessful fit).

Thetgap, c, and pvalues can be refined using the final results of a first, tentative aftershock determination, to do a second one.

The genetic plane search for the plane, characterized by its azimuth ϕand dipδ, which minimizes the L1 norm of the perpendicular distances from the fault plane to all events, is as follows. An initial set of equispaced ϕjand δjvalues, covering the acceptable ranges, is built and the fit to the aftershock candidates is evaluated for each pair of values.

To estimate the error of fit for each candidate plane, for each azimuth ϕjand dip δjpair, event coordinates (x,y,z)are rotated as:


and in the (x,y,z)coordinate system the equation of the plane isx=0, so that the L1 fit error for the j’th parameter pair is easily computed as


where Nais the number of remaining aftershock candidates.

The parameter pairs corresponding to the best Npfits are chosen as the parents of the next generation; they are sorted by the absolute value of their respective errors, so that parent number 1 is the best fit and parent number Nphas the largest errorεNp. The standard deviations of the parent’s parameters, σϕandσδ, are evaluated.

Next, Np(Np1)/2children are constructed with parameters which are averages of those of each pair of parents; other Ncchildren are constructed for each parent jby randomly modifying one or the other of the parameter values, according to


where Ndesignates the normal distribution, σϕm=max{σmin,σϕ}, σδm=max{σmin,σδ}, and σminis a minimum allowable value that ensures significant variations.

Errors are computed for the children and the Npbest fits among the whole population, parents plus children, are chosen as the parents for the next generation. The process is repeated until the goodness of fit criterion is met (and the process ends) or until a preset number of generations is attained.

For the current iteration the best fit is for the plane corresponding to parent number one, ϕ=ϕ1andδ=δ1; using these values event coordinates are rotated as in (4) and the standard deviation, σ, of the xdistances is evaluated. Those events with


where fσis a damping factor (default is 1.25) and uiis the location uncertainty, adjusted for the plane dip as


are eliminated as outliers.

This method has been implemented as a Matlab program, aftplane.m, which allows the user to interactively set most parameters, 3D plots the input hypocenters, identifies the aftershocks of one mainshock, optionally plots the aftershocks and/or independent shocks (in different colors), and optionally outputs the corresponding catalogs. A definite advantage of using Matlab for this algorithm is that both data and trial models are handled as matrices, so that rotating, sorting, and identifying values is done more efficiently and with less lines of code than would be possible in other programming languages like FORTRAN or BASIC.

A variation of this method is used as a function by program cleancat.mwhich identifies and eliminates the aftershocks of all mainshocks in a catalog. For each event not previously identified as an aftershock, aftershocks are identified following the same steps described above, except that, after fitting the plane a search is made for events with magnitudesMMmainΔM, which will not be considered as aftershocks, and if any are found, the aftershock list is cut so as to exclude all events beginning with the first of these non-aftershocks. Thus, for mainshocks occurring at timesti, instead of the total Omori number of aftershocks (Utsu, 1970; Ogata, 1983)


where H(t)is the Heavyside function, we are actually evaluating




Π(t)is the boxcar function, and tktiis the time of the succeeding event occurring along the same fault plane and cluster or the time when the gap criterion is met (tk=if no such event or gap occur). An event at tkmay “inherit” some of the aftershocks from the one atti, but they will be identified as aftershocks and be duly eliminated. Båth’s law indicates that ΔMshould be 1.2, but this criterion results, for most catalogs, in too many mainshocks; our default value is thusΔM=1.0, which means that we accept as aftershocks those with energy ratioE(M)/E(Mmain)~0.016, but this ΔMvalue can be easily changed by the user. The largest events in the catalog are plotted vs. time above the initial and resulting event densities, to illustrate which aftershocks were eliminated. The program iterates the whole process, as many times as needed, until no more aftershocks are found

In cleancat, parameters are not set interactively, but can be easily adjusted in a list of adjustable parameters at the beginning of the code.

The program optionally outputs a catalog excluding identified aftershocks.


3. Application

We will now show some examples of the application of the method. The aftplaneprogram will be used to identify fault planes and aftershocks from three mainshock-aftershock sequences from two different parts of the world featuring different tectonic environments. The cleancatprogram will be used to clean the catalog of a fault system.

3.1. Aftplane: transcurrent regime, Joshua Tree and Landers earthquakes

In 1992, two large earthquakes, the Joshua Tree, April 23 MW6.2, 33.9°N, 116.3°W and the Landers, June 28 MW7.3, 34.2°N, 116.5°W, occurred close together on a line previously unrecognized as a potential throughgoing seismogenic fault (Nur et al., 1993). We chose these events as illustration because, although both events have mainly strike-slip mechanisms, they have slightly different strikes and dips, so that we wanted to test whether the method could identify these small differences.

Figure 1 shows the location of the study area in souther California, USA, and its recent seismicity; the faults ruptured during the Joshua Tree and Landers earthquakes are located within the red diamond.

For both events we used maximum allowable horizontal and vertical location uncertainties of uhhmax=0.2 km anduvmax=0.5 km, respectively; a priory cutoff timesnyr=4years; maximum distance for spatial clustering r=1km; maximum permissible time gap tgap=10daysestimated as the average of ngap=10inter-event times; Omori weighting parameters c=2.0days andp=1.0; aftershock magnitude criterionΔM=1.0; fit criterion maximum error εmaxkm.

Figure 1.

Seismicity map of southern California showing the location of the Joshua Tree and Landers faults (both within the red diamond.) (Modified from Lin et al, 2007.)

3.2. Aftplane: Joshua Tree earthquake

The catalog for the Joshuea tree earthquake contained 5075 events spanning 66.30 days (~0.182 yr). Figure 2 (left) shows the 3497 remaining events after the first rough elimination by acceptable uncertainties and by an estimated expected fault length of ~7.97 km corresponding to a critical distance of 11.96 km. Figure 2 (right) shows as blue circles the 3379 shocks identified as clustering with the main event.

Figure 2.

Joshua Tree mainshock plus 3497 acceptable aftershock candidates (left) and Joshua Tree mainshock plus 3379 clustered aftershock candidates (right).

Figure 3.

Joshua Tree mainshock (red asterisk) plus 1094 aftershocks (blue diamonds); left: plan view showing 171.6º faultplane azimuth; right: view along faultplane azimuth clearly showing the 86.6º dip.

The main Joshua Tree shock and the identified 1094 aftershocks are shown in figure 3, both in a plan view (left) which clearly shows the resulting 171.6º faultplane azimuth, and a cross section along the fault plane azimuth (right) which shows the resulting 86.6º faultplane dip. The values found by aftplane agree extremely well with those estimated by Velasco et al. (1994) of strike 171º, dip 89°. Figure 4 shows a cross section parallel to the fault plane, illustrating aftershock concentrations.

Figure 4.

Joshua Tree mainshock (red asterisk) plus 1094 aftershocks (diamonds), cross section seen along azimuth 81.6º, perpendicular to the 171.6º faultplane azimuth.

3.2.1. Aftplane: Landers earthquake

The catalog for the Landers earthquake contained 49,605 events spanning 4,932.52 days (~13.514 yr). Figure 5 (left) shows the 17,553 remaining events after the first rough elimination by acceptable uncertainties and by an estimated expected fault length of ~76.78 km corresponding to a critical distance of 115.17 km. Figure 5 (right) shows as blue circles the 12,834 shocks identified as clustering with the main event.

Figure 5.

Landers mainshock plus 17553 acceptable aftershock candidates (left) and Landers mainshock plus 12834 clustered aftershock candidates (right).

The main Landers shock and the identified 3,225 aftershocks are shown in figure 6, in a cross section seen along the determined 340.6º fault plane azimuth, which shows the resulting 70.1º faultplane dip. The values found by aftplane agree extremely well with those estimated by Velasco et al. (1994) of strike 341°º, dip 70°. Figure 7 shows a cross section parallel to the fault plane, illustrating aftershock concentrations.

Figure 6.

Landers mainshock (red asterisk) plus 3225 aftershocks (blue diamonds), view along 340.6º faultplane azimuth. The location of the mainshock hypocenter is obscured by those of the aftershocks. Dip 70.1º

Figure 7.

Landers mainshock (red asterisk) plus 3225 aftershocks (blue diamonds), view along azimuth 70º, perpendicular to 340.6º faultplane azimuth.

3.3. Aftplane: subduction regime, Armería (Tecomán, Colima) earthquake

The Armería MW7.4 earthquake, also known as the Tecomán or Colima 2003 earthquake, occurred on January 22, 2003, on the subduction zone along the boundary of the Northamerican and Pacific plates in central-western Mexico.

Figure 8 shows the location of the study area, the mainshock epicenter (red star) and the subsequent seismicity recorded and located by the Colima Seismic Network (RESCO).

Nuñez et al (2004) and Yagi et al (2003) estimated a fault plane with a 300° strike and a quite shallow 20° dip, which agrees with the 20° to 30° dip of the subduction zone determined by Andrews et al (2010).

For aftplane we used the RESCO catalog with the same parameter values mentioned above, except for horizontal and vertical location uncertainties of uh=0.075 km anduv=0.50 km, and maximum distance for spatial clustering r=4km;

The catalog for the Armería earthquake contained 11,475 events spanning 1,529.9 days (~4.192 yr). Figure 9 (left) shows the 10,275 remaining events after the first rough elimination by acceptable uncertainties and by an estimated expected fault length of ~92.73 km corresponding to a critical distance of 139.09 km. Figure 9 (right) shows as blue circles the 7,109 shocks identified as clustering with the main event.

Figure 8.

Location of the Armería, 22 January 2003, MW 7.3 earthquake; the star indicates the epicenter of the main shock, circles are located events following the mainshock, located by the RESCO network.

The main Armería shock and the identified 460 aftershocks are shown in figure 10, in a cross section seen along the determined 86.2º fault plane azimuth, which shows the resulting 33.2º faultplane dip. The values found by aftplane agree extremely well with the above mentioned estimates strike 300° and 20° to 30° (Nuñez et al., 2004; Yagi et al., 2004; Andrews et al., 2010). Figure 11 shows a cross section parallel to the fault plane, illustrating aftershock concentrations.

Figure 9.

Armería mainshock plus 10868 acceptable aftershock candidates shown as black crosses (left) and Armería mainshock plus 7109 clustered aftershock candidates shown as blue circles (right).

Figure 10.

Armería mainshock (red asterisk) plus 460 aftershocks (blue diamonds), view along 86.2º faultplane azimuth.

Figure 11.

Armería mainshock (red asterisk) plus 460 aftershocks (blue diamonds), view along azimuth 356.2º, perpendicular to 86.2º faultplane azimuth.

3.4. Cleancat: whole Joshua Tree-Landers fault system

To illustrate the use of program cleancatwe chose the catalog covering the whole Joshua Tree-Landers fault system (Figure 1), because this is a system with many close-lying, subparallel, faults, which gives scope to the iterative aftershock recognition scheme of the program.

The parameters used are uhhmax=0.2 km anduvmax=0.5 km; a priory cutoff timesnyr=4years; maximum distance for spatial clustering r=2km; maximum permissible time gap tgap=30daysestimated as the average of ngap=10inter-event times; Omori weighting parameters c=2.0days andp=1.0; aftershock magnitude criterionΔM=1; fit criterion maximum error εmax=0.35km.

Figure 11 shows all events in the catalog (black crosses), and identified aftershocks as yellow circles. Total processing consisted of 10 iterations which identified and eliminated 11,665, 4,212, 1,702, 86, 94, 30, 80, 18, 49, and 1 aftershocks, respectively, for a total of 17,937 aftershocks.

Figure 12.

Joshua Tree- Landers fault system seismicity (black crosses) and identified aftershocks (yellow circles).

Figure 12 shows the occurrence times and magnitudes of the largest events in the catalog (top), and below them (middle) is plotted the ocurrence density (per Δt=46day) vs. time (blue line); peaks in the occurrence rate after large shocks aftershocks are clearly seen. The bottom panel of figure 12 shows the occurrence density after aftershock elimination using ΔM=1.0(red line); although densities are much lower than before filtering, the peaks coinciding with the occurrence of the Joshua Tree and Landers earthquakes indicate that many aftershocks are not being identified.

Use of ΔM=0.9(Fig.13) effectively diminishes the troublesome above mentioned peaks to background seismicity level; the seismicity rate peak aroundt~7.25105days is associated with distributed seismicity with events about the same size. Thus, we see that Båth’s law is not appropriate for the mainshock-aftershock relationships in the seismicity of the Joshua Tree-Landers fault system; a small (0.1 unit) change in the magnitude criterion can make a large difference in the aftershock recognition capability.

Figure 13.

Joshua Tree-Landers fault system seismicity versus time. The top panel shows the largest events in the period (blue circles with vertical lines). The middle and bottom panels show seismicity rates, for 46 day-long time intervals, before (middle) and after (bottom) processing by cleancat, respectively; note the different vertical scales.

Figure 14.

Joshua Tree-Landers fault system seismicity rates, for 46 day-long time intervals, after processing by cleancat, withΔM=0.9.


4. Conclusions

We present a simple method for identification and/or elimination of aftershocks, based on the generally accepted assumption that aftershocks are related to the fault rupture of the mainshock. The method has been tried on various catalogs with good results and, when aftershocks are numerous enough, good estimates of rupture planes that agree very well with those reported in the literature.

A variation of the method used for eliminating all aftershocks from a seismicity catalog (“catalog cleaning”) uses, iteratively, a variation of the same principle. Using the seismicity occurrence time rate as illustration and criterion of the effectiveness of the method, indicates that the required difference between mainshock and aftershocks ΔMis a key parameter for correct aftershock identification, and that ΔM=1.2(Båth’s law) may be too strict for some geographic areas and/or seismo-tectonic settings.



Many thanks to José Frez, Juan García A., and María Luisa Argote for useful criticism and comments. We are grateful to RESCO and Gabriel Reyes, and to the SCEC for the use of their catalogs. We also thank Mr Davor Vidic of Intech for the kind invitation to participate in the present book.

© 2011 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution-NonCommercial-ShareAlike-3.0 License, which permits use, distribution and reproduction for non-commercial purposes, provided the original is properly cited and derivative works building on this content are distributed under the same license.

How to cite and reference

Link to this chapter Copy to clipboard

Cite this chapter Copy to clipboard

F.A.Nava, V.H.Márquez and J.F.Granados (August 1st 2011). Aftershock Identification Through Genetic Fault-Plane Fitting, Scientific and Engineering Applications Using MATLAB, Emilson Pereira Leite, IntechOpen, DOI: 10.5772/23733. Available from:

chapter statistics

2458total chapter downloads

1Crossref citations

More statistics for editors and authors

Login to your personal dashboard for more detailed statistics on your publications.

Access personal reporting

Related Content

This Book

Next chapter

Sea Surface Temperature (SST) and the Indian Summer Monsoon

By S. C. Chakravarty

Related Book

First chapter

Tips and Tricks for Programming in Matlab

By Karel Perutka

We are IntechOpen, the world's leading publisher of Open Access books. Built by scientists, for scientists. Our readership spans scientists, professors, researchers, librarians, and students, as well as business professionals. We share our knowledge and peer-reveiwed research papers with libraries, scientific and engineering societies, and also work with corporate R&D departments and government entities.

More About Us