Earthquake Source Characterization Essay

1 Introduction

The regional patterns of macroseismic intensity, I, contain information about earthquake sources. Over the years, the problem of extracting this information has been addressed in different ways. Most authors have used empirical relationships between the areas of isoseismals, or sparse observations, and their barycenters to infer magnitudes and epicenters [e.g., Frankel, 1994; Johnston, 1996; Bakun and Wentworth, 1997; Gasperini et al., 1999; Musson and Jiménez, 2008] and treated the geographical distribution of high I degrees to infer the directions of the fault sources at depth [Gasperini et al., 1999, 2010]; the approach by Gasperini et al. was also implemented in the catalogue of seismogenic sources in Italy that were larger than M 5.5 [DISS Working Group, 2010].

Suhadolc et al. [1988], Zahradnik [1989], and Molchan et al. [2004], however, attempted to include the source kinematics in a model. This approach is interesting but requires reliable I patterns, which are rare for historic earthquakes. The nonlinear intensity-based automatic inversion technique KF niching genetic algorithm (NGA) [Sirovich and Pettenati, 2004], which is used in the present work, was developed to exploit the extraordinary database of high-quality macroseismic observations in the Mercalli-Cancani-Sieberg (MCS) scale [Locati et al., 2011]. These observations were made by seismologists and historians, who interpreted the reports that were written by the officers of various kingdoms of the time, the gazettes, and other sources. This KF-NGA algorithm roughly identifies the principal geometric and kinematic characteristics of the causative sources of destructive earthquakes from their regional I patterns (see section 2.1).

In previous papers, we performed geophysical nonlinear inversions by using a simple kinematic-function (KF) model [Sirovich, 1996] and a grid-search approach [Pettenati and Sirovich, 2003]. Since 2004, this kind of source inversion used a NGA algorithm. Because both the KF model and the NGA have been previously presented in this journal [Sirovich and Pettenati, 2004; Sirovich et al., 2013], only a brief summary of these methods is given here. In particular, the KF-NGA technique was verified by using data from earthquakes in California [Pettenati and Sirovich, 2007], the Fennoscandian Shield-Caledonian Range transition [Bungum et al., 2009] and Greece [Pettenati et al., 2010]. In California and Greece, the fault sources of earthquakes that were retrieved by KF-NGA inversion were confirmed by those that were obtained from instrumentation. Recently, the source mechanism of the 1980 Ms 6.9 Irpinia earthquake was inverted, and the results agreed well with instrumental data and field observations of the rupture [Sirovich et al., 2013]. In the same paper, the fault source of the 1694 Mw 7.0 earthquake was retrieved in the same region, and the result was confirmed by the paleoseismological study by Galli et al. [2014]. The vivid descriptions of the surface rupture on 8 September 1694 and its changes during the following weeks, which were recorded by anonymous high officers from the Kingdom of Naples (printed in 1694), prove the correctness of the inversion results [see Sirovich et al., 2013].

This paper shows that the KF-NGA source inversion of the I data by Guidoboni et al. [2007] for the Ferrara earthquake on 17 November 1570 (Mw 5.5 [Rovida et al., 2011]) retrieves a fault source that is compatible with the regional seismotectonic setting and present evolution of the outer front of the Apenninic chain. This result also explains the role of the earthquake in the definitive diversion of the Po River course at Ficarolo (20 km NW of Ferrara). The study area is shown in Figure 1, and the point intensities that are used in this work are shown in Figure 2. This diversion dried the former principal branch of the river's delta and moved it from the area of “Valli di Comacchio” to its present position 40 km to the north.

First, the basic kinematic information of the fault source of the destructive earthquake that hit the southern Po Plain on 17 November 1570 is retrieved (including the double-couple, or DC, mechanism; see Figure 3). This earthquake had similar kinematics to the 20 May 2012 (M 6.1) and 29 May 2012 events (M 5.9; both magnitudes from Ganas et al. [2012] [see Anzidei et al., 2012]). The three aforementioned events occurred at similar distances (20–30 km) en echelon (see Figures 4-6).

In a second step, the present paper reviews the existing geomorphological, archaeological, and historical information (Figure 4) on the regional uplift of the southern flank of the Po Plain over the last millennia. These pieces of data are intricate because the sedimentation and erosion rates, eustatic history, tectonic uplift, natural subsidence that is caused by consolidation and compaction and recent anthropogenic subsidence are superimposed onto the Po Plain. In this context, the hydrographic net is the most reliable reference, showing that the river course over the last millennia has shifted to the north by 20–40 km. The 1570–1591 seismic sequence appears to have been the last step in the tectonic process that caused the definitive diversion of the course of the Po River.

