Analysis of epistemic uncertainty associated with GMPEs and their weight within the logic tree for PSHA : Application to Taiwan

In probabilistic seismic hazard analysis (PSHA), the standard practice is to select a set of appropriate ground motion prediction equations (GMPEs) and assign weights on the logic tree, especially for regions where strong motion data are sparse and where no indigenous GMPE exists. Subjectively assigning weights to a set of models usually has the disadvantage of not obtaining mutually exclusive and collectively exhaustive models because of sparse or unavailable data. Therefore, the development of logic tree weightings in PSHA remains a major challenge. In this study, a distance metric measure for GMPEs is first analyzed to show how a set of GMPE’s prediction models can be partially reconciled by using high-dimensional information visualization techniques. Visualization of a large suite of GMPEs onto a 2-D graphical map provides a powerful theoretical framework that can guide the selection of a set of representative models. These models are considered mutually exclusive and collectively exhaustive, and have the ability to represent the center, body, and range of ground-motion distribution in a logic tree analysis. Second, determination a set of weights for PSHA are estimated based on the residuals, likelihood and EDR-index which are known as data – driven weight. The other weight type is non-data driven which is calculated based on the probability of a random model generated on a Sammon’s map. The methods presented here, that improve consistency in the weight assignment, can help to reduce overall epistemic uncertainties and offer a way of assigning weight on the logic tree. Finally, an example for assigning the weights on the logic tree for PHSA are discussed. Article history: Received 8 May 2018 Revised 11 August 2018 Accepted 13 August 2018


InTroducTIon
The island of Taiwan is located on the complex boundary between the Eurasian Plate and the Philippine Sea Plate (Wu et al. 2013).Figure 1a shows the tectonic setting of Taiwan.The Philippine Sea Plate is moving towards the northwest and subducting beneath the Eurasian Plate along the Ryukyu trench.South of Taiwan the Eurasian Plate is being subducted beneath the Philippine Sea Plate along the Manilla trench.The movement was studied at a rate of 80 mm yr -1 of the Philippine plate (Seno et al. 1993;Yu et al. 1997), and thus has been caused high seismicity on the island (Tsai 1986).According to the Central Weather Bureau (CWB), more than 18000 seismic events annually strike the island with approximately 1000 felt earthquakes.On average, 50 earthquakes occur per day; most of which are shallow crustal earthquakes with depth range from 0 to 30 km.More than 425 destructive earthquakes on the island have recorded in the seismic areas of Taiwan since 1983 as recorded in the database.These include shallow crustal earthquakes and subduction zone earthquakes.Considering the seismic activity in this area, probabilistic seismic hazard analysis becomes a very important issue in most of the engineering constructions.
Probabilistic seismic hazard analysis (PSHA) (Cornell 1968;McGuire 2008) is widely used for evaluating the seismic risk at a specific site of an engineering project, and for estimating seismic loads.This is typically done using ground motion prediction equations (GMPEs), which can generally be expressed as mathematical functions relating the peak ground motion (PGA), or the response spectral acceleration, to earthquake related parameters such as magnitude (M), distance (R), etc.For a given site, one can select GMPEs published in the literature that satisfy the minimum criteria proposed by Cotton et al. (2006) and later revised by Bommer et al. (2010) for use in PSHA.In general, the predictions of GMPEs will differ from each other and these differences can be interpreted as epistemic uncertainty in median prediction, which is a major contributor to the overall epistemic uncertainty of the seismic hazard at a low exceedance rate (Toro 2006).
In this study, for the analysis of epistemic uncertainty associated with GMPEs in Taiwan area, we consider only shallow earthquakes relevant to the development of GMPE by selecting the events with magnitude range from 3.5 to 7.65 and with the depth range from 0 to 50 km.The poorly sampled data was not selected; as a result, a total number of 13423 recordings obtained from 158 earthquakes can be considered depending on the model developers.78 of which are reverse faulting, 23 of which are normal faulting, and 57 of which are strike-slip faulting.Figure 1b shows the magnitude-distance distribution of data that was used for the adjustment/development process for general case.The adjustment/development process was performed for each of the original GMPE using the same database, but the number of records used to adjust each model is substantially different from each other.
The adequate quantification of epistemic uncertainties is an important task in PSHA concerning how to obtain the center, body, and range of the technically defensible inter-pretations (the CBR of the TDI) of both the seismic source model and the GMPEs.This can usually be done by employing a logic-tree framework in which different weights will be given to each branch of the logic tree where the candidate GMPEs exist.There has been considerable progress recently regarding various aspects of logic tree analysis and GMPEs for PSHA.For example, the likelihood (LH)-based method, the average sample log-likelihood (LLH) proposed by Scherbaum et al. (2004Scherbaum et al. ( , 2009)), and recently the Euclidean distance-based ranking (EDR) method suggested by Kale and Akkar (2013) have become useful tools to help experts to judge the applicability of different GMPEs to a specific region when the observed ground motion data are known.Atkinson et al. (2014) introduced the so-called "scale-backbone approach" in which a single GMPE is selected and scaled up/down by changing its magnitude and distance to cover the range of epistemic uncertainty.This approach provides a method to develop a set of GMPEs that are mutually exclusive and collectively exhaustive in terms of the ground motion amplitude of a given magnitude and distance, but that do not capture the range of different scaling with magnitude and distance (GeoPentech 2015).
Whether hazard analysts choose to use different GMPEs or a scaled back-bone approach, they are still tasked with assessing the center, body, and range of epistemic uncertainty.Scherbaum et al. (2010) proposed using high-dimensional visualization techniques to provide a graphical representation of different GMPEs.They used Sammon's map (Sammon 1969) as a tool to project median GMPE prediction onto a two-dimensional map, where each GMPE can be visualized as a single point.A set of GMPEs on the map can be thought of as a two-dimensional projection of the GMPE model space.The Sammon's map approach can also help to develop a large suite of new models from existing GMPEs, and thus lead to a continuous distribution of the median GMPEs that comprise the space where the assessment of center, body, and range should be based.
In this paper, after reviewing the use of mapping GMPEs using the work of Scherbaum et al. (2010), GeoPentech (2015), Kuehn et al. (2015) for PSHA, the ground motion space of the CBR of the TDI of the GMPEs is assessed through a set of candidate GMPEs and applied to Taiwan.A metric to measure the distances between GMPEs is described and how these distances can be interpreted in a manageable way will be presented.Finally, logic tree weights for median GMPEs based on either the available data or the common form models are developed.The weights are carefully examined and judged so that they are adequately explainable with regards to the center, body, and range of the technically defensible interpretation of the GMPE prediction.The weighting scheme provides guidance and offers a method of assigning weights to the logic tree, and the average weight is suggested for use in PSHA.

