A Refined Vs 30 Map for Taiwan Based on Ground Motion Attenuation Relationships

Seismic hazard evaluations require an estimate of the expected ground motion at the site of interest usually by using attenuation relationships. The mean shear-wave velocity over the top 30 m (Vs30) is incorporated in the ground motion attenuation relationships in this study. By comparing the standard deviations of the residuals between the observed and predicted values before and after incorporating the site effect term Vs30, the reduction in standard deviation for the peak ground velocity (PGV) is significantly reduced by about 11%. Clearly, the refined attenuation relationships will be more useful for engineering purposes. Analyzing the site effect term using the amplification factor (relative to a site with Vs30 = 760 m s-1), has revealed that the Changhua Plain, Chianan Plain, Pingtung Valley, Ilan Plain, and Taipei Basin have high values, implying large ground motion amplification. Following a disastrous earthquake, quick assessment and timely peak ground acceleration (PGA) and PGV map reporting will be critical for effective emergency response operations. After an earthquake we can combine the simple attenuation relationships, as determined from Model 1, to provide near real-time estimation and reporting of the PGA and PGV values for the Taiwan area. We can also use the relations between the intra-event site residual and the Vs30 to estimate the Vs30 for stations that have recorded strong motions, but do not yet have Vs30 information. Our approach including sites with estimated Vs30 has resulted in a refined Vs30 contour map that can be used for more realistic seismic hazard assessment for Taiwan. This approach is especially applicable to the foothill and mountain areas.