Incidentally, this geomorphological evidence solves the ±180° ambiguity of the KF-NGA algorithm in the calculated rake angle because the uplift of the southern Po riverbank is only compatible with a compressive mechanism of the Apenninic type. Therefore, we could draw the beach ball diagram in black and white; however, the diagram remains white to emphasize the limitation of the standard procedure that was used. This work continues a multiyear project to obtain quantitative seismological information on preinstrumental earthquakes from their regional I patterns.

2 Methods

Two algorithms are used: (i) the KF geophysical inversion of intensities, which applies the NGA genetic approach and (ii) the natural-neighbor (n-n) bivariate interpolation scheme, which contours point I data and draws isoseismals [Sirovich et al., 2002]. For the present purposes, we only recall here that the interpolation algorithm uses the n-n coordinates for weighting, interpolating, and contouring the I points. The interpolation is local because the weight of an observed I brought to a new neighbor point is proportional to the area of the intersection of their Voronoi polygons. In the n-n approach, the interpolant (a) fits the data exactly at the observation sites, (b) is isoparametric and bounded by the data values, and (c) is continuously differentiable at all points except the data sites. Moreover, the n-n isoseismals do not require (contouring) parameters and act as a compromise between the crude objectivity of the Voronoi tessellation and intuitive appeal of the somewhat subjective classical isoseismals. Thus, the values of the point intensities in this paper are given by the colors of the isoseismal areas.

2.1 The KF-NGA Inversion Procedure

The 11-parameter KF inversion procedure considers the total horizontal component of the body wave radiation from a rupture plane of unit width in an elastic half space within a distance range from approximately 10 to 80–100 km from the source. Incidentally, the nonparametric regression analyses of the attenuation of the peak ground acceleration (PGA) from moderate-sized earthquakes in the Po Plain by Bragato et al. [2011] demonstrated that the PGA in the present study area is systematically enhanced for distances between 80 and 200 km by postcritically reflected S waves and multiples from the Moho discontinuity. This situation has been suggested to hold for the PGA and pseudo-acceleration at various frequencies and for the intensities of strong earthquakes [Sugan and Vuan, 2013].

The 11 parameters that are inverted are the hypocentral latitude and longitude, the DC (strike, dip, and rake angle), the seismic moment (M0), the depth (H) of the line source, the shear wave velocity (Vs), the rupture velocities (Vr) along strike and antistrike, and the along-strike percentage of the total rupture length (L). The L and source width parameters are derived from the M0 by using the empirical relationships given by Wells and Coppersmith [1994]. Moreover, L is the sum of the absolute values of the along-strike (considered positive, L+) and the antistrike source segments (L−). The signs of the Vr/Vs ratios (i.e., Mach numbers) follow the same logic (Mach+ and Mach−).

The nondimensional values produced by the KF are transformed into intensities through the empirical relationship that was developed by Sirovich et al. [2001, equation (2)] from 1720 I values of Californian earthquakes. These 1720 observations were more often taken in sedimentary basins where people tend to live; consequently, the aforementioned empirical correlation between the KF and I values is likely to implicitly incorporate site effects that are similar to those in the Po Plain.

The inversion is performed by an NGA sharing technique with four independent subpopulations, which are used to search for the absolute minimum variance model in the 11-source parameter model space. Each initial subpopulation has 1000 random sources and evolves independently from the others. One selection and three evolutionary (stochastic) steps are applied within every generation. First, the best individuals (90% of the initial source models) are chosen in the Selection step and are then “married” in the Crossover step (with 90% fertility). In the Mutation step, 6% of the 500 new individuals (450 “sons” plus 50 nonfertile parents) receive a new value for one parameter by chance. In the Crowding step, identical individuals are deleted and each subpopulation is prevented from quickly converging to a false minimum. The genetic process continues for N generations to ensure that the virtual space of the parameters is properly explored. In particular, the best source of each generation that is transmitted to the new generation is included in the evolution logic. The NGA also involves the sharing step, which allows subpopulations of sources to survive within parameter subspaces (niches) so that each source of each subpopulation is not in competition with the sources of other subpopulations that live in other niches. In this step, each subpopulation of sources evolves independently from the others, and the normalized distance between each source of a subpopulation and each source of all the other subpopulations obeys a certain condition [Koper et al., 1999, equation (1)].

The objective function during the inversion is , where rs is the intensity calculated at a site minus the field intensity (obtained by historians and seismologists; the suffix denotes the sites). Each parameter space is explored with the sampling steps and within the ranges given in Table 1. The resolutions of the uncertainties of the parameters obtained by the inversion coincide with the widths of the steps shown in the table. The width of the fault source that is used by the KF is also obtained with the aid of Wells and Coppersmith [1994] (using empirical relationship number 4 in Table 2A on p. 990 and “Slip type” “All”) from the M0 via the relationship Mw = 2/3(log M0 − 9.1), which was suggested by International Association of Seismology and Physics of the Earth's Interior [2005], where Mw is the moment magnitude. The inversion errors are calculated by inverting randomized sets of data and following the standard by Pettenati and Sirovich [2007] and Sirovich et al. [2013] (see later).

Nucleation latitude N (deg)44.50–45.200.01
Nucleation longitude E (deg)11.20 – 12.000.01
Strike angle (deg)0 – 3591
Dip angle (deg)20 – 901
Rake angle (deg)0 – 1791
Nucleation depth (km)4 – 500.1
Vs (km/s)3.50 – 3.960.01
Mach in along-strike direction0.50 – 0.950.01
Mach in antistrike direction0.50 – 0.950.01
M0 (N m) 10171.0 – 8.00.01
Percentage of L along-strike1 – 1001
Nucleation latitude N (deg)44.97 ± 0.0744.97 ± 0.08
Nucleation longitude E (deg)11.55 ± 0.1011.63 ± 0.09
Strike angle (deg)121 ± 16127 ± 16
Dip angle (deg)26 ± 628 ± 7
Rake angle (deg)73 ± 1877 ± 16
Nucleation depth (km)35.3 ± 4.834.9 ± 5.4
Vs (km/s)3.96 ± 0.06


In this thesis, I present a series of works on the characterization of source properties and physical mechanisms of various small to moderate earthquakes through both observational and numerical approaches. From the results, we find implications on a broader scheme of topics relating to larger earthquakes, shear zone structure, frictional properties of faults, and seismic hazard assessment.

Part I consists of two studies using waveform modeling. In Chapter 2, we present an in-depth study of a series of intraslab earthquakes that occurred in a localized region near the downdip edge of the 2011 Mw Tohoku-Oki megathrust earthquake. By refining source parameters of selected events, simulating their rupture properties and comparing their mechanisms to stress changes caused by the main shock in the region, we are able to identify the true rupture plane and the reactivation of a subducted normal fault, enhancing our understanding on the downdip shear zone. In Chapter 3, based on similar techniques, we further develop a systematic methodology to perform fast assessments on important source properties as an earthquake occurs. For two Mw 4.4 earthquakes in Fontana, moment magnitude and focal mechanism can be accurately estimated with 3 to 6 s after the first P-wave arrival, while focal depth can be constrained upon the arrival of S waves. Rupture directivity can also be determined with as little as 3 seconds of P waves. This study opens the opportunity to predict ground motions ahead of time and can potentially be useful for Earthquake Early Warning.

Part II involves the modeling of seismic source properties and physical mechanisms of interacting earthquakes in dynamic rupture simulations. In particular, we focus on small repeating earthquake sequences that trigger one another. In Chapter 4, we quantify the relative importance of physical mechanisms that contribute to earthquake interaction and identify that the stress change caused by post seismic slip is the dominating factor. Our findings introduce the possibility to constrain frictional properties of the fault based on earthquake interactions. We further apply this working model in Chapter 5 to reproduce the actual interacting repeating sequences in Parkfield. We are able to identify possible physical mechanisms that cause the inferred high stress drops of these repeating events, as well as reproduce their synchronized seismic cycles. Results from our simulations are consistent with the observed scaling relation between the recurrence time interval and the seismic moment of these events. Our findings indicate that the difference between the observed and the theoretical scaling relations can be explained by the significant aseismic slip in the rupture area.

Item Type:Thesis (Dissertation (Ph.D.))
Subject Keywords:Seismology, earthquakes, waveform, seismic rupture, numerical simulation
Degree Grantor:California Institute of Technology
Division:Geological and Planetary Sciences
Major Option:Geophysics
Thesis Availability:Public (worldwide access)
Research Advisor(s):
  • Lapusta, Nadia (co-advisor)
  • Helmberger, Donald V. (co-advisor)
Thesis Committee:
  • Avouac, Jean-Philippe (chair)
  • Clayton, Robert W.
  • Tsai, Victor C.
  • Ampuero, Jean-Paul
  • Lapusta, Nadia
  • Helmberger, Donald V.
Defense Date:26 October 2016
Record Number:CaltechTHESIS:11232016-044435496
Persistent URL:
Related URLs:
Default Usage Policy:No commercial reproduction, distribution, display or performance rights in this work are provided.
ID Code:9983
Deposited By:Ka Yan Semechah Lui
Deposited On:29 Nov 2016 17:03
Last Modified:12 Dec 2016 23:31


PDF - Final Version
See Usage Policy.

Repository Staff Only: item control page


Leave a Reply

Your email address will not be published. Required fields are marked *