*Article* **Use of Macroseismic Intensity Data to Validate a Regionally Adjustable Ground Motion Prediction Model**

### **Yuxiang Tang 1,\*, Nelson Lam 1,2, Hing-Ho Tsang 2,3 and Elisa Lumantarna 1,2**


Received: 2 September 2019; Accepted: 27 September 2019; Published: 30 September 2019

**Abstract:** In low-to-moderate seismicity (intraplate) regions where locally recorded strong motion data are too scare for conventional regression analysis, stochastic simulations based on seismological modelling have often been used to predict ground motions of future earthquakes. This modelling methodology has been practised in Central and Eastern North America (CENA) for decades. It is cautioned that ground motion prediction equations (GMPE) that have been developed for use in CENA might not always be suited for use in another intraplate region because of differences in the crustal structure. This paper introduces a regionally adjustable GMPE, known as the component attenuation model (CAM), by which a diversity of crustal conditions can be covered in one model. Input parameters into CAM have been configured in the same manner as a seismological model, as both types of models are based on decoupling the spectral properties of earthquake ground motions into a generic source factor and a regionally specific path factor (including anelastic and geometric attenuation factors) along with a crustal factor. Unlike seismological modelling, CAM is essentially a GMPE that can be adapted readily for use in different regions (or different areas within a region) without the need of undertaking any stochastic simulations, providing that parameters characterising the crustal structure have been identified. In addressing the challenge of validating a GMPE for use in an area where instrumental data are scarce, modified Mercalli intensity (MMI) data inferred from peak ground velocity values predicted by CAM are compared with records of MMI of past earthquake events, as reported in historical archives. South-Eastern Australia (SEA) and South-Eastern China (SEC) are the two study regions used in this article for demonstrating the viability of CAM as a ground motion prediction tool in an intraplate environment.

**Keywords:** seismic hazard; attenuation; GMPE; crustal model; MMI; intraplate region

#### **1. Introduction**

A ground motion prediction equation (GMPE) is a set of algebraic functions of earthquake magnitude, distance, and a site parameter, and may also include parameters to identify the style of faulting [1–5]. GMPEs are used to define the characteristics of ground motions for specific regions and earthquake scenarios (expressed in terms of magnitude-distance, or M-R combinations). In tectonically active regions (e.g., Western North America) where plenty of strong motion data can be captured by a network of densely distributed recording instruments, GMPEs are typically developed from regression analysis of recorded strong motion data. Empirical GMPEs, which were derived mainly from field recordings, should be capable of capturing regional specific earthquake ground motion characteristics [6,7]. However, empirical GMPEs may give results that are very sensitive to data from isolated records and the type of regression techniques adopted, and more so when data are scarce [8].

For tectonically stable regions of low-to-moderate seismicity, the conventional approach based on the regression of empirical data may not be feasible for developing GMPEs. To overcome the challenge of paucity in strong motion data, stochastic simulations of a seismological model expressed in the form of Fourier amplitude spectrum (FAS) can be undertaken to generate synthetic accelerograms [9,10]. This approach basically makes use of a seismological model to define the amplitude of the individual sinusoids constituting the acceleration time histories on the ground surface, whilst assigning random phase angles to the sinusoids. This approach of synthesizing ground motions is known as stochastic simulations. Standard calculation procedures may then be applied to derive the intensity of ground motion for any given earthquake scenario. A set of GMPE for the target regions can then be developed using the synthetic data obtained from simulations applied in a repetitive fashion to cover for different M-R combinations.

Complexity in the earthquake generation process and the associated uncertainties that are embodied in empirical GMPEs cannot possibly be captured completely by seismological modelling, nor by stochastic simulations. To address these intrinsic deficiencies of seismological modelling, a new class of (semiempirical) modelling approaches, namely the hybrid empirical method (HEM) [11,12] and the referenced empirical method (REM) [13,14], have been developed. However, the application of either HEM or REM requires a suitable "host region" containing abundant strong motion recordings and well-developed GMPEs.

In intraplate regions like South-Eastern Australia (SEA) and South-Eastern China (SEC), it is difficult to find a suitable "host region", or alternative well-established GMPEs, for applying HEM or REM. However, very useful macroseismic intensity information expressed in terms of the modified Mercalli intensity (MMI) has been recorded on isoseismal maps for earthquakes that occurred in this region for over hundreds of years. This type of data has been used in investigations for studying the attenuation behaviour of earthquake ground motions in Australia [15,16] and in many other regions of low-to-moderate seismicity [17]. The analysis of MMI data can result in the development of a GMPE that can be representative of the strong motion transmission properties of earthquake-affected areas based on observations over a long-time span (without requiring recording instruments to be placed close to the epicentre of any earthquake). For this reason, MMI data that have been recorded to date have much better coverage of strong motion conditions than instrumented data in regions of low-to-moderate seismicity, both in terms of time and space. Many current GMPEs developed in regions of low-to-moderate seismicity for the modelling of seismic hazards were derived from datasets of small magnitude earthquake events (for example, GMPEs developed in SEA were developed from datasets with a maximum moment magnitude M = 5.4 [18,19]). These GMPEs can be heavily biased to ground motion behaviour typifying small magnitude earthquakes. Their applicability in modelling large magnitude earthquake events occurring in the future is in doubt.