INTRODUCTION
Local site effects play an important role in modifying seismic motions.This fact is evidenced by the observed damage distribution patterns in numerous earthquakes and is reflected in some seismic codes and ground-motion prediction equations (GMPEs), which prescribe stronger motions for sites located in soft sediments compared to bedrock sites (Borcherdt 1970).Local site conditions at an accelerograph station can dramatically affect the strong motions recorded (Douglas 2003).In particular, soft sediments with low seismic velocities relative to the underlying bedrock can lead to large ground motion amplification at the surface (Pratt et al. 2003).Some studies have found that shallow soil sites have significantly higher ground motions than rock or stiff soil sites (Campbell 1981).Site-effects estimation has been a major issue in engineering seismology for the last 20 years.This is due to the fact that seismic hazard is strongly influenced by site effects because site conditions strongly affect the frequency content and also the amplitude of ground motions.This is also because most populated areas are located in sedimentary basins (Drouet et al. 2008).
The average shear-wave velocity of the upper 30 m of a soil profile (Vs30) is used as a quantitative parameter for most earthquake ground-motion site-effect studies (Boore et al. 1993;Anderson et al. 1996;Bolt and Abrahamson 2003;Huang et al. 2007Huang et al. , 2009; Lee and Tsai 2008).The Vs30 is adopted as an international standard for soil classification since it was proposed by the United States National Earthquake Hazard Reduction Program (NEHRP) (Gallipoli and Mucciarelli 2009).Starting from the work of Borcherdt (1994), the Vs30 has become a standard parameter for quantifying response at a site.Accordingly, successive NEHRP regulations (NEHRP 1994(NEHRP , 1997(NEHRP , 2001) ) included Vs30 to define their respective local site categories (Cadet et al. 2010).For the development of ground-motion models, developers systematically evaluated a list of predictor parameters to consider for predicting earthquake shaking intensity.The most significant decision made by all developers was to use Vs30 as a parameter for characterizing soil-stiffness effects on ground motions (Power et al. 2008).For example, the "Next Generation Attenuation of Ground Motions Project (NGA)" of the Pacific Earthquake Engineering Research Center (PEER) has directly used Vs30 for the ground-motion attenuation models (Power et al. 2008).In Taiwan, Huang et al. (2007Huang et al. ( , 2009) ) used the Vs30 data in Central Taiwan and in the Taipei Basin to estimate the high frequency site amplification using the quarter-wavelength method.Lee et al. (2012) developed a new empirical Arias intensity attenuation relationship for shallow crustal earthquakes in Taiwan considering Vs30.Their results show that the incorporation of Vs30 can significantly reduce regression error.
Several studies have been made to evaluate the site effects based on accelerographic recordings.For example, examination of the residuals for sites with different soil categories was shown to be a useful method using sets of records where site information was not complete (Abrahamson and Litehiser 1989).Liu and Tsai (2005) analyzed the residuals of peak ground acceleration (PGA) and peak ground velocity (PGV) with respect to site conditions.They found the contour maps of site residuals for the PGA and PGV data, especially for the PGV, were highly correlated with the regional geology and topography of Taiwan.They further pointed out that almost all major metropolitan areas coincided with high residual areas that would require special attention in structural seismic design.Choi and Stewart (2005) evaluated the ground motion amplification factors from residuals between accelerations from recordings and attenuation relationships for active seismic regions.They also found site amplification to increase with decreasing Vs30.
In this study we first develop refined attenuation relationships by incorporating a site effect term related to Vs30, the average shear-wave velocity in the upper 30 m of sediments, aiming to reduce the standard deviation of the predicted ground motion in large earthquakes.We will also investigate the variation in residuals with Vs30 in order to construct a refined Vs30 map for the entire Taiwan area.The results of this study will provide valuable information for site evaluation of critical facilities in high earthquake hazard regions, as well as for land-use planning.

DATA
Beginning in 1991 the Central Weather Bureau Seismology Center (CWBSC) embarked on a six-year seismic strong-motion instrumentation program, known as the Taiwan Strong Motion Instrumentation Program (TSMIP) (Liu et al. 1999).The TSMIP accelerographic network dataset used in this study consists of 617 free-field stations, primarily covering densely populated areas.The locations of these stations are given in Table 1 and shown in Fig. 1.These free-field stations are densely spaced approximately 5 km apart, and only about 3 km apart in urban areas.Concurrently, the Institute of Earth Science, Academia Sinica installed another strong motion network with similar instrumentation in the Central Mountain Strong Motion Array (CMSMA) to provide strong motion data for the mountainous areas (http://www.earth.sinica.edu.tw/~smdmc/cma/cma.htm).The locations of these 10 stations are given in Table 2 and shown in Fig. 1.Strong seismic ground motion data obtained from a total of 627 stations from the TSMIP and CMSMA networks are used in this study to derive new ground motion attenuation relationships.
The Central Weather Bureau (CWB) and National Center for Research on Earthquake Engineering (NCREE) have undertaken a free-field strong-motion station drilling project to construct an engineering geological database for TSMIP (EGDT) since 2000.A total of 468 free-field strong motion stations were surveyed by 2010, with 439 of these stations drilled and their P-and S-wave velocities measured using a suspension PS-logging system (Kuo et al. 2011(Kuo et al. , 2012)).Measured Vs30 records are available at 426 of the 617 TSMIP free-field strong-motion stations.Records from these stations are used to develop the attenuation relationships which incorporate a site effect term, aimed to reduce the standard deviation of the predicted ground motion for engineering applications.The crosses in Fig. 1 represent the TSMIP stations with available Vs30 records.Relevant localities are also indicated.
Over 7900 digital accelerograms were recorded from 51 crustal earthquakes by the TSMIP and CMSMA networks, as shown in Fig. 1, and selected for this study (Liu and Tsai 2005).The M w magnitudes range from 4.0 -7.1.These records are analyzed to investigate the attenuation relationship dependence on magnitude, regional earthquake clusters and site effects.Figure 2 shows the data distribution of these events in terms of magnitude and source-to-site distance.

METHODOLOGY
The strong-motion attenuation relationships express earthquake ground motion parameters as functions of simple parameters characterizing the earthquake source, the propagation path between the earthquake source and the site, and the geological conditions beneath the site.The following equation form is used in this study (Liu and Tsai 2005) where PGA and PGV are ground motion parameters, X is the closest distance to the rupture surface or hypo central distance.We characterize the source-to-site distance in terms of the closest distance to the rupture surface, r rup .If the rupture surface is not defined for an event, the hypo central distance is then used as the source-to-site distance, M w is moment  †Determined by Lee et al. (2001), on the basis of surface geology and borehole data.‡Determined according to UBC1997.†Determined by Lee et al. (2001), on the basis of surface geology and borehole data.‡Determined according to UBC1997.
magnitude, a is the geometric spreading coefficient, b is the inelastic attenuation coefficient, c is magnitude coefficient, d is a constant, h 1 and h 2 are close-in distance saturation coefficients.Vref is the reference shear-wave velocity.The coefficients a, b, c, d, dl, e, h 1 , h 2 , and Vref are to be determined by regression from the data.In Eq. ( 2) the constant changed from d to dl while we specify Vref, the reference velocity to be 760 m s -1 .In this case the soil amplifications are specified relative to motions that would be recorded on a NEHRP B/C boundary site condition.
The coefficients in the equation for predicting ground motion were determined using a two-stage regression procedure.The least square method was used in the regression.A similar approach was used previously by Joyner and Boore (1993) and Liu and Tsai (2005).
The two models are adopted to account for different situations.Model 1, uses recordings from all 627 strongmotion sites, is especially useful for early warning systems to make quick assessment and timely reporting of the PGA and PGV maps.Their results will be critical for effective emergency response operations.Model 2, using 426 strongmotion sites to incorporate a site effect term with available Vs30, is aimed to reduce the standard deviation of the predicted ground motion.This approach emphasizes direct use of strong ground motion recordings for seismic ground shaking estimation for engineering applications.
The residual value, i.e., site response factor, is defined as the difference between logarithms of the observed and the predicted ground motion, and is expressed by the following equation: where Y s is the observed value, Y r is the predicted value from Eq. ( 1), h is the earthquake inter-event errors with standard deviation equal to x, and f is intra-event errors with standard deviation equal to v. The h and f are assumed to be normally distributed independent variants with variances 2 x and 2 v .The amplification factors of site effect can be calculated from exponent of (r).The standard deviation of total residual T v is given by the equation: The residuals due to regression were decomposed into interevent (earthquake-to-earthquake) and intra-event residuals.
We analyzed the relations between the intra-event residuals f and the average shear-wave velocity in the upper 30 m of sediments, Vs30 using the following equation (Boore et al. 1997;Choi and Stewart 2005;Liu et al. 2013): where the coefficient f and Vref are to be determined by regression from the data.

Attenuation Relationships for Vertical and Horizontal PGA and PGV
Ground motion characteristic studies in Taiwan require ground-motion attenuation models.Attenuation relationships, or "GMPEs", provide an efficient means for predicting the level of ground shaking and its associated uncertainty at any given site or location, as well as for use in seismic hazard analyses (Bolt and Abrahamson 2003).An attenuation relationship is a mathematical equation that relates a specific strong-motion ground shaking parameter to a number of earthquake seismological parameters and the recording site.The seismological parameters quantitatively characterize the earthquake source, the wave propagation path between the source and the site, and the soil and geological profile beneath the site (Campbell 2003).
Regressions on the dataset for Model 1 without differentiating site conditions and Model 2 with site conditions have resulted in the coefficients of the attenuation relationships, as given in Tables 3 and 4, respectively, for the vertical and horizontal components of PGA and PGV in Taiwan area.In Tables 3 and 4, σ 1 and σ 2 are standard deviations on ln(PGA, PGV).

Analyses of Site Total Residuals for PGA and PGV
Residual examination for sites with qualitative soil categories is a useful method for sets of records where site information is not complete, and hence cannot be included explicitly within the equation (Abrahamson and Litehiser 1989).In this study we analyze the residuals to investigate variations in PGA and PGV with respect to site conditions.
The mean and standard deviation of total and intra-event residuals for 627 stations for the vertical and horizontal PGA and PGV are given in columns 5 -8 of Tables 1 and 2. The corresponding contour maps of mean of total residuals for the horizontal PGA and PGV are shown in the Figs. 3 and 4, respectively.The total residual patterns, especially those for PGV, agree reasonably well with the regional geology and topography patterns.Notably, the PGV residual is more  sensitive to the site class than the PGA residual, because the PGA is a high-frequency parameter which is less affected by local site conditions.Local site conditions can dramatically affect the strong motions recorded (Douglas 2003).Furthermore, in order to understand the relation between the ground motion parameter site residuals and Vs30, we plot the Vs30 contour map based on 426 measured data as shown in Fig. 5.We also found that the PGV residual contour patterns in major plain areas are similar to the Vs30 contour patterns.The purpose of Model 1 is early warning systems to make quick assessment and timely reporting of the PGA and PGV maps.Hence, the equation form is so simple that is does not include the fault type (focal mechanism), intra-or inter-event term and site effect term such as Vs30.Their results will be critical for effective emergency response operations.Accordingly, following a disastrous earthquake, quick assessment and timely reporting of PGA and PGV maps will be critical for effective emergency response operations.Thus, after an earthquake we can combine the simple attenuation relationships, as determined from Eq. ( 1), and the total residuals, as determined from Table 1, to provide near realtime estimation and reporting of the PGA and PGV values for the Taiwan area.

Comparisons of Model Predictions Including Vs30
In early times most attenuation relationships used broad site categories such as "rock", "stiff-soil", and "softsoil".There has recently been a move toward using quantitative site classifications based on the shear-wave velocity measured at the strong-motion site.The most commonly used parameter is the average shear-wave velocity over the top 30 m (Vs30) (Boore et al. 1997;Bolt and Abrahamson 2003).For ground-motion model development, developers systematically evaluated a list of predictor parameters to consider for predicting earthquake shakings.One of the most significant decisions made by all developers was to use the average shear-wave velocity in the upper 30 m of sediments, Vs30, as the parameter for characterizing soilstiffness effects on ground motions (Power et al. 2008).
A major limitation in using quantitative site descriptions is that the Vs30 information is not available for most strong-motion recording sites.This situation has been greatly improved for our study because measured shear-wave velocity profiles are available at 439 strong-motion recording sites in the Taiwan area.Of these, 426 stations were used in this study as given in Table 1.To incorporate a site effect term based on the average shear-wave velocity over the top 30 m (Vs30) into the attenuation relationship, the analytical form used in this study is given in Eq. ( 2).
As our understanding and modeling of attenuation relationships improve, there will be a trend toward reducing the modeling variability.In empirical attenuation models the modeling variability given for the model is the standard deviation (Bolt and Abrahamson 2003).In addition to the median ground motion, the standard deviation of the ground motion is also important for seismic hazard analyses.After incorporating a site effect term, Vs30 in the attenuation relationships, the standard deviations between the observed and predicted values are reduced from 0.683 -0.647 for horizontal PGA, and from 0.663 -0.587 for horizontal PGV, respectively.We found only minor reduction in standard deviation for PGA.In contrast, the PGV standard deviation is significantly reduced by about 11%.
We analyzed the site effect term using the amplification factor (relative to a site with Vs30 = 760 m s -1 ): exp[-e ln(760/Vref) + Res_av], where Res_av is the average intra-event residual from Eq. ( 1) for each site and the e and Vref of horizontal PGA and PGV are given in Table 4.The corrected site amplification factor contour maps relative to a NEHRP B/C boundary site condition (Vs30 = 760 m s -1 ) are plotted in Fig. 6 for horizontal PGA and in Fig. 7 for horizontal PGV, respectively.The amplification factors are contoured at 0.2 intervals.They range from 3.4 -0.4 for horizontal PGA and from 4.0 -0.4 for horizontal PGV, respectively.Both total residual and amplification factor contour maps of horizontal PGV have similar patterns, revealing that the Changhua Plain, Chianan Plain, Pingtung Valley, Ilan Plain, and Taipei Basin have high values, implying large amplification of ground motions.

Construction of a Refined Vs30 Map
Correlations between the intra-event residual and  5. Reading from Fig. 11 the value of R was 0.712 (R 2 = 0.5070) indicating that the horizontal PGV residuals are much more dependent on Vs30 than the residuals for the other three ground motion parameters, probably due to their frequency content which is less affected by local site conditions.
Since the Vs30 information is not available for a significant number of TSMIP and CMSMA stations, especially for those having a NEHRP B and C site conditions, here we use the relations between the intra-event site residual and the Vs30 for horizontal PGV, as given by Eq. ( 6), to estimate the Vs30 for these stations that have recorded strong motions.The estimated Vs30 values, denoted as Vs30N, are given in the last column of Tables 1 and 2 We compared our refined Vs30 results with studies by the United States Geological Survey (USGS) (http://earthquake.usgs.gov/hazards/apps/vs30/predefined.php) and Lee and Tsai (2008).The refined Vs30 pattern in this study is similar to that in USGS.Previously, Lee and Tsai (2008) mapped the Vs30 distribution in Taiwan using 230 P and S wave velocity (PS) logging at soil and soft rock strong-motion station sites and 4885 engineering boreholes.The Vs30 map supplied important knowledge for each strong-motion station, and for sites between stations in Taiwan.However, the accuracy of such mapping is inevitably dependent on the amount and quality of data.In eastern Taiwan and the southern tip of Taiwan, where boreholes are few and scattered, the map accuracy is relatively poor.The mapping results use the 2000 -2005 Vs measurements and may be further refined based on new data.Kuo et al. (2012) reclassified the 439 drilled free-field TSMIP stations and compared the results with that of Lee and Tsai (2008).It was found that 80 of the 436 stations in common were misclassified by Lee and Tsai (2008).It is probably that Lee and Tsai (2008) focused on mapping the soil sites, and thus tentatively assigned all of the unmeasured rock site stations a value of 760 m s -1 .
We compared the Vs30 estimations from Lee and Tsai (2008) with the Vs30N from this study.The above relation is plotted in Fig. 13.The regression equations and the coefficient of determination R 2 are also shown in the figures.R is

CONCLUSIONS
From above results and discussion, we can summarize our findings as follows: (1) The PGV residual contour pattern is highly similar to the Vs30 contour pattern.The horizontal PGV residual as found more dependent on Vs30 than the vertical and horizontal PGA, as well as vertical PGV, probably due to their high-frequency content which is less affected by local site conditions.(2) By comparing the standard deviations between the observed and predicted ground motion values before and after incorporating the site effect term Vs30, the reduction in standard deviation for PGA is only moderate.In contrast, the PGV standard deviation is significantly reduced by about 11%.Evidently, the refined attenuation relationships will be more appropriate for engineering applications.
(3) After analyzing the local site effect in terms of the amplification factor (relative to a site with Vs30 = 760 m s -1 ), the Changhua Plain, Chianan Plain, Pingtung Valley, Ilan Plain, and Taipei Basin were revealed to high values, implying large ground motion amplification.Large parts of the Central Mountain Range have low values, implying potential ground motion de-amplification.(4) Following a disastrous earthquake, quick assessment and timely reporting of PGA and PGV maps will be critical for effective emergency response operations.Thus, after an earthquake, we can combine the simple attenuation relationships, as determined from Eq. ( 1), and the site residuals, as given in Tables 1 and 2, to provide near real-time estimation and reporting of the PGA and PGV values for the Taiwan area.(5) Finally, using the correlations between the intra-event residual and Vs30 according to Eq. ( 6), we can estimate Vs30 for stations that recorded strong motions, but whose Vs30 information is not available.Our approach including sites with estimated Vs30 has resulted in a detailed Vs30 contour map to facilitate more realistic seismic hazard assessment for Taiwan, especially applicable to the foothill and mountain areas.
TotalRes = Total Residual in ln unit.(5) IntraRes= Intra-event Residual in ln unit.(6) PGAh = The mean of residual for horizontal component of peak ground acceleration.(7) PGVh = The mean of residual for horizontal component of peak ground velocity.(8) sd = standard deviation.(9) R = Number of Records.(10) Vs30 = The average shear-wave velocity in the upper 30 m of sediments in m s -1 .(11) Vs30N = The estimated Vs30 values from the relation between the intra-event residual and the Vs30 in m s -1 .

Fig. 1 .
Fig. 1.Distribution of the earthquakes, Taiwan Strong Motion Instrumentation Program (TSMIP) and Central Mountain Strong Motion Array (CMSMA) free-field stations used in this study.Relevant localities and topography are also indicated.

Fig. 2 .
Fig. 2. Magnitude-distance distribution of the earthquakes used in this study.

Fig. 3 .
Fig.3.The total residual contour map (in ln unit) for horizontal PGA.
. A refined Vs30 contour map that includes the 201 Vs30N sites and 426 Vs30 sites is shown in Fig. 12.By comparing the Vs30 contour maps, as shown in Fig. 5 based on 426 Vs30 sites with Fig. 12, we can find the following results: (1) Fig. 12 offers a more detailed Vs30 contour map for using data from more Vs30 sites.(2) The class D contour pattern is similar.(3) The class E area appears in the Changhua Plain in Fig. 12 needed to pay more attention.(4) The rocky areas of classes B and A coincide with the foothill and mountain areas distribution in Fig. 12.

Fig. 6 .
Fig. 6.The corrected site amplification factor contour map relative to a NEHRP B/C boundary site condition (Vs30 = 760 m s -1 ) for horizontal PGA.

Fig. 7 .
Fig. 7.The corrected site amplification factor contour map relative to a NEHRP B/C boundary site condition (Vs30 = 760 m s -1 ) for horizontal PGV.

Fig. 8 .
Fig. 8.The relation between the intra-event residual of vertical PGA and the average shear-wave velocity in the upper 30 m of sediments Vs30.

Fig. 9 .
Fig. 9.The relation between the intra-event residual of horizontal PGA and the average shear-wave velocity in the upper 30 m of sediments Vs30.

Fig. 10 .
Fig. 10.The relations between the intra-event residual of vertical PGV and the average shear-wave velocity in the upper 30 m of sediments Vs30.

Fig. 11 .
Fig. 11.The relations between the intra-event residual of horizontal PGV and the average shear-wave velocity in the upper 30 m of sediments Vs30.

Fig. 12 .
Fig. 12.The refined Vs30 contour map based on combined measured and estimated data.

Fig. 13 .
Fig. 13.The relation between the Vs30 estimated from Lee and Tsai (2008) and the Vs30N from this study.

Table 1 .
Station code, location, site classifications, ground motion total and intra-event residuals, measured and estimated Vs30 of the TSMIP stations.

Table 3 .
Coefficients for the vertical and horizontal components of PGA and PGV from Eq. (1).

Table 4 .
Coefficients for the vertical and horizontal components of PGA and PGV from Eqs. (2) and (3).

Table 5 .
Coefficients for the vertical and horizontal components of PGA and PGV from Eq. (6).