orIGInAl GMPES And cAndIdATE GMPE ModElS
The most important works in PSHA are the selection of a number of GMPEs, their adjustment to a specific target region and/or site, and potentially testing them against existing observations.In this study, a set of GMPE models is selected under the framework of Taiwan SSHAC level 3 projects.Then, these GMPE models were adjusted by modifying their key modeling parameters and calibrating them to match the empirical data in Taiwan.These adjusted/developed GMPEs are listed as follows: (1) Adjusted GMPE of Abrahamson et al. (2014), referred to as ASK14adj.
The base-case function forms of these GMPEs can be found in Appendix A. Among these seven adjusted GMPEs, four of the NGA-West2 models (ASK14, BSSA14, CB14, and CY14), originally known as the global models were developed for multiple regions, which are more complicated than the other three models (Id14, ASB14, and Bi14).The adjusted GMPEs were performed using the events with Mw less than 6.5.The minimum number of model's parameters were sought to analyze against Taiwan data, while other model's parameters associated with the terms controlling, e.g., large magnitude scaling and non-linear site effect were held fixed to the values determined from the original GMPEs.
Table 1 summarizes the important features of the adjusted coefficients used in the adjusted GMPEs.Chao17 and Phung17 models are the newly developed using the same database with the adjusted models.However, their regression approaches were different from each other.Chao17 and Phung17 were regressed for many more model's parameters as compared to the adjusted models.Moreover, Phung17 was developed by using a data set mainly from Taiwan earthquake events and partially from foreign events (for M > 6.5), while Chao 17 used only Taiwan ground motion data.The major differences in data subset used for the development of candidate GMPEs can be briefly summarized below.
For adjusted GMPEs, Taiwan data was used with following criteria: (1) Consider main shock events only with at least 15 number of recordings per event.

GMPE constant Style of faulting Magnitude
Table 1.Adjusted coefficients used in the candidate GMPEs.
(3) Consider all data in the distance range 0 -R max in which R max is maximum usable distance available in the database.
For Phung17 model, a combined NGA-West2 and Taiwan data set was used with following criteria: (1) Consider all crustal events including main shock and aftershock event from Taiwan.(2) Select magnitude events with M range from 6.5 to 7.9 and R rup range from 0 to 150 km from NGA-West2 database.
(3) Select the events with at least 5 recordings.
(4) Consider all data in the distance range 0 -R max in which R max is maximum usable distance available in the databse.
For Chao 17 model, Taiwan data was used with following criteria: (1) Consider all crustal events including main shock and aftershock event from Taiwan.
(2) Consider all data in the distance range 0 -R max in which R max is maximum usable distance developed by the developer.
The adjusted/developed GMPEs from now onwards are called as the candidate GMPEs and are used for the case of study.Comparisons between the original GMPEs and the adjusted/developed models are shown in Fig. 2. In Fig. 2, the plots with a 5% damped response spectrum for the magnitude scaling (R rup = 30 km), the distance scaling (M = 7), the period of T = 0.01 sec (or PGA), and for vertical strikeslip earthquakes and engineer rock site conditions (V S30 = 760 m s -1 ) are shown.