Another modelling challenge is reconciling modern GMPEs, which typically make use of peak ground acceleration (PGA), peak ground velocity (PGV), or response spectral accelerations (RSA) as intensity measures, with a GMPE derived from MMI data. Converting macroseismic intensity MMI data into ground motion parameters such as PGA or PGV has been studied by many scholars in the past few decades. A simple, and well-known, MMI-PGV conversion relationship was recommended by Newmark and Rosenblueth [20], Gaull et al. [21], and Lam et al. [15]. More rigorous conversion relationships have also been developed more recently (e.g., [22–24]). Regression analysis of MMI data, as obtained from historical archives, can therefore be transformed into attenuation relationships that are expressed in terms of PGA, PGV, or RSA.

In the absence of detailed spectral properties of a historical earthquake (which could only be derived from analyses of instrumented records), PGV is the preferred ground motion parameter, as opposed to PGA or RSA, for characterising the intensity of the earthquake for reasons outlined by Bommer and Alarcon [25]. First, PGV is a reliable, simple, and measurable intensity metric, which has a good correlation with damage distribution [26]. The accuracy of "shake maps" relies on good correlations between PGV and MMI values [22]. Second, PGV data can be employed for

estimating the risk of damage to buried pipelines because of the good correlations between horizontal PGV and material strains [27]. Thus, fragility functions for buried pipelines expressed in terms of PGV can be found in open sources [28]. Third, PGV data can also be used for assessing the risks of liquefaction [29–31]. Fourth, PGV is one of the three parameters for scaling an elastic response spectrum model for engineering applications [32,33]. The seismic action model stipulated by the current Australian standard is based on scaling design PGV values that were derived from probabilistic seismic hazard analysis (PSHA), employing GMPEs found on MMI data [34]. This is evident of PGV being recognised as the key ground motion intensity measure.

This paper is aimed at introducing a generic, regionally adjustable GMPE (known as the component attenuation model with acronym: CAM) for the prediction of PGV for given earthquake scenarios in low-to-moderate seismicity regions. The use of macroseismic historical MMI data collected from SEA and SEC regions (the distribution of magnitude and distance for collected data can be found in Figure 1) to demonstrate the accuracy of predictions by CAM is also presented.

**Figure 1.** M-R combinations for selected regions. (**a**) South-Eastern Australia (SEA); (**b**) South-Eastern China (SEC).

#### **2. The Component Attenuation Model (CAM) of PGV**

The framework of CAM as a GMPE is defined in a generic functional form, which decouples the source and attenuation effects into different components, as shown by Equation (1).

$$\text{PGV} = \Delta \times \alpha \times \beta \times \mathbf{G} \times \gamma\_c \times \mathbf{C} \tag{1}$$

Equation (1) can also be presented in the logarithmic (base 10) format as shown by Equation (2):

$$
\log \text{PCV} = \log \Delta + \log \alpha + \log \beta + \log \mathbf{G} + \log \gamma\_c + \log \mathbf{C}, \tag{2}
$$

where Δ is the referenced PGV measured on reference hard rock sites (typifying Eastern North America) for the reference scenario of M = 6 and R = 30; α is the source factor, which is a function of the moment magnitude (M) and Brune stress drop (Δσ); β is the regional anelastic attenuation factor, which is responsible for the attenuation effects, excluding geometric attenuation (G); G is the geometric attenuation factor; γ*<sup>c</sup>* is the crustal modification factor, which accounts for the amplification and attenuation phenomena within the upper 4 km of the rock crust, as well as amplification at the "mid-crust", which is in the surrounding of the location of the fault plane; and C is the calibration factor (which is used to minimise discrepancies between the model predictions and empirical recordings). The value of α is set as unity (and value of β and G close to unity) for a stress parameter value of 200 bars at M = 6 and R = 30 km. Seismological parameter values of stochastic simulations that CAM is

based upon are summarised in Table 1. The listing of the values of model coefficients associated with each component factor in CAM are presented in the next section.

**Table 1.** Parameter values used in stochastic simulations for component attenuation model (CAM) modelling.


<sup>1</sup> The equation is only used when the exponent value is unknown.

CAM, as introduced in the foregoing, is similar in methodology and philosophy to a generic ground motion prediction model (GMPE) introduced in reference [37]. A unique feature of CAM, which is not found in any GMPE, is the use of a geology-based crustal modelling approach, wherein crustal shear wave velocity profiles of bedrock to depths of tens of kilometres are used to derive modification factors of the upper crust (not to be confused with the site factor for modelling the effects of surficial sediments down to tens of metres only). Description of the crustal modelling methodology, complete with case studies, can be found throughout the rest of the article.

#### *2.1. Generic Source Factor*

The generic source factor is expressed in the form of Equation (3).

$$
\alpha \log \alpha = a\_1 \times \mathcal{M}^{a\_2} \times \Delta \sigma^{a\_3} + a\_{4\prime} \tag{3}
$$

where *a*1–*a*<sup>4</sup> are model coefficients. Equation (3), which was derived from the simulated data, has been normalised at M = 6. The model covers the range of M4–M8 and Δσ = 30–300 bar.

#### *2.2. Regional Whole Path Anelastic Attenuation Factor*

The regional anelastic attenuation factor (β) is used to account for the attenuation of ground motions along the entire transmission path of the seismic wave from source to site. Results derived from stochastic simulations of the seismological model with increasing distance have been normalised at R = 30 km and the effects of geometric attenuation effect (G) have been removed (and accounted for separately). The anelastic attenuation factor so derived is expressed in the form of Equation (4):

$$\log \beta = (b\_1 \times \mathcal{M} + b\_2) \times Q\_0^{b\_3} \times (\log \mathcal{R})^{b\_4} + b\_5. \tag{4}$$

where *b*1–*b*<sup>5</sup> are model coefficients; M is moment magnitude, *Q*<sup>0</sup> is the regional dependent quality factor for wave transmission, and R is the hypocentral distance.

Higher frequency waves (with a larger number of wave cycles for a given distance) are more susceptible to anelastic attenuation along the wave travel path than lower frequency waves. Given that the frequency characteristics of an earthquake are magnitude-dependent, the extent of anelastic attenuation is accordingly magnitude-dependent. Hence, earthquake magnitude is a parameter in the predictive expression of Equation (4). The reliability (ability) of this expression to reproduce results of stochastic simulations of the upper crustal amplification phenomenon reasonably accurately is demonstrated in Section 3.

Predictions derived from Equations (3) and (4) have been subject to residual analysis, wherein the residual value is defined as δ = log Ysim/Ypred , where Ysim and Ypred are PGV measurements obtained from stochastic simulations and CAM predictions, respectively. δ > 0 refers to underestimation of CAM and δ < 0 refers to overestimation of CAM. A fourth order polynomial expression (Equation (5)) for defining the value of the adjustment factor βadjustment\_factor is used to adjust, through multiplication, predicted values of β from Equation (4) to minimise the values of the residuals, along with any other systematic modelling errors.

$$
\beta\_{\text{adjusted\\_factor}} = b\_6 \mathbf{R}^4 + b\_7 \mathbf{R}^3 + b\_8 \mathbf{R}^2 + b\_9 \mathbf{R} + b\_{10} \tag{5}
$$

where *b*6–*b*<sup>10</sup> are model coefficients.

#### *2.3. Crustal Modification Factor*

The crustal factor γ<sup>c</sup> is to account for the combined effects of the amplification and attenuation of the upper earth crust, which is mainly dependent on the shear wave velocity profile in the upper 3–4 km of the earth crust. The literature refers to those phenomena as crustal modifications. The modification factor can be resolved into three components: (i) amplification of the upper crust; (ii) attenuation of the upper crust; and (iii) modification of the mid-crust. The respective factors representing each of these components are combined in a multiplicative manner, as represented by Equation (6).

$$
\gamma\_c = \gamma\_{am} \times \gamma\_{an} \times \gamma\_{mc\prime} \tag{6}
$$

where γ*am*, γ*an*, and γ*mc* are factors representing amplifications of the upper crust, attenuation of the upper crust, and modification of the mid-crust, respectively.

Details of the derivation of each of these component factors contributing to crustal modifications are described in the rest of this section under separate subheadings.

#### 2.3.1. Upper-Crustal Amplification

In modelling upper crustal amplification, information presented in reference [38], which is abbreviated herein as BJ97, and more recent updates of the model [39], have been incorporated for constructing the shear wave velocity profiles, which have also been used to infer the crustal density profiles [40]. The two profiles were then called up jointly to derive the upper crustal amplification factors by use of the square-root-impedance (SRI) method. The shear wave velocity profiles and the corresponding frequency-dependent amplification factors, so derived for the study regions, are presented in Figures 2 and 3.

The upper-crustal amplification factor γ*am* was derived accordingly as a function of VS30 (time-averaged shear wave velocity in the upper 30 m of the earth crusts [41]) by curve-fitting simulated data. The general form of the function is shown by Equation (7):

$$\log \mathbf{y}\_{\text{am}} = \boldsymbol{\gamma}\_1 \times \mathbf{M}^{\boldsymbol{\gamma}\_2} \times \mathbf{V}\_{\text{S30}}^{\boldsymbol{\gamma}\_3} + \boldsymbol{\gamma}\_4 \times \mathbf{V}\_{\text{S30}} \tag{7}$$

where γ1–γ<sup>4</sup> are the regression coefficients.

**Figure 2.** Shear-wave velocity profile: (**a**) SEA; (**b**) SEC, in which CRUST1.0 is a global crustal database (https://igppweb.ucsd.edu/~gabi/crust1.html, last accessed in January 2019), the shear wave velocity data can be found in the "Supplementary Materials".

**Figure 3.** Frequency-dependent crustal modification factor. For both conditions, <sup>κ</sup><sup>0</sup> <sup>=</sup> 0.057/Vs0.030.8 <sup>−</sup> 0.02 [42,43].

Upper crustal amplification of seismic waves is also frequency-dependent, given that the amplification of higher frequency waves (of shorter wave lengths) are more sensitive to the shear wave velocity gradient of the earth crust than lower frequency waves (of longer wave lengths). This phenomenon can be explained by reference to the quarter wavelength principles in the analysis for crustal amplification. Further explanations can be found in reference [38]. The extent of upper crustal amplification is accordingly dependent on the frequency characteristics of the upward propagating seismic waves. Given that the frequency characteristics of an earthquake are magnitude dependent, the extent of upper crustal amplification is, therefore, also magnitude dependent. Hence, earthquake magnitude is a parameter in the predictive expression of Equation (7). The validity of using VS30 to characterise upper crustal conditions has been demonstrated in reference [42]. The reliability (ability) of Equation (7) to reproduce results of stochastic simulations of the upper crustal amplification phenomenon reasonably accurately is demonstrated in Section 3.

#### 2.3.2. Upper Crustal Attenuation

The Kappa factor (κ0) is another key parameter to be considered in upper crustal modelling. The attenuation behaviour of the ground motion in the upper crust at high frequencies (in addition to the whole path anelastic attenuation covered in Section 2.2) is controlled by the value of κ0. Numerous studies that were targeted at modelling the value of κ<sup>0</sup> have been reported in the literature [41–46]. In simulations undertaken in this study, the value of κ<sup>0</sup> accordingly ranges from 0.001 s to 0.1 s. The functional form for the expression for determining the values of γ*an* to account for the effects of upper crustal attenuation is represented by Equation (8).

$$
\log \gamma\_{an} = \gamma\_5 \times \mathbb{M}^{\gamma\_6} \kappa\_0^{\gamma\_7} + \gamma\_{8\prime} \tag{8}
$$

where γ5–γ<sup>8</sup> are regression coefficients.


**Table 2.** Summary of the selected NGA-East seismological models.

<sup>1</sup> source models presented in the original references have been replaced by the more updated source model presented in the table based on recommendations in reference [35].2 R is hypocentral distance in km.

As for whole path anelastic attenuation, upper crustal attenuation is also magnitude-dependent. Thus, earthquake magnitude is also a parameter in Equation (8) for predicting the extent of upper crustal attenuation. The significance of the magnitude term is well demonstrated in Section 3.

Residual analysis has been conducted for CAM in totality, incorporating the generic source, regional path (including anelastic and geometric attenuation factors), and crustal modification factor. An apparent trend with increasing distance R has been identified. A fourth order polynomial expression (Equation (9)) for defining the value of the adjustment factor γadjustment\_factor is used to adjust, through multiplication, predicted values of PGV from Equation (3), (4), (5), (7), and (8) to minimise values of the residuals, along with other systematic modelling errors.

$$\gamma\_{\text{adjusted\\_factor}} = \gamma\_{\text{g}} \mathbf{R}^4 + \gamma\_{10} \mathbf{R}^3 + \gamma\_{11} \mathbf{R}^2 + \gamma\_{12} \mathbf{R} + \gamma\_{13\star} \tag{9}$$

where γ9–γ<sup>13</sup> are regression coefficients.

#### 2.3.3. Mid-Crustal Modification

Another crustal factor to consider is the mid-crustal modification factor γ*mc*, which is used to account for the effects of the density and shear wave velocity of the earth crust at the depth of the source of the earthquake. Mid-crustal amplification is purely a source phenomenon, as it represents the increase in the amplitude of shear waves (generated at the source of the earthquake) with decreasing shear wave velocity and density of the earth crust at the depth of the source from where seismic waves are emitted. In calculating the amplification factor a high shear wave velocity value of 3.8 km/s and a crustal density value of 2.8 g/cm3, both of which are characteristics of hard rock conditions in the shield regions of Eastern North America, are used as the "benchmark" conditions, for which the amplification factor is set as unity. These benchmark parameters for shear wave velocity and crustal density are denoted as β0sim and ρsim, respectively.

The value of γ*mc* can be found using Equation (10), which has been derived by the authors based on the source model defined in reference [35].

$$\gamma\_{mc} = \frac{\rho\_{\rm sim}}{\rho\_{\rm S}} \times \left(\frac{\beta\_{\rm Osim}}{\beta\_{\rm OS}}\right)^{k} \tag{10}$$

where *<sup>k</sup>* = <sup>−</sup>0.273 <sup>×</sup> <sup>M</sup> + 3.278, <sup>ρ</sup>sim = 2.8 g/cm3, <sup>β</sup>0sim = 3.8 km/s, and <sup>ρ</sup><sup>S</sup> and <sup>β</sup>0S are the density and shear wave velocity of the earth crust at the depth of the source (i.e., the mid-crust).

#### *2.4. Verification of CAM Using Various Seismological Models*

To verify CAM as a tool for translating seismological models into GMPEs for predictions of PGV, four NGA(Next Generation Attenuation)-East seismological models (as listed in Table 2), namely Atkinson and Boore (1995) [1], Silva, Gregor, and Darragh (2002) [47], Atkinson (2004) [48], and Boatwright and Seekins (2011) [49], have been used to verify CAM. The acronyms of these four models are AB95, SGD02, A04, and BS11. The PGV values, so derived from stochastic simulation of the (selected) seismological models, were compared with PGV values predicted by use of CAM for the respective model. Results of the verification analyses are presented in Section 3.2.

#### *2.5. Validation of CAM Using MMI Data*

Historical macroseismic MMI data collected from South-Eastern Australia (SEA) and South-Eastern China (SEC) has also been adopted to validate CAM as a regionally adjustable generic GMPE. The magnitude and distance distributions are shown in Figure 1 for the two selected study regions. The well-known MMI-PGV correlation proposed by Atkinson and Kaka [23] was adopted in this study (with residual corrections) to transfer the predicted PGV values (as derived from CAM) into MMI values for comparison purposes. The correlation relationship without residual corrections is expressed by Equation (11), and that with residual corrections is expressed by Equation (12). The historical MMI data can be found in this paper in the section presenting "Supplementary Materials".

$$\text{MMI} = \begin{cases} 4.37 + 1.32 \times \log \text{PGV} & \log \text{PGV} \le 0.48\\ 3.54 + 3.03 \times \log \text{PGV} & \log \text{PGV} > 0.48 \end{cases} \tag{11}$$

$$\text{MMI} = \begin{cases} 4.37 + 1.32 \times \log \text{PGV} + 0.47 - 0.19 \times \text{M} + 0.26 \times \log \text{R} & \log \text{PGV} \le 0.48\\ 3.54 + 3.03 \times \log \text{PGV} + 0.47 - 0.19 \times \text{M} + 0.26 \times \log \text{R} & \log \text{PGV} > 0.48 \end{cases} \tag{12}$$

where M and R are the moment magnitude and hypocentral distance, respectively, and the unit of PGV is cm/s.

For recordings collected from SEA, some of the events were recorded in terms of local magnitude (ML). The conversion between M (moment magnitude, same as MW) and ML would need to be undertaken in the first place to obtain correct estimates of the event magnitude. In this study, the bilinear conversion relationship developed for Australian conditions [50] as defined by Equation (13) was adopted.

$$\mathbf{M} = \begin{cases} 2/3 \mathbf{M}\_{\mathrm{L}} + 1.2, & \mathbf{M}\_{\mathrm{L}} \le 4.5 \\ \mathbf{M}\_{\mathrm{L}} - 0.3, & \mathbf{M}\_{\mathrm{L}} > 4.5 \end{cases} \tag{13}$$

As the earthquake recordings for SEC were derived originally from the ancient yearbook, some records dated back to as early as the 15th century [51,52]. Thus, magnitude conversion is filled with uncertainties. The magnitude recordings compiled in this study (which can be found in the section titled "Supplementary Materials") are expressed in terms of moment magnitude.

Another important component in CAM is the upper crustal modification factor. Regional shear wave velocity (SWV) profiles can be used for deriving the upper crustal amplification factor. The shear wave velocity profiles in this study were constructed based on the use of a geology-based modelling approach, as recommended in reference [53]. With this modelling approach, a compressional (P) wave velocity profile is converted into a shear (S) wave velocity profile using relationships that have been developed by regression analysis of recorded data collated from multiple sources. More details in relation to the process of modelling the shear wave velocity profile for the two study regions can be found in the "Supplementary Materials" section. The parameter values used for constructing SWV profiles are listed in Table 3. In Table 3, the values of ZS (depth of the upper sedimentary crustal layer) and ZC (combined thickness of the soft and hard sedimentary crustal rock layers) for each region were obtained from the average estimates of the thickness of soft sediment layer and total sediment layers, respectively, in CRUST1.0 database (https://igppweb.ucsd.edu/~gabi/crust1.html, last accessed in January 2019). VS0.03 (shear-wave velocity values at the depth of 0.03 km) values were identified by curve-fitting to minimise discrepancies (defined as sum of squares errors) between the modelled and recorded SWV values. The SWV value at 8 km depth (VS8), representing conditions at the source of the earthquake, has also been determined for the two study regions based on information of the regional crustal structure. The detailed information about the recorded SWV value can be found the references [54–59]. The velocity profile showing the upper bound of 2.78 km/s in Figure 2 is the shear wave velocity profile for generic hard rock conditions, as recommended in reference [38]. The velocity profile showing the lower bound of 0.76 km/s in Figure 2 was derived from interpolation between the shear wave velocity profiles for generic rock and generic hard rock conditions, as recommended in reference [38]. The modelling approach introduced in reference [39] has been adopted. The frequency-dependent modification factors for the study regions are shown in Figure 3.


**Table 3.** Parameter values for modelling shear wave velocity (SWV) profiles for South-Eastern Australia (SEA) and South-Eastern China (SEC).

1. VS0.03 is the shear wave velocity value at the depth of 0.03 km; 2. VS8 is the shear wave velocity value at the depth of 8 km; 3. VS0.2 is the shear wave velocity value at the depth of 200 m; 4. ZS is the depth of the upper sedimentary crustal layer; 5. ZC is the combined thickness of the soft and hard sedimentary crustal rock layers.

The parameter values used in CAM and the selected GMPEs for comparison purposes for SEA and SEC regions are summarised in Table 4.


**Table 4.** Parameter values used in CAM for SEA and SEC alongside the selected GMPEs.

<sup>1</sup> Both Δ values for SEA and SEC were determined from the calculated PGV at M6R30 from program GENQKE for hard rock conditions [9,15]; <sup>2</sup> Vs30 values are based on Figure 2. <sup>3</sup> PGVS and PGVR refer to PGV value on an average soil site and rock site, respectively, the value of 1.5 for the rock to the soil site conversion is based on stipulations by the earthquake loading standard for sites on shallow soil sediments [63] as explained in reference [64]; <sup>4</sup> SGC refers to Somerville et al. (2009) [65]; 5A12 refers to Allen (2012) [18]; <sup>6</sup> CB08 refers to Campbell and Bozorgnia (2008) [66]; <sup>7</sup> CY08 refers to Chiou and Youngs (2008) [2].

#### **3. Results of Verification Analyses**

#### *3.1. PGV Modelling*

Values for each of the component factors in CAM are presented in, Figures 4–8, alongside results generated from stochastic simulations of the respective seismological model. Figure 4 shows the source factor (α) as a function of M and Δσ.

Figure 5 shows the path factor (β) as a function of moment magnitude M, R, and wave transmission quality factor (Q0). Figure 6 shows the upper crustal amplification factor (γ*am*) as a function of M and VS30. Figure 7 shows the upper crustal attenuation factor (γ*an*) as a function of M and κ0. Figure 8 shows the overall PGV obtained both from stochastic simulations and CAM predictions. The regression coefficients together with the regression goodness (R2) values for different component factors demonstrating excellent agreement between the two sets of results are summarised in Table 5.

**Figure 4.** Simulated values of peak ground velocity (PGV) (symbols) alongside predictions by CAM (lines) for varying moment magnitudes (M4.5–M7.5) and stress drops (30–300 bar), at R = 30 km.

**Figure 5.** Anelastic attenuation factor (β). (**a**) Simulated PGV values normalised at Q0 = 680 for varying M values (M4.5–M7.5) and distances (4–800 km) shown alongside predictions by CAM (lines); (**b**) simulated PGV values normalised at M = 6 for varying Q0 values (150–800) and distances (4–800 km) shown alongside predictions by CAM (lines).

γ

γ

**Figure 6.** Regression analysis results of upper crustal amplification factor (γ*am*).

**Figure 7.** Regression analysis results of upper crustal attenuation factor (γ*an*).

**Figure 8.** Simulated PGV values (symbols) for varying M values (M4.5–M7.5) and distances (4–800 km) shown alongside predictions by CAM (lines) for Δσ = 200 bar, VS30 = 0.76 km/s, and κ<sup>0</sup> = 0.025 s.


**Table 5.** Coefficients of CAM for modelling PGV.

This paper focuses on the use of CAM for predictions of PGV. CAM also provides predictions for response spectral accelerations. Refer to "Supplementary Materials", which provides the link to access CAM for response spectrum modelling.

#### *3.2. Translating Seismological Models into GMPEs in Terms of PGV*

Four well-known NGA-East seismological models (as listed in Table 2) have been used as examples to verify the accuracy of CAM. The selected seismological models are namely Atkinson and Boore (1995) [1], Silva, Gregor, and Darragh (2002) [47], Atkinson (2004) [48], and Boatwright and Seekins (2011) [49], with acronyms: AB95, SGD02, A04, and BS11. Each of these models has its own attenuation properties (encompassing geometric attenuation and whole path anelastic attenuation). In this study, the same generic source factor (of the generalised additive double-corner frequency form with stress drop of 200 bar) and site conditions (VS30 = 0.76 km/s and κ<sup>0</sup> = 0.025 s) have been input into the seismological models for defining the Fourier amplitude spectrum (FAS) and for making predictions of PGV through stochastic simulations. All the selected seismological models have been translated into PGV predictive models with a reasonable level of accuracy (as shown in Figure 9).

**Figure 9.** Comparison between predictions by CAM (lines) and simulations of seismological models (symbols). (**a**) AB95 model, with trilinear geometric spreading and Q0 = 680; (**b**) SGD02 model, with magnitude-dependent geometric spreading and Q0 = 351; (**c**) A04 model, with trilinear geometric spreading and Q0 = 893; (**d**) BS11 model, with bilinear geometric spreading and Q0 = 410.

The process of transforming a seismological model into a ground motion prediction model is summarised as follows: (a) stochastic simulations of a seismological model based on a given set of seismological parameters for generating an ensemble of accelerograms, each of which has its own array of random phase angles; (b) calculation of the response spectrum for every accelerogram that has been generated; (c) statistical analysis of the calculated response spectral ordinates for determining their mean values; and (d) developing a GMPE for providing median ground motion predictions by collation of the mean response spectral values (across the natural period range of engineering interests) for different combinations of seismological parameters including M and R.

#### *3.3. Comparing with Historical MMI Data and Existing GMPEs*

In regions of low-to-moderate seismicity where instrumented strong motion data are lacking, ground motion models can be verified by use of macroseismic data. The most commonly used metrics of this type are modified Mercalli intensity (MMI) data. A more recently developed alternative to the MMI scale is the environmental seismic intensity (ESI) scale. Further descriptions of the ESI scale can be found in Section 4. At present, only data presented in the MMI scale are currently available in the two study regions for verifying CAM.

Four candidate GMPEs have been selected for use in the comparative study for evaluating the accuracies of GMPEs in terms of their level of agreement with MMI data. For modelling PGVs in SEA, GMPEs that have been compared are namely: (i) SGC09 (non-cratonic condition) [65], (ii) A12 (shallow earthquake) [18], and (iii) CAM-SEA (this study). For modelling PGVs in SEC, GMPEs that have been compared are namely: (i) CB08 [66], (ii) CY08 [2], and (iii) CAM-SEC (this study). CB08 model and CY08 model are selected because that PGV predictions by these two models have been suggested in the literature to be appropriate for use in South China [61,62]. Details of seismological parameters that have been identified for the two regions for input to CAM can be found in Table 4. Results of evaluations for the two regions are presented in Figures 10 and 11, respectively. In both figures the *x*-axis is the historical recorded MMI values based on observations on average soil sites (as presented in "Supplementary Materials"). The *y*-axis shows the MMI values that are predicted from the four candidate GMPEs and have incorporated an average site factor of 1.5 for transforming predictions from rock sites to average soil sites. The key reference for sourcing intensity information from isoseismal maps was AGSO (1995) [67]. Information presented in website: https://earthquakes.ga.gov.au/ has also been used to identify locations of epicentres of the historical earthquake events.

**(a)** SGC09

**Figure 10.** *Cont.*

**(c)** CAM (this study)

**Figure 10.** Comparison between historical recorded and model-predicted modified Mercalli intensity (MMI) for SEA.

**(a)** CB08

**Figure 11.** *Cont.*

**(c)** CAM (this study)

**Figure 11.** Comparison between historical recorded and model-predicted MMI for SEC.

#### **4. Discussion**

CAM essentially decouples the effect of earthquake source, path, and crust into separate components, thereby enabling it to be regionally adjustable in order that predictions for SEA and SEC can be covered in one model.

In Figure 9, all the simulated results (symbols) and predictions by CAM (lines) are shown to match well for the selected magnitude ranges, indicating that CAM can accurately represent the four selected seismological models: AB95, SGD02, A04, and BS11. The minor mismatch, as displayed in the figures, can be explained by the different relationships that have been used to determine the value of the exponent factor *n* in the Q0 function (*Q* = Q0 *f <sup>n</sup>*), as recommended in the literature; refer to reference [36]. A prominent feature of the A04 model is that the function for determining the *Q* factor is not of the typical form involving parameters Q0 and *n,* but is instead defined as: *Q*(*f*) = max(1000, 893 *f* 0.32). The strategy of setting Q0 = 893 for all Tn < 1.0 s, and Q0 = 1000 for all Tn ≥ 1.0 s has been adopted in predicting response spectral values. According to Bommer and Alarcon [25], the ratio between response spectral acceleration (RSA) and PGV (RSA/PGV) is nearly constant at Tn = 0.5 s, thus, Q0 = 893 was adopted in CAM for predicting the value of PGV in this study.

By a rough glance at Figure 10, predictions by CAM for SEA (Figure 10c) are shown to be in better agreement with the recorded values for MMI than the other two candidate GMPEs (Figure 10a,b). Adopting the geometric attenuation factor of R−1.33, as per recommendations by A12 [18], at short distance as opposed to the conventional factor of R−<sup>1</sup> is controversial, whilst good match between the model predictions and field recorded data has been demonstrated in references [18,19]. It is noted that when calibrating seismological parameters (to achieve agreement between predictions from a seismological model and field recorded data) there are trade-offs of the assumed stress drop values with the assumed rate of geometrical attenuation. Stress drop behaviour of earthquakes, as assumed by different groups of investigators (based on calibration), can be very inconsistent. The geometrical factors adopted by the two study groups can accordingly be very inconsistent too (R−1.33 versus R<sup>−</sup>1), whilst achieving good agreement between predictions and recorded data in their respective studies.

Specific studies on ground motion characteristics of the SEC region have not been reported in the literature. In this study, the geometric attenuation form was adopted from AB95 model to represent the attenuation characteristics in the SEC region. CB08 [66] and CY08 [2] models were selected for comparison purpose, given that PGV predictions by these two models have been suggested in the literature to be appropriate for use in South China [61,62]. There appears to be under-predictions by CAM for SEC (Figure 11c), warranting further investigations whilst overpredictions by the other two models are also shown (Figure 11a,b).

Discrepancies of results presented in this study are considered to have been resulted from the following causes:


Another matter of consideration is the selection of the best metrics for quantifying the intensity of earthquake ground shaking (and ground deformation) in historical, and prehistorical, earthquake

events. The most commonly used macroseismic intensity metrics is the MMI scale that has been used for verifying the accuracy of CAM, as described in the previous section of the article. A review of metrics presented in the MMI scale for quantifying historical earthquake hazards can be found in Ref. [24]. The alternative environmental seismic intensity scale (ESI) [69] is a new macroseismic metric that is based on traces left on the landscape that were caused by earthquake activities in the natural environment. ESI allows the intensity of ground shaking in an earthquake event to be post-dicted in situations where no information on damage to buildings is available or when diagnostic damage-based elements have saturated [70]. Thus, prehistorical events that are not within the scope of any historical archives can be covered. The basic idea of the ESI is to make use of traces of geological and/or geomorphological nature that have been left behind by primary and secondary surface ruptures and mass movements, so generated by large magnitude earthquakes to post-estimate the intensity of the hazard and the magnitude of the event [71]. Example applications of the ESI scale can be found in references [71–75]. Results reported in reference [73] indicate that incorporating ESI data into probabilistic/deterministic seismic hazard analysis can result in significant changes to the modelled PGA values. Another type of macroseismic information is paleo liquefaction data, which can also be employed to provide very valuable information in support of seismic hazard analyses [76].

#### **5. Conclusions**

CAM was introduced as an engineering tool for providing realistic ground motion predictions for intraplate regions. The modelling principle is based on utilising up-to-date (and generic) knowledge on the frequency behaviour of seismic waves that radiated from the source of a local (small and medium magnitude) earthquake in combination with knowledge on local geophysical and crustal conditions, which control frequency modifications along the wave travel path. The seismological model provides the framework for integrating both knowledge bases into the predictions. This decoupling approach to modelling is more scientific, and rational, than simply applying the logic tree procedure to assign weighing factors on existing ground motion models (that are typically biased to conditions and areas where instrumented data are abundant).

This article focused on PGV as a ground motion intensity metric that can be translated into MMI data (which can be compared with data extracted from isoseismal maps of historical earthquake events). The use of algebraic expressions to make predictions of the PGV is an important feature in CAM. Essentially, an adaptable seismological model is presented as a GMPE, thereby waiving away the need for undertaking stochastic simulations along with response spectral computations. An important, and significant, achievement of this article was to have the algebraic expressions in CAM verified. The outcome from the MMI comparative study undertaken for SEA also shows a great deal of promise for CAM, but there are also scopes of improving the match with SEC data in a future study.

The user-friendly setting of CAM (featuring the use of algebraic expressions) serves to facilitate engineering professionals to become more involved with ground motion modelling, thereby gaining a good perspective of the modelling rationale and the underlying assumptions. In summary, this article represents a contribution towards improving transparencies in seismic hazard modelling and in the selection, and scaling of, accelerograms for engineering applications.

The success of CAM in the future relies on the investment of resources into studying crustal and geophysical conditions in intraplate regions for the users of CAM become better informed. This is important given that the quality of the output from any predictive model can only be as good as the quality of input into the model. Further investigations on the MMI-PGV conversion relationship are also warranted for comparison across GMPE models involving the use of MMI data to become more robust.

**Supplementary Materials:** Supplementary materials can be found at http://www.mdpi.com/2076-3263/9/10/422/s1. The detailed information about the proposed geology-based shear wave modelling approach is available online at: https://www.dropbox.com/s/jp4o7j08fe3gqxn/Geology-based-SWV-modelling.pdf?dl=0. The **MMI** recordings for the SEA region can be found online at: https://www.dropbox.com/s/v3eqh3w686ho97j/MMI% 20events%20SEA.xlsx?dl=0. The **MMI** recordings for the SEC region can be found online at: https://www. dropbox.com/s/8vr4wq9mbsee5c0/MMI%20Events%20SEC.xlsx?dl=0. The vs. data can be found online at: https://www.dropbox.com/s/u3v127wjt10cwwk/Shear%20wave%20velocity%20data.xlsx?dl=0. The link to access **CAM** for response spectrum modelling: https://www.dropbox.com/s/jjfbfc8cm2srub3/CAM-Response-spectralacceleration.pdf?dl=0.

**Author Contributions:** Conceptualization, N.L., H.-H.T., and Y.T.; methodology, Y.T., H.-H.T.; software, Y.T.; validation, Y.T. and N.L.; formal analysis, N.L., Y.T.; investigation, Y.T., N.L.; resources, Y.T.; data curation, N.L., H.-H.T., and Y.T.; writing—original draft preparation, Y.T.; writing—review and editing, N.L., H.-H.T., and E.L.; visualization, Y.T.; supervision, N.L., H.-H.T., and E.L.; project administration, N.L.; funding acquisition, N.L.

**Funding:** The first author was funded by China Scholarship Council, CSC grant number 201506030102. The project was funded by Commonwealth of Australia through the Cooperative Research Centre program.

**Acknowledgments:** The extensive technical support given by Yiwei Hu at the University of Melbourne is gratefully acknowledged.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).