SAMMon'S MAPPInG for PSHA
Sammon's mapping (Sammon 1969) is a non-linear dimensional reduction algorithm to detect data structure and is a tool to project a vector of median ground motion prediction from a high-dimensional ground motion space onto a two-dimensional ground motion space.The median GMPE prediction for each model on a high-dimensional ground motion space is represented by a vector whose component is the ground motion intensity (e.g., PGA) estimated with many different predictor variables such as magnitude and distance.The median GMPE prediction for two models on two-dimensional ground motion space is represented by two points, where the relative distance between them is nearly the same as those distances in the high-dimension ground motion space.The Sammon's mapping configuration can be obtained by minimizing the misfit function (E) that measures the error between the Euclidian distance of the GMPE on the high-dimensional space and the Euclidian distance of the GMPE on the two-dimensional space.The iterative process is continued until the objective function reached 10e -5 or less so as to provide a sufficiently accurate projection that maintains GMPE structure from the high-dimension to the 2-D map.The misfit function is described as follows in Eq. ( 1): where ij map D corresponds to the GMPE-distance between GMPEs i and j on the Sammon's mapping configuration with two dimension, and ij GMPE D is the GMPE-distance between GMPE model-i and model-j in high-dimension ground motion space which corresponds to the Euclidian distance (the L2 norm) and can be expressed below in Eq. ( 2) where w k denotes the weight given to the k th scenario.If the GMPE-distances treat each magnitude/distance scenario equally with uniform weight (w k = 1), the unweighted GMPEs map is generated even though some of them may not contribute to the hazard.For PSHA, we are interested in hazard space, and thus a map resulted by weighting GMPEdistances is considered for the contribution of magnitude/ distance scenarios to the hazard.For this, w k is obtained from the de-aggregation matrix for the hazard site, and the GMPE-distances are now weighted by their contribution to the hazard.The weight, w k according to GeoPentech (2015) is defined in Eq. (3).
. , w DEAGG M R NS 0 5 1 where NS is the total number of scenarios, M k and R k are the magnitude and distance of the k th scenario, and ^h is the hazard de-aggregation matrix for the bin containing the M k and R k .Figure 3 shows an average hazard de-aggregation for T = 0.01 sec and 10000 return period, and this is resulted from average of seven adjusted models and two newly developed models with SSC model utilizing in Taiwan for the hazard sites.
For the purpose of comparing the differences among GMPEs, each GMPE is evaluated at the hazard scenarios (magnitude-distance pairs) that is shown in Fig. 3 for T = 0.01 sec.The GMPEs are evaluated at moment magnitude M = 5.1, 5.3, …, 8.3; rupture distance R rup =1, 3, 5, 7, 9, …, 294 km; and depths to the top of the rupture Z tor = 1, 3, 5, 7, 9, …, 32.5 km.The depth to the top of the rupture Z tor is set to the limit of 32.5 km because the shallow crustal events that have occurred in Taiwan are mostly within a depth of 35 km.Hence, each GMPE is evaluated at 16 magnitudes, 19 distances, and 11 depths, leading to a total of 3344 scenarios.Consequently, one can represent each GMPE as a point in a 3344-dimensional ground motion space, where each coordinate corresponds to the prediction for one particular magnitude-distance pair.For the R jb -based GMPEs, the predictor variable R jb can be easily calculated from R rup and Z tor since the dip angle is known, while sites with negative R jb are not included in the calculation.The reverse, normal and vertical strike -slip faulting styles are randomly assigned for 3344 sites, and this is approximately equivalent to one third of the computed sites for each style of faulting.All sites are located on the footwall side with an assumed dip angle of 45° for the reverse and normal fault types.The site condition is set to V S30 = 760 m s -1 .All other predictor variables, such as the fault rupture width (W) and the depth to a shear wave velocity horizon of 1000 m s -1 (Z1.0), are set to the default values for the GMPEs as required.
In Fig. 4, the resulting Sammon's map is shown for comparison between the original GMPEs (seven models) and the candidate GMPEs (nine models) using unweighted GMPE-distances and weighted GMPE-distances with the same scenarios.It can be seen that when using the weighted GMPE-distances, the relative distances between the GMPEs become smaller as compared to that used unweighted GM-PE-distance.This is due to significant contribution of large magnitude and short distance de-aggregation scenarios to the hazard, which can already be seen in Fig. 2. In other words, GMPEs are now compared using more focused range -the range that higher de-aggregation weights are given to large magnitude-short distance scenarios but not given to small magnitude -long distance scenarios.In overall, the changes of the candidate models are relatively minor with relative configuration of the models similar for both maps.However, the changes of the original models are seen to be more significant because one possible reason is that the de-aggregation weights were computed from the candidate models but not from the original models.Another observation is that the directions in which the scaling of the GMPEs changes are not the same between two maps any more.In particular, the average model over the whole magnitude and distance range is not necessarily the average model over the weighted range.Nevertheless, the map calculated based on the weighted GMPE-distances provides valuable information about the relative differences between the GMPEs over hazard relevant magnitude/distance ranges.
It is believed that the map calculated based on the weighted GMPE-distances provides valuable information about the relative differences between the GMPEs over relevant hazard magnitude/distance ranges.From the description of Sammon's mapping, the relative positions (or distance) between two GMPEs on the map does have a meaning which corresponds to the relative positions in high-dimensions.Therefore, they can be rotated and mirrored without changing the ij map D .

AnAlySIS of EPISTEMIc uncErTAInTy THrouGH THE InTErPrETATIon of GMPES on SAMMon'S MAP
In this section, Sammon's map is generated for the original candidate GMPEs (seven models) and the candidate GMPEs (nine models) using the weighted GMPE-distance.The weighted map, which considers magnitude/distance scenarios that contribute to the hazard from hazard disaggregation, is used to quantify the epistemic uncertainty.To facilitate the interpretation of the map, reference models are also added together with the set of thirteen GMPEs.These reference models are defined as (Kuehn et al. 2015): (1) The average of all models with uniform weights, Aver N w GMPE models is shown in Fig. 5; the S scaled models differ by a constant, the M scaled models are shown with the change in magnitude scaling rate, and the R scaled models are shown with the change in distance attenuation rate.These reference models can be used to orientate the map in an explainable way.First, the average model is set at the center of the map, i.e., point {0, 0}, and then the map is rotated or mirrored so that the up/down scaling is aligned horizontally along the x-axis with the model S++ on the right (along positive direction).Second, the positions of the positive M scaled is oriented along the positive y-axis.Finally, the R scaled models on the map are automatically arranged depending on the positions of the M and S scaled models.
Figure 6a shows the Sammon's map for a set of GMPEs including the candidate models and the original models together with reference models.The map reveals that models CB14adj, BSSA14adj, Id14adj, and Id14 are separated from the other models.Comparison with the reference models shows that CB14adj and BSSA14adj are close to the S--and the S-model, which is the average model down-scaled by 0.67 and 0.8, respectively.This can already be seen in Fig. 1 where the CB14adj is consistently lower than the other models.In general, the reference models help with the interpretation of directions in which GMPEs changed in a systematic way.This allows us to quickly assess the difference between GMPEs over the wide magnitude/distance range used to generate the map.
For the purposes of assessing the changes between the original models and the adjusted models from the map, one can observe in Fig. 6a that pair models CB14adj -CB14, BSSA14adj -BSSA14, and ASK14adj -ASK14 have roughly appeared along the S direction, which indicates that they differ by a constant.The pair model of CY14adj -CY14 has appeared generally in the R direction, which is demonstrated that the model of CY14 is adjusted to scale in distance attenuation.The CY14adj model appeared to be close to the R--model implying that it has stronger attenuation with distance than its original model.The Id14adj has similar attenuation with its original model, but has stronger magnitude scaling as it moves closer the location of M+.On the other hand, the model pairs ASB14adj -ASB14 and Bi14adj -Bi14 have appeared in all scaling references, but their relative distance is small compared to the other pairs.
Another aspect is that one can also identify three GMPEs, i.e., BSSA14, CB14adj, and CY14adj, which are very different in scaling and are oriented in different reference directions on the Sammon's map.Their magnitude scaling and distance scaling is portrayed in Figs.6b and c.The BSSA14 and CB14adj models are more or less parallel in both magnitude and distance scaling because they are approximately up/down scaled versions of each other.That is why they align roughly in the same directions as the S++, S--models.On the other hand, the BSSA14 and CY14 have less similar overall ground motion, but their magnitude scal-ing is very different.Hence, in Figs.6b and c, they differ in their y-coordinate, which corresponds to differences in magnitude scaling, as evidenced by the reference models.The Sammon's map shown in Fig. 6a can help us to discover regional differences in associated median GMPE predictions.Under the scenarios considered in Fig. 2, one can construct three ellipses by using: (1) the candidate models, (2) the original models, and (3) combined the candidate models and the original models.The ellipses are subjectively constructed using the fit to the polygons (closed-dashed lines) from individual group of GMPEs.The range shown in Fig. 6d can cover all the models from each individual group.For example, considering an ellipse that encloses all the candidate models, the ellipse covering all candidate GMPEs is used to quantify the range of epistemic uncertainty associated to median prediction for Taiwan PSHA.It defines the mean to be at the point {-0.14, 0.1}, and the range to be an ellipse with a half axis of length 0.47 in the major direction and 0.22 in the minor direction.On the other hand, Fig. 6d also reveals that the epistemic uncertainty range is not captured adequately by the candidate GMPEs because of the existing gaps between the GMPEs.Hence, to fully obtain the range of epistemic uncertainty, new models need to be selected or the candidate models need to be combined, adjusted, or even newly generated.The Sammon's map is a good tool that can be used to identify missing models or areas that are not captured by existing GMPEs.
Moreover, the range shown in Fig. 6d using the candidate models, a possible partition within the selected uncertainty range can be generated.Three concentric ellipses can be constructed and a total of 17 sub-regions can be developed.Such a partition could be used to define the number of models needed to adequately capture the range of epistemic uncertainty.Since the generated sub-regions cover the range deemed appropriate, one objective was to generate or select a set of GMPEs for those sub-regions such that each sub-region is represented by one model.This could ensure that a good coverage of the epistemic uncertainty is associated with GMPEs, since the map is a representation of the continuous GMPEs space.In the context of Fig. 6d, conceptually, the inner ellipse corresponds to the technically defensible GMPEs; if one wants to assign the center of the ellipse with higher weight, one could define a function f that produces larger values in the center and falls off towards the outer range.Then, the weight for sub-region i becomes the expectation of the function over its area A i .
In other words, the function f could define the relative importance of the center, body, and range, which will be discussed in the following section.However, a scheme based on Sammon's map and Eq. ( 4) provides a way to weight GMPEs, with a clear understanding that each of these models represents a certain part of the range of epistemic uncertainty in GMPE predictions.

IdEnTIfy THE rEPrESEnTATIVE SuITE of GMPES
Since the CBR of TDI range provides a basis for assigning weights on the logic tree, therefore, method to evaluate the CBR of TDI range on Sammon's map needs to be discussed.It has already been discussed from the previous section that there is a need for generating new models to identify a suite of models.In this study, the development of the representative suite of GMPEs is based on a continuous distribution of the median GMPE models in which these new models can be sampled from the mean and covariance of a set of coefficients.This set of coefficients are the result from the fit to the candidate GMPEs using a common form model.The common form used in this study is described below.
, , ) Where, SA ref refers to prediction at a reference rock site condition with V S30 = 760 m s -1 , M is the moment magnitude, R rup is the rupture distance, and Z tor is the depth to the top of the rupture.Cases for the strike-slip, reverse, and normal events are also considered.The styles of faulting factors are simply identified as the magnitude-independent parameters 10 i and 11 i even though some of the GMPEs provide magnitude-dependent styles of faulting factor.The constraint condition was imposed on the fit to remove non-physical models.The regression details are described briefly in the following context.
The first step in the derivation of the common for model was to generate a set of ground motion prediction values from candidate GMPEs.These prediction values are estimated based on the magnitude -distance scenarios that are important to Hazard.Therefore, for each of the candidate GMPEs, median predictions are calculated for M5 to M8 with 0.2 magnitude increment, and also for M = 5.5, 6.5, and 6.5 due to capturing a break in magnitude scaling of some candidate GMPEs.The distance bins relevant to hazard is short distance ranged from 0 to 70 km.For this reason, we sampled short distance (below 70 km) with 27 distance scenarios while 4 distance scenarios were set for distance greater than 70 km.The maximum distance in the fit was 200 for capturing the curvature of candidate GMPEs.
Due to the fact that there much more reverse and strike slip events than normal events occurred in Taiwan.This condition guides us to consider a proper set of prediction values necessary for estimating style of faulting effect in the common form model.We set 80 percent number of scenarios resulted from reverse and strike slip fault with magnitude from M5 to M8, while 20 percent of which was used for normal faulting event.Also, due to Taiwan's geology, crustal earthquakes occurred in Taiwan with Z tor -depth up to 50 km or more, and this might cause different ground motion characteristic of Taiwan region as compared with other regions  while worldwide events were recorded with Z tor -depth less than 20 km.The consideration of deeper Z tor -depth for Taiwan was carried by using scaled factor from Z tor -Magnitude relationship of Chiou and Youngs (2014).The scaled factors were 0.5, 1, 2, 3, and 5.The Z tor -M relationship for the reverse faulting and for the combined strike-slip and normal faulting are described as equations followings: (1) For the reverse faulting model: All other predictor variables, such as the depth to a shear wave velocity horizon of 1000 m s -1 , are set to the default values for the GMPE that requires.For each value of R jb and Z tor , the R rup can be computed from the fault geometry created for each of the scenarios described above.The scenarios are defined from the footwall site for the general case: Hence, each GMPE is evaluated at 19 magnitudes, 31 distances, 5 Z tor leading to 2945 scenarios for the fit.
All the coefficients ( i i ) are determined by fitting common form function to each of candidate GMPEs using least square regression method.The regression results were heavily influenced by initial values of all coefficients and the scenarios used to carry the fits.Therefore, the fits were carried out with a number of steps by resetting the initial values, and the scenarios are generated with heavier samples for short distances than long distance.The coefficients 6 i and 7 i are first determined using scenarios with distance less than 50 km, and rest of the coefficients are determined by the whole scenarios after constraining coefficients 6 i and 7 i .It is not necessary to determine coefficients 8 i for ASB14adj model, so we only use scenarios with distance less than 70 km.
Figure 7 shows the results of fitting to the common form model for ASK14adj model and Phung17, respectively.The results were shown to be good fit implying that the main feature of candidate GMPEs was captured by the common-form model.Through the fit, the candidate GMPEs used different functional forms are thus represented by a set of common form models creating a distribution of coefficients population.Sampling a set of new coefficients is carried out using the multivariate normal density function with the calculated mean ( n i ) and covariance matrix (R i ) multiplied by a factor of 2 which can be used to broaden the range of the common form models rather than only filled the gap among the candidate GMPEs.
Figure 8 shows the prediction of the candidate GMPEs (for the PSA at T = 0.01 sec) in comparison with the range of the 2000 sampled common form models for magnitude scaling (V S30 = 760 m s -1 and R rup = 30 km) and for the distance scaling (V S30 = 760 m s -1 and M = 7).It is observed that the generated common form models can also cover all the candidate GMPEs.Figure 9 shows the corresponding map for 2000 common form models mapping together with the candidate GMPEs.The range of the common-form models is seen to be capture the range of candidate GMPEs on the Sammon's map space.It is thus created a continuous distribution of the median prediction as it can be visualized on the Sammon's map.The issue is then raised of how to adequately quantify the epistemic uncertainty that represents the center, body, and range of technically defensible interpretations (the CBR of the TDI) of the GMPEs.It is believed that, as described in the preceding section, the ellipse was an adequate representation of the general shape of the continuous model distribution and can be used to partition the Sammon's map space into a small number of sub-regions.The concept of developing the ellipse is based on: (1) the ranges that are created by the candidate GMPEs and the candidate GMPEs with their epistemic uncertainties (Atik and Youngs 2014) on Sammon's map.(2) Each of the reference directions on which the scaling of GMPE models should be based is assumed to be followed by a line.(3) The ellipse can be constructed based on the position of the candidate GMPE and their uncertainty on each reference direction (details are shown in Appendix B). Figure 9 shows the candidate GMPEs (that call GMPEs with big color dots) and the candidate GMPEs plus/minus two sigma epistemic uncertainties (that call GMPEs ± 2σ AY14 with small color dots).The use of ±2σ AY14 corresponds to 95% confident interval for the uncertainty in the median ground motion.One can select different sigma level depending on their own judgments.As shown in Fig. 9, two convex-hulls are represented by the dash lines.Each of which was constructed by connecting the outermost points (blue asterisks) which are the projection of the outermost models on each reference directions (M-dir, R-dir, and S-dir).The inner and outer convex-hull is for the GMPEs and the GMPEs ± 2σ AY14 , respectively.The inner convex-hull leaves the range mostly covering the GMPEs, but the outer convex-hull can cover all the GMPEs and some degree of uncertainty.The shape of convex-hulls can guide to find the CBR of TDI range, and an ellipse was created as it goes through the pairs (blue asterisks) corresponding to the GMPEs and the GMPEs ± 2σ AY14 projected on the reference directions.The ellipse was shown in Fig. 9a.The minor axis is oriented by the R-dir.This is deemed valid because the R-dir indicates the change of model scaling in attenuation rate at far distances, while the contribution of far distance scenario to hazard is not significant.On the other hand, the major axis is dominated by the M-dir and the S-dir since the change GMPE models in magnitude and up/down scaling are most contributed to the hazard.The CBR range developed such the way leads to a suite of common form models that captures the range of all the candidate GMPEs with additional epistemic uncertainty in terms of their position in the Sammon's map space.
In the previous section, the ellipse range was divided into 17 sub-regions in which a single representative common-form model is given to each sub-region.Figure 9b shows the CBR range portioned by the same approach.In Fig. 9b, each ellipse plays a different role: the center ellipse corresponds to the center (C), the middle ellipse corresponds to the body (B), and the outer ellipse corresponds to the range (R) of the technically defensible interpretations (TDI) of GMPEs.A single representative model for each cell could be (1) the model whose hazard curve is closest to the mean hazard curve for a cell, (2) the average model in a cell, (3) the model is closest to the centroid of a cell, and (4) the model having either the smallest mean between event residual or the largest log-likelihood is represented in a cell.A single representative model, which is closest to the centroid of that cell, is selected from the common form models in each cell, as shown in Fig. 9b.In this way, the suite of seventeen representative models can represent a set of mutually exclusive and collectively exhaustive models, since the weights on logic tree branches are usually treated as probabilities.It is believed that the current approach could be one possible choice, which differs from that of SWUS GMC SSHAC LEVEL 3 (GeoPentech 2015).Figure 10 shows the comparison of magnitude scaling and distance scaling between the 17 representative common form models (gray lines) and the candidate GMPEs.The selected representative common form models are tightly clustered inside the candidate GMPEs.This is the desired outcome, which fills the existing gap that the candidate GMPEs model left on the Sammon's map of the ground motion space.

dETErMInATIon of WEIGHTS for A SuITE of rEPrESEnTATIVE ModElS
The weights for the median models can be estimated using the observed ground motion data (data-driven weight) and/or the existing models on the GMPE space (nondatedriven weight).The log-likelihood, the mean between event residuals and the EDR index are the different metrics used to measure the goodness of fit between a set of models and observed ground motion data (Scasserra et al. 2009;Kale and Akkar 2013;Arroyo et al. 2014).The likelihood is actually the probability that indicates the variability of the data among the prediction from the GMPEs, while the mean between-event residuals represent the mean bias of the data with respect to the median prediction.The EDR index considers the combination of aleatory variability and model bias.The bias between the median ground motion estimation and observed data is taken into account by a scheme similar to residual analysis.These concepts make it different with respect to the log-likelihood and residual analyses.Unlike the residual weight and likelihood weight, the prior weight are non-data driven weight, which is derived from the normal variate distribution since the mean and covariance of common form coefficients are known.The prior value is actually the probability density of the common form parameters distribution.The prior value explains how a common form model is likely to be generated on the GMPE space.
These metrics are the informative basis for evaluating the weights on the logic tree.As for data-driven weights, the residuals, the likelihood and the EDR are used.On the other hand, the distribution of common-form model parameters is used to calculate non-data driven weight.

computation of residual, log-likelihood and Edr-Index
Computation of the log-likelihood, mean betweenevent residuals and EDR-index requires observed data that is collected from the sites relevant to the hazard.Figure 11a shows the magnitude-distance distribution of data which is made up of the global (the NGA-West2 database) and Taiwan data.The selected data is then corrected to the reference site V S30 of 760 m s -1 ; this is necessary because all the common form models were generated from the refitted candidate GMPEs with a specific V S30 of 760 m s -1 .We performed the correction of earthquake records for site response to the reference V S30 = 760 m s -1 using two candidate GMPEs that are ASK14adj and CY14adj.These two models were selected because their Hanging-wall model can be used to correct the sites subjected to hanging-wall effect.
In addition, considering multiple GMPEs can be captured epistemic uncertainty associated with site effects.The procedure for site correction was then performed by estimating two sets of prediction values relevant to measured V S30 and reference V S30 of 760 m s -1 given M, R rup , and hanging-wall flag.The site amplification factor was determined between two sets of prediction value, and thus the ground motion at reference site condition with V S30 = 760 m s -1 is equivalent to the observed ground motion differed by the site amplification factor.Figure 11b shows the plot of lnPSA(T = 0.01 sec) -R rup for the observed data (hollow circles) compared with the corrected data at reference site (color dots) after filtering the criteria described in Table 2.The ground motion at reference site condition is equivalent to the observed ground motion differed by the site amplification factor.Since data and models are available, the mean between-event residual and the log-likelihood can be calculated using Eqs.( 7) and ( 10) of Abrahamson and Youngs (1992) as described below.
(2 ) ( ) ( ) , and ij n is the predicted value corresponding to observation y ij .
The calculation was employed a fixed inter-event standard deviation ( x) and fixed intra-event standard deviation (z) obtained from CY14adj.The EDR-index is estimated using Eq. ( 11) of Kale and Akkar (2013).Figure 12 respectively shows the contour plot for the log-likelihood, the mean between-event residuals and the EDR-index of the selected dataset that are generated by each common form model.The contours from the log-likelihood, as shown in Fig. 12a, are distinguished by ten distinctly colored spans from the smallest to the largest likelihood value (-2509 to -2076).The contours from the residuals, as shown in Fig. 12b, are separated by a number of bands having a constant width of 0.15 unit and varying from -0.8 to 0.9. Figure 12c shows the contours for the EDR-index with a span color of 0.1 unit varying from 0.82 to 1.77.Figure 12d plots distribution of the common-form models with no contours, indicating that the weight is derived from the models themselves and not from the data.

determination of Weight Associated with Each representative Model
The weight associated with each of the representative model (or the weights for median), used in GeoPentech (2015) is defined in Eq. ( 10).The weight is proportional to its area and weighted by the mean value of the selected metric (L ij ) as follows: where A i is the area of the i th cell and N j is the number of the model in the i th cell.The metrics L ij are considered following: (1) One over the absolute mean between event residual plus a small constant c (e.g., c = 0.0075 to avoid a possible singularity).This is equivalent to (3) The shifted EDR-index values by a constant (4) The prior values, is the probability density of the common form parameters.
( , , ) It is noted that the direct use of the log-likelihood value or the EDR value in the computation may lead to an indistinguishable weight; that is, all seventeen representative models can receive a portion of equally about 5.88%, and this is not desirable.In other words, if the likelihood weight (wLL) computed from the original log-likelihood value (without C), the result was not significant importance of the center body and range of TDI of GMPE prediction.It is therefore suggested that both the log-likelihood and the EDR values are subtracted by a constant C to shift log-likelihood/the EDR values to those that can be used to refine the weights.Such a constant C is the average metric value and can be determined from the models outside the outer polygon.As a result, we propose to use constant C to shift the original log-likelihood values as it can be used to normalize the weights giving either adequately distinguishable or more informative.
Figure 13 shows the weighting scheme from different metric.The residual weight (wR), tends to distribute a higher portion to the sub-regions close to the zero contour (thickblack contour) where the models fit the data better; those are the sub-regions S9, S1, S8, S17, S7, S5, S6, and S13 contributing up to 68.57% of the total weight.The log-likelihood Table 2. Data selection criteria for crustal source.
Note: Select events with at least 5 recordings.Select events with at least 1 recording within 20 km.Remove Chi-Chi main shock and aftershock events from the NGA-West2 database event.Fig. 13.Comparison on the weighting schemes calculated using the mean between event residual, the log-likelihood, the EDR, and the prior value.
weight (wLL) and the EDR weight (wEDR) are comparable with the similar trend.The weighting schemes produce larger values in the center and falls off towards to the outer range.
The results show that these weighting schemes conform with the log-likelihood or the EDR contours and gives higher portions to the inner ring: 8 cells including S17, S1, S2, S3, S4, S5, S6, S7, and S8 (up to about 63% of the total weight).Unlike the residual weight and likelihood weight, the prior weight is a non-data driven weight, which is derived from the normal variate distribution since the mean (μ θ ) and covariance (Σ θ ) of common form coefficients are known.The prior weight reflects how often the generated common form models appear in a particular location.Indeed, the cells S8, S17, and S16, where the models are more densely sampled, have higher portions of 12.41, 10.75, and 10%, respectively.The cells S11, S14, S13, and S12, conversely, have few models and thus received lower weights with portions of 1.01, 1.48, 2.35, and 2.52 %, respectively.The contour maps as shown in Fig. 12 for mean between event residual, log-likelihood and EDR-index are consistent with each other.Given a number of GMPEs, it is suggested that a GMPE model receive a relatively high weight if it was examined to fit the observed data; in other words, a model whose log-likelihood is larger or its residual (or EDR) value is smaller in comparison with others is considered as a good prediction and thus given a higher weight.

An Alternative Way to Evaluate the Weight for the Median
The idea of method is based on the distribution of calculated metrics of the continuous models on the Sammon's map.The mean between event residual, log-likelihood, EDR and the prior values were calculated from the common form models and dataset, as presented in Fig. 14a.The color dots indicate the values of candidate GMPEs, and the black dots are the average metric value over a sub-zone.Corresponding to the bar chart, the cumulative density function can be obtained and shows in Fig. 14b.Each of the black dots takes a value from the cumulative density function, and thus results in a set of seventeen numbers, which are normalized to one to obtain the weight.Figure 15 shows the proper weight scheme on which each of the component of the proposed approach is being compared with that calculated based on the formula Eq. (10) (i.e., wR, wLL, wEDR, and wPri).The result can be used to check, as compared to the weight calculated based on the formula in Eq. ( 10), whether the weights are properly assigned to each sub region.The alternative method leaves an aspect that is especially beneficial in the context of probability.On the one hand, the method does provide a cumulative distribution function (CDF) curve that quickly help to quantify the weight given to candidate models.For example, on Fig. 14b, the residual weight of Chao17 and Bi7adj models could be given higher portions as com-pared to representative models as well as the other candidate models.On the other hand, the CDF for the log-likelihood shows that the model of CY14adj, ASK14adj, Phung17, and Chao17 appear up to above 80% and higher than the rest models indicating the dominant weights of those models if one considers the weights possibly given to the set of candidate models.Therefore, the present approach offers clear advantages over previous attempt related to the weights.

concluSIonS
In this paper, the work of GeoPentech (2015) in the SWUS 2015 report was carefully studied, re-examined, and applied to the Taiwan region to develop the representative GMPE as well as the weighting schemes for PSHA study.In the analysis, the selected candidate GMPEs were reviewed with great care, and then adjusted or developed for use in the Taiwan region.The Sammon's mapping technique, developed by Scherbaum et al. (2010), was then applied to assess the spread of median GMPE prediction and correlation through 2-dimensional visualization.Such maps can qualitatively provide insight on the epistemic uncertainty associated with GMPEs, by weighting the GMPE-distances through informative map distances.In the generation of the map, the hazard de-aggregation was assigned to weight GMPE distances.A set of GMPEs models including both the original models and candidate models together with the reference models were projected together on the same 2-D Sammon's map to assess differences between GMPEs over an important magnitude-distance range.It is shown that the changes in GMPE distance (either uniform weighted or de-aggregation weighted GMPE distance) for each pair of GMPEs (original GMPE vs. candidate GMPEs) are visually explained by looking at the reference models.PSHA is interested in the hazard space; therefore, the de-aggregation weighted GMPE distance is suggested to have the final hazard outcome.
To generate the maps in this paper, a large suite of new models (2000 random GMPEs) were created from the candidate GMPEs, which both interpolate between and extrapolate beyond the range of the existing GMPEs.The ground motions of the candidate GMPE models and 2000 random models were combined into a high-dimensional vector of the ground motion values, which was then projected into the 2-D Sammon's map.The median GMPE predictions in the Sammon's map approach created a continuous distribution values as a proxy to assign the weights on the logic-tree.
The representation of the candidate GMPEs (GMPEs) and additional epistemic uncertainty around the median (GMPEs ± 2σ AY14 ) provide a widely enough range on 2-D Sammon's map, and an ellipse was adequately considered as a general representative shape to enclose the outermost points of all GMPEs.The ellipses were guided to define the range and body of the distribution of the common-form models and discrete the range into a small number of subregions in which a number of representative common form models can be selected for each sub-region.Different weighting schemes were presented, based on comparisons to ground motion data and on the distribution of the covariance matrix for the common form coefficient weights, to assign weights to the selected representative common form models developed from each sub-region.The weights assigned to the inner cells are equally likely to be stable, but those assigned to the outer cells vary significantly from one cell to another.Finally, the developed weight as well as the representative GMPE can be used for PSHA.
The weights formula used in GeoPentech (2015) for the median GMPEs were verified through computation for a number of metrics which either the likelihood, the EDRindex and residuals derived from models and data or the distribution of model parameters based on mean and covariance of model's coefficients.Regarding to the likelihood weight, the weight cannot be directly obtained from the raw log-likelihood.It is therefore a constant C was introduced to find a shifted value so that the log-likelihood and/or the EDR weight can be meaningful and justifiable.The results showed that the weight scheme does lead to adequately explainable with regards to center body and range of technically defensible interpretation of GMPE prediction.We propose an alternative way to re-evaluate the weighting schemes based on the cumulative distribution function curve, and this method provide the benefit enable to quick assess the weights not only for the suite of representative models but also for the candidate models.The weighting schemes obtained by this new approach can also be compared with the weight using the formula in Eq. ( 10) and thus achieved the goal of PSHA study.Besides, the EDR index was tested on the Global + Taiwan dataset using 2000 common-form models.The derived weight using EDR index performs comparable with the log-likelihood weight and the residual weight.Therefore, we concluded that the EDR index is a good metric for estimating the logic tree weight, and this appropriate metric could consequently be used as an additional branch on the logic tree.The method is not the only option to correctly define the weighting of GMPEs for PSHA; however, it offers a systematic way to combine different sources of knowledge, such as recorded data and prior information for PSHA. Figure 16 shows an example of logic tree weight assignment for the median base models utilizing the presented method.In fact, when performing PSHA in a practical way, a sensitivity analysis discussing different choices for weight assignment to each branch of the logic tree and how those changes affect the result of the hazard curve needs to be presented.

data and resources
This study used the ground-motion database compiled for the Taiwan SSHAC Level 3 PSHA project.The strong ground-motion data in the SSHAC GMC Ground-Motion database version 4-20161228 is provided by NCREE and SINOTECH.All calculations are carried out using the MATLAB 2017 program.
Fig. 1.(a) Tectonic setting of Taiwan (Wu et al. 2013), (b) general data used for this study.

Fig. 2 .
Fig. 2. Plots showing for (a) magnitude scaling and (b) distance scaling at the R rup of 30 km and M of 7.0 for the periods of T = 0.01 sec (or PGA) under vertical strike slip earthquake and engineer rock site conditions (V S30 = 760 m s -1 ).The model abbreviations given in the legend are used henceforth.
Fig. 4. Sammon's map showing for comparison of GMPEs: (a) based on un-weighted GMPE-distance, (b) based on weighted GMPE-distance.The original GMPEs are denoted by colored diamonds, and the candidate GMPEs are denoted by colored dots.

Fig. 5 .
Fig. 5. Difference in scaling of reference models that help the interpretation of the Sammon's map: (a) reference models for up-down scaling, (b) reference models for magnitude scaling rate, and (c) reference models for distance attenuation scaling.
Fig 6.(a) Sammon's map of GMPEs based on weighted GMPE-distances, together with thirteen reference models.The original GMPEs are denoted by colored diamonds, and the candidate GMPEs are denoted by colored circles.(b) and (c) Magnitude scaling and distance scaling of three GMPEs for a detailed interpretation of the relative differences.(d) The possible partition of selected range of epistemic uncertainty shows by the gray lines.

Fig. 7 .
Fig. 7. Comparison between the refit to common form model and the candidate model: (a) ASK14adj and (b) Phung17 for T = 0.01 s.
Fig 9. (a) The ellipse range used to define the CBI of TDI.(b) Selected seventeen representative common form models (black dots) for each cell.
Fig. 10.Comparison of (a) magnitude scaling and (b) distance scaling between the representative common form model and the candidate GMPEs.

Fig. 12 .
Fig. 12.The contour lines underlying the Sammon's map indicate; (a) the value of the mean between-event residual, (b) the likelihood value, and (c) the EDR value computed for each model.(d) The distribution of common form models.

Fig. 15 .
Fig. 14.(a) Bar chart showing histogram of log-likelihood of the models on the Sammon's map.(b) Cumulative density function obtaining from the histogram.The color dots indicate the likelihood value of each candidate GMPEs, and the black dots are the average metrics over one sub-region and given to that sub-region.