- Email: [email protected]
Updating a probabilistic seismic hazard assessment with instrumental and historical observations based on a Bayesian inference
Nuclear Engineering and Design 350 (2019) 98–106
Contents lists available at ScienceDirect
Nuclear Engineering and Design journal homepage: www.elsevier.com/locate/nucengdes
Updating a probabilistic seismic hazard assessment with instrumental and historical observations based on a Bayesian inference
T
E. Vialleta, , N. Humbertb, P. Mottierc ⁎
a
EDF, Nuclear New Build Engineering Division, France EDF, Hydropower Engineering Division, France c Polytechnique Montreal, Canada b
ARTICLE INFO
ABSTRACT
Keywords: Probabilistic seismic hazard assessment PSHA Bayesian inference Bayesian updating
Since the basic work of Cornell, many studies have been conducted in order to evaluate the probabilistic seismic hazard (PSHA) at a given site or at a regional scale. In general, results of such studies are used as inputs for regulatory hazard maps or for risk assessments. Such approaches are nowadays considered as well established and come more and more used worldwide, in addition to deterministic approaches. Nevertheless, some discrepancies have been observed recently in some PSHA, especially from studies conducted in areas with low to moderate seismicity. The lessons learnt from these results lead to conclude that, due to epistemic uncertainties inherent to such a domain, some deterministic choices have to be made and, depending on expert judgments, may lead to strong differences in terms of seismic motion evaluation. In this context, the objective of this paper is to present a methodology that can be used to take into consideration instrumental and historical observations in order to reduce epistemic uncertainties in a PSHA. The method developed here is based on a Bayesian inference technique used in order to quantify the likelihood of the prior estimation, and finally update the PSHA. The period of observation under consideration is the completeness period for each set of observation used. The updating process is developed at a regional scale, over a significant number of stations. The potential correlation between points of observation is also discussed and accounted for. Finally, a case of application is proposed on the French metropolitan territory to demonstrate the efficiency of this updating method and draw perspectives for further applications.
1. Introduction and background
consideration, cf. Eq. (1).
Bayes theory on conditional probabilities was developed in the second half of the 18th century by Reverend Thomas Bayes and was strongly promoted by Pierre-Simon de Laplace in the early 19th century. The use of Bayesian theorem of conditional probability (also known as Bayesian inference), which has drastically increased during the second half of the 20th century, is nowadays largely disseminated in a wide range of fields and applications such as biology, pharmaceutical statistics, economics, finance, environmental and earth science, industrial statistics, mechanics (ISBA, 2018) … The principle of Bayes theorem of conditional probabilities is simple to present (but not always simple to solve!). Based on a probabilistic assessment of a given parameter A (usually called as the prior estimation), Bayesian theorem consists in calculating the likelihood of the prior estimation using actual observations B of the parameter under
P (A |B ) =
⁎
P (B |A). P (A) P (B )
(1)
where P(A) is the prior estimation based on a probabilistic assessment P(B|A) is conditional probability of the observed event B according to the prior assessment P(B) is the total probability of the observed event B according to the prior assessment P(A|B) is the expected updated (posterior) estimation
Corresponding author at: EDF, Nuclear New Build Engineering Division, Direction Technique, 19 rue Pierre Bourdeix, 69007 Lyon, France. E-mail address: [email protected] (E. Viallet).
https://doi.org/10.1016/j.nucengdes.2019.04.034 Received 29 December 2018; Received in revised form 22 April 2019; Accepted 24 April 2019 Available online 20 May 2019 0029-5493/ © 2019 Elsevier B.V. All rights reserved.
Nuclear Engineering and Design 350 (2019) 98–106
E. Viallet, et al.
2. Application of Bayes theorem in the field of nuclear industry and earthquake engineering
3.2. Prior estimation: PSHA complemented by PGA-to-Intensity relationships
In the field of nuclear engineering, Bayesian approaches have been used for a long time, especially in the field of reliability analyses (Colombo, 1985) and mechanical and structural engineering (Madsen, 1985; Heinfling et al., 2005). In the field of earthquake engineering, Bayesian approaches are also largely used, in various aspects, such as seismology and structural dynamics. In seismology, the Bayesian inference is used in elicitation techniques (Runge et al., 2013), for comparison of PSHA hypothesis (Humbert et al., 2008b; Secanell et al., 2017) or for direct updating of GMPEs on local data (Ordaz et al., 2007; Wang et al., 2009; Stirling et al., 2010). In the last decade, many developments and applications were published in the field of PSHA testing and Bayesian updating, especially in the framework of OECD/NEA (OECD/NEA, 2008) and SIGMA Project (Keller et al., 2014; Tasan et al., 2014). PSHA testing and Bayesian updating are now considered as state of the art and should be systematically implemented. This process is not antagonist with expert judgments, which are of course necessary to guaranty the robustness of the prior PSHA and the correct propagation of random and epistemic uncertainties. PSHA testing and Bayesian updating have to be considered as an additional step that allows actual observations at a given site to be taken into consideration to reduce epistemic uncertainties (i.e. uncertainties that could not be reduced by experts), taking into account random and epistemic uncertainties they intrinsically incorporate (OECD/NEA, 2015).
Theoretically, any kind of observation and any period of observation, including instrumental seismicity, historical seismicity and paleoseismicity data can be used in such a process. In this section, we’ll focus on historical seismicity. The objective is to determine an exceedance rate of a given macroseismic Intensity based on the hazard curve under consideration. The principle of the method was described previously (Secanell et al., 2017) and extended in this paper as described below. In order to use historical observations, PGA-to-Intensity relationships are used to transform PGA to macroseismic intensities, based on Eq. (2).
I = a . Log (PGA) + b ( + /
)
(2)
Random and epistemic uncertainties are taken into consideration through the full propagation of the standard deviation of the PGA-toIntensity relationship and through the use of different PGA-to-Intensity relationships. The goal is to determine the rate of exceedance of a given intensity from a given hazard curve expressed in term of PGA. The calculation is performed based on Eq. (3). Amax
NIi =
(NA +
A
NA ). PIi | A
Amin
(3)
where
NIi is the annual rate of exceedance of a given Intensity Ii NA is the annual rate of exceedance of a given PGA A PIi | A is the probability of a given PGA to produce the given intensity Ii according to Eq. (2), including σ (i.e. percentage of the contribution of the given PGA level to the considered class of intensity) ΔA is a discretization step small enough in order to get stable results of NIi (in this study, the range of PGA from 1 to 1000 cm/s2 is discretized into 40 steps)
3. The innovating method developed in this paper The method developed in this paper is based on the Bayesian updating process described in previous papers (Viallet et al., 2008; Humbert et al., 2008b). In these previous references, instrumental seismicity is taken into consideration as observations to update the prior estimation. In addition, a previous work (Viallet et al., 2017) demonstrated that historical seismicity can also be taken into consideration through the same Bayesian process. In the current paper, previous works are extended and the Bayesian inference is applied based on both historical and instrumental seismicity, in order to update the prior PSHA. The objective is thus to demonstrate that such an approach can be applied combining multiple sources of observation (instrumental and historical seismicity in this case) and to show that the proposed process has an improving impact on the weighting process of the prior PSHA branches.
Finally, the annual exceedance rate is multiplied by the completeness period of the historical catalogue in order to determine the prior rate of exceedance of a given hazard curve in terms of macroseismic Intensity. The illustration of this process is given in Fig. 1. This allows to express the prior PSHA estimation in term of rate of exceedance of a given Macroseismic Intensity. Through this process, the full range of uncertainties is propagated in order to guaranty that none of the possibilities are disregarded.
3.1. Principle for incorporating historical observations
3.3. Historical observations consideration
First of all, the approach needs a PSHA following the state of the art approach as input data. This prior estimation has to take into consideration all sources of random and epistemic uncertainties. This PSHA has to be carried out though a well-documented procedure, which should explicitly describe the elicitation of experts’ judgments, such as the so-called SSHAC procedure. Hazard curves expressed in term of PGA (Peak Ground Acceleration) rate of exceedance will be considered as the prior estimation for the implementation of Bayes theorem of conditional probability. The hazard curve is used to determine an exceedance rate of a given macroseismic Intensity based on Intensity to PGA relationships. Random and epistemic uncertainties are fully propagated, as described in the following sections. Historical observations are expressed in term of macroseismic Intensity observed at the site under consideration. Finally, the Bayesian theorem of conditional probability is applied for each single hazard curve of the prior estimation in order to calculate its likelihood and though to update its weight accordingly. The likelihood function becomes the posterior weight of the hazard curves.
The historical observations have to be well documented. This means that a historical database has to be available and its quality has to be high enough to allow to define completeness period and to be sure that it is possible to quantify the number of event that were felt on a given location over this completeness period. Once the previous requirements are fulfilled, the historical observations are the number of historical events that produced a Macroseismic Intensity higher or equal to a given one. An application is presented in the next section. We would like to mention that the historical observation which are used here are also usually used in the PSHA as input data (for earthquake recurrence parameters’ determination for instance). This means that these data could be used twice in the process, which may be confusing. However, we are dealing here with a reduced number of historical observations compared with the huge number that are considered in the prior PSHA. Then, it can be easily shown that removing these observations from the set of data used in the prior PSHA would have no impact on the PSHA result. Then, there is no risk of giving them an arbitrary high weight in the overall process (by using them twice). 99
Nuclear Engineering and Design 350 (2019) 98–106
PGA to Intensity relationship : Illustration
1.E+01 1.E+00
Annual rate of exceedance
1.E-01 1.E-02 1.E-03 1.E-04 1.E-05 1.E-06
Rate of exceedance (over the completeness period : XXX y)
E. Viallet, et al.
35
30
N Ii =
25
Amax
∑ (N
A + ∆A
VII-VIII / 265y
VIII / 515y
Amin
− N A ).PI i
A
20
15
10
5
0 V / 135y
V-VI / 135y
VI / 165y
VI-VII / 165y
VII / 265y
VIII-IX / 515y
IX / 715y
Intensity range (above) + Illustration of PGA contribution (below)
1.E-07 10
100
1000
PGA (cm/s²)
1.0
1.2
1.4
1.6
1.9
2.3
2.6
3.1
3.6
4.3
5.0
5.8
6.9
8.0
9.4
11.0
12.9
15.2
17.8
20.8
24.4
28.6
33.5
39.3
46.1
54.0
63.3
74.2
87.0
101.9
Fig. 1. Process used to transform hazard curve in PGA into macroseismic Intensity. Left: Discretized hazard curve – Middle: PGA-to-Intensity relationship – Right: Intensity rate of exceedance.
3.4. The updating process: accounting for random occurrence of events
multiple sites of observation, the Poisson’s occurrence model Eq. (6) is generalized through a negative binomial distribution that allows to account for correlation between stations of observation, as described in detail in a previous paper (Humbert, 2015) and shown in Eq. (7).
The updating process (Eq. (1)) applied to our case is shown in Eq. (4).
P (Model|Observation) =
P (Observation|Model) . P (Model) P (Observation)
P (Observation|Ai ) = P (n, t ) = exp ln
(4)
where
ln (1 + n) .
P(Model) is the prior estimation, direct output of the PSHA based on a probabilistic assessment P(Observation|Model) is the probability of the observed event according to the prior assessment P(Observation) is the total probability of the observed event according to the prior assessment P(Model|Observation) is the expected updated (posterior) estimation
P (Observation|Ai ) . P (Ai ) P (Observation)
(Ai ). n!
(A ) 1
n
(7)
P (Observation|Ai )
(8)
3.4.3. Likelihood function and posterior estimation Finally, the likelihood function, Eq. (9), becomes the updated weight of the hazard curve Ai under consideration, which allows to recalculate the range of (reduced) epistemic uncertainties of the PSHA, and becomes the posterior estimation.
P (Observation|Ai ) P (Observation)
3.4.1. Calculation of P(Observation|Ai) Based on the rate of exceedance calculated as described in Eq. (3), P (Observation|Ai) is calculated assuming a Poisson’s occurrence model, Eq. (6).
e
1 K
N
P (Ai ) is the prior weight of the hazard curve under consideration, which results directly from the weighting process of the PSHA under consideration, including expert judgment and epistemic uncertainties propagation P (Observation|Ai ) is the probability of the observed event according to the prior assessment P (Observation) is the total probability of the observed event B according to the prior assessment
P (Observation|Ai ) = P (n, t ) =
1
K
3.4.2. Calculation of P(Observation) P(Observation) is simply calculated through Equation (8).
(5)
P (Observation) =
t )n
1
ln
lnΓ is the logarithmic gamma function λ(A) is the cumulative rate of exceedance of a given acceleration A for all the stations under consideration n is the number of events with an acceleration higher that A K is the average number of sites impacted by one earthquake (this parameter is indicative of the correlation among stations. If K tends to 1, the negative binomial distribution goes towards a Poisson distribution).
where
.t (
A
(A ) +n 1
where
This equation is applied to each single hazard curve of the prior PSHA (Ai) so that the updating process will allows to calculate a likelihood of each single hazard curve of that prior estimation, Eq. (5).
P (Ai |Observation) =
1K K
K
(9)
4. Application: French metropolitan territory 4.1. Description of the PSHA used for this application
(6)
In order to illustrate the efficiency of this method, the application is performed based on the French Metropolitan territory. The prior PSHA is based on a previous study (Humbert et al., 2008a). This study is not described in detail in this paper in order to concentrate on the Bayesian updating process. However, this PSHA includes the propagation of a large range of random and epistemic uncertainties, especially on the GMPEs and the random model as shown by Eq. (10).
where λ is the annual rate of exceedance of a given acceleration Ai t is the period of time under consideration n is the number of event with an acceleration higher that Ai In order to perform this updating process at a regional scale using 100
Nuclear Engineering and Design 350 (2019) 98–106
E. Viallet, et al.
Fig. 2. Location of the 19 EDF NPPs on the French Metropolitan territory (green diamonds). NB: Yellow triangles are French accelerometric networks stations not used in this study. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
( ) (±a ) + (±a ). M . a a R + ( ± ) + (± ) + (± )
( )
a1 a2 Log (PGA) = ± .M+ ± . M2 + 1 2 Log
2
2
5
5
3
relationships from different authors were used. The only criterion that was used to select PGA-to-Intensity relationships is the domain of validity in term of PGA that should be large enough to allow to estimate Intensities without any bias (especially boundary effects). These selected PGA-to-Intensity relationships are presented in Table 1 and plotted in Fig. 4.
4
3
4
6
6
(10)
where PGA is the Peak Ground Acceleration, M is the magnitude, R is the hypocentral distance, Every parameter (a1 to a6) has its own uncertainty, even the random part σ.
4.3. Completeness periods for historical seismicity Based on SisFrance historical seismicity database (SisFrance, 2018), the period of completeness of historical seismicity is known for the French Metropolitan territory, as shown in Table 2. This assumption has to be consistent with the one that is made for the prior PSHA, in term of period and uncertainties, if any. This consistency was confirmed in our case. In addition, it should be pointed out that the completeness of observation data is crucial and should be checked and confirmed in order to guaranty that the updated process does not introduce any bias (if observation database is incomplete, the Bayesian approach will artificially increase the weight of branches that are consistent with this biased observation set).
This process allows to obtain hazard curves for a large number of locations. In our case, we are using the results obtained for 19 EDF Nuclear Power Plant locations, in order to use historical seismicity as observation, which is well documented for these sites. The location of EDF NPPs is shown in Fig. 2. For this paper, 42 hazard curves have been generated among the whole logic tree that contains much more than this number. This choice (only 42 hazard curves) was made in order to ease the interpretation and to better illustrate the updating process. In this purpose, the 42 selected hazard curves were randomly chosen among the whole logic tree of the prior PSHA. Fig. 3 shows the 42 hazard curves obtained at one site.
4.4. Observations: historical earthquakes felt on sites Based on SisFrance database and Seismic Hazard Assessment of the EDF’s NPP, it is possible to quantify the number of historical events felt on the 19 EDF sites. This number of observed events, quantified per range of magnitude, is presented in Table 3.
4.2. PGA-to-Intensity transformation In order to propagate a large range of uncertainties, different 101
Nuclear Engineering and Design 350 (2019) 98–106
E. Viallet, et al.
1.E+01
1 592 466 564 626 688 5 143 189 551 557 593 1551 1557 1593 2551 2557 2593 3551 3557 3592
1.E+00
Annual rate of exceedance
1.E-01 1.E-02 1.E-03 1.E-04 1.E-05 1.E-06
62 746 139 363 259 781 30 176 550 554 573 1550 1554 1573 2550 2554 2573 3550 3554 3573 3593
1.E-07 10
100
1000
PGA (cm/s²) Fig. 3. Illustration of the 42 hazard curves resulting from the logic tree under consideration at a given site.
4.5. Logic tree developed for Bayesian updating based on historical observation
is quite comparable for all hazard curves. This illustrates that for such high Intensities, the updating process is not able to rank the hazard curves, simply because at this level of Intensity (which corresponds to long return periods), all hazard curves do not predict any observation, as actually observed. When there is not enough observation, the updating process does not modify prior weights, as also observed in previous studies (Secanell et al., 2015). Finally, the updated weight of the hazard curves resulting from the historical seismicity logic tree is shown in Fig. 9. This figure clearly shows the efficiency of the updating process. The updating process, which does not modify any of the parameters of the prior PSHA, only modifies the weight of the hazard curves. One can finally observe that none of the prior hazard curves are completely eliminated. Their weights are reduced according to the likelihood that is calculated through Bayes theorem of conditional probabilities, according to observations.
The logic tree developed to update the weight of the 42 hazard curves is described in Fig. 5. This process allows to account for different PGA-to-Intensity relationships and to account for different Intensity ranges (i.e. different return periods). 4.6. Results The updated weight (Bayesian likelihood function) of the hazard curves are presented in the next figures. The likelihood function obtained for each PGA-to-Intensity relationship is shown in Fig. 6. In this figure, the height of the colored bar corresponds to the value of the likelihood of the hazard curve under consideration, for the PGA-to-Intensity relationship under consideration. The likelihood function obtained for each Intensity range included in the logic tree is shown in Fig. 7. In this figure, the height of the coloured bar corresponds to the value of the likelihood of the hazard curve under consideration, for the Intensity range under consideration. Fig. 8 “radar view” illustrates how the likelihood evolves with Intensity. This figure shows that for high Intensity (I = IX), the likelihood
4.7. Discussion First of all, we can conclude that the Bayesian updating procedure is successfully implemented and allows to update a PSHA based on actual historical seismicity, as observed on the sites under consideration. We also observe that if observations are not sufficient, the updating process
Table 1 PGA-to-Intensity relationships used in this study. Author
Domain of application
Relationship (PGA in cm.s−2)
σ
Atkinson and Kaka (2006)
Missouri et California I ∊ [II;IX] PGA ∊ [0.4;260] Centre US et California I ∊ [II;IX]
I = 2.315 + 1.319*log(PGA) + 0.372*log(PGA)2
0.93
If log(PGA) ≤ 1.69, I = 2.65 + 1.39*log(PGA)
1.01
If log(PGA) > 1.69, I = −1.91 + 4.09*log(PGA) I = 2.58*Log(PGA) + 1.68
0.35
I = 2.3*Log(PGA(g)) + 10 If log(PGA) ≤ 1.57, I = 1.78 + 1.55*log(PGA)
0.3 0.73
If log(PGA) > 1.57, I = −1.60 + 3.70*log(PGA) If I ≤ V, I = 7.58 + 2.2*log(PGA(g))
1.7
Atkinson and Kaka (2007) Faenza and Michelini (2010) Marin et al. (2004) Worden et al. (2011) Wald et al. (1999) and Atkinson and Sonley (2000)
Italy I ∊ [I;VIII] France California I ∊ [II;VIII-IX] PGA ∊ [0.07;800] California-France
If I > V, I = 10.18 + 4.35*log(PGA(g))
102
Nuclear Engineering and Design 350 (2019) 98–106
E. Viallet, et al.
10 9 Atkinson et Kaka (2007)
Macroseismic Intensity
8 7
Marin & All (2004)
6 5
Faenza & Michelini (2010)
4
Worden & All (2011)
3
Wald, Atkinson&Sonley (BRGM-2000)
2
Atkinson et Kaka (2006)
1 0 1
10
100
1000
PGA (cm/s²) Fig. 4. PGA-to-Intensity relationships used in this study.
Table 2 Completeness periods for historical seismicity. Zone
Intensities IV to IV-V
Intensities V to V-VI
Intensities VI to VI-VII
Intensities VII to VII-VIII
Intensities VIII to VIII-IX
Intensities IX to IX-X
Entire France Alpes Britany North Provence Pyrenees Rhine Other areas
1920 1950 1920 1920 1925 1930 1885 1960
1880 1880 1880 1880 1880 1910 1870 1850
1850 1830 1750 1850 1850 1910 1750 1750
1750 1800 1700 1750 1750 1750 1630 1750
1500 1500 1500 1500 1500 1500 1500 1500
1300 1300 1300 1300 1300 1300 1300 1300
Table 3 Total amount of historical events felt on all EDF NPP sites. Intensity range Completeness period
V 135y
V-VI 135y
VI 165y
VI-VII 165y
VII 265y
VII-VIII 265y
VIII 515y
VIII-IX 515y
IX 715y
IX-X 715y
X 715y
Total number of events (whole 19 EDF’s sites)
24
12
3
2
0
0
0
0
0
0
0
PGA to Intensity relationships Weight Atkinson&Kaka (2006) Weight = 1/6 Atkinson&Kaka (2007) Weight = 1/6
Bayesian updating Logic tree => Implemented for each hazard curve of the PSHA
Faenza & Michelini, 2010 Weight = 1/6 Marin, Avouac, Schlupp, & Nicolas, 2004 Weight = 1/6 Worden, Gerstenberger, Rhoades, & Wald, 2011 Weight = 1/6 Wald et Atkinson & Sonley (BRGM-2000) Weight = 1/6
Intensity Weight
I = VI Weight = 1/4 I = VII Weight = 1/4
Updated weight
I = VIII Weight = 1/4 I = IX Weight = 1/4
Fig. 5. The logic tree developed to update the weight of the 42 hazard curves based on historical seismicity.
103
Nuclear Engineering and Design 350 (2019) 98–106
E. Viallet, et al.
Baysian updating - Historical seismicity Effect of Intensity - PGA relationship Atkinson et Kaka (2006)
0.6
Updated weight
0.5
Wald, Atkinson&Sonley (BRGM-2000) Worden & All (2011)
0.4 0.3
Faenza & Michelini (2010) Marin & All (2004)
0.2 0.1 0 1
259 564 688
5
Atkinson et Kaka (2007)
176 551 573 1551 1573 2551 2573 3551 3573
Hazard curve id Fig. 6. Likelihood function depending on PGA-to-Intensity relationships.
Baysian updating - Historical seismicity Effect of Intensity range 0.4
Updated weight
0.35 0.3 0.25
Weights (I=IX)
0.2
Weights (I=VIII)
0.15
Weights (I=VII)
0.1
Weights (I=VI)
0.05 0 1
259
564
688
5
176
551
573 1551 1573 2551 2573 3551 3573
Hazard curve id Fig. 7. Likelihood function depending on range of Macroseismic Intensity used.
35923593 3573 3557 3554 3551 3550
1
Baysian updating - Historical seismicity
62 139 259 363 466 564 592
2593
Effect of Intensity range
626
2573
688
2557
746
2554
781
2551
5
2550 1593 1573 1557 1554 1551 1550593
30
573
143 176 189 550 551 557 554
Weights (I=VI) Weights (I=VII) Weights (I=VIII) Weights (I=IX)
Fig. 8. Likelihood function depending on range of Macroseismic Intensity used – Radar view.
does not affect the prior estimation (typically the case for Intensity IX in our study). This means that there is no risk of performing updating: the process is self-reliable. This also means that even with few (even no) observation, Bayesian updating can be implemented. In addition, it is important to point out that the updating process does not modify any of the assumptions or input data of the prior estimation. It only updates the weight of the hazard curves (output of the PSHA), based on observations. This process allows to reduce the
epistemic uncertainties that are contained in the prior model, based on their intrinsic (and objective) likelihood. Finally, the procedure developed in this paper, based on historical seismicity, takes into consideration all sources of uncertainties (random and epistemic) in the PGA-to-Intensity transformation process and in the random occurrence of events, which avoid any deterministic or inappropriate choices…
104
Nuclear Engineering and Design 350 (2019) 98–106
E. Viallet, et al.
Updated weight 0.1
Updated weight
0.075
Historical seismicity logic tree
0.05
0.025
0 1
259
564
688
5
176
551
573 1551 1573 2551 2573 3551 3573
Hazard curve id Fig. 9. Historical seismicity logic tree (6 PGA-to-Intensity relationship + 4 ranges of Intensity): Updated weight.
Instrumental seismicity - PGA threshold = 0.01g – Weight = 0.2 Bayesian updating Logic tree => Implemented for each hazard curve of the PSHA
Historical seismicity (see dedicated logic tree in fig. 5) Total weight = 0.8 - Intensity = VI Weight = 0.2 - Intensity = VII Weight = 0.2 - Intensity = VIII Weight = 0.2 - Intensity = IX Weight = 0.2
Updated weight
Fig. 10. The full logic tree developed to update the weight of the 42 hazard curves based on both instrumental and historical seismicity (see dedicated logic tree in Fig. 5 for historical seismicity).
Updated weight 0.100
Updated weight
0.075
Full logic tree updated weight (likelihood function)
0.050
0.025
0.000 1
259
564
688
5
176
551
573
1551 1573 2551 2573 3551 3573
Hazard curve id Fig. 11. Full logic tree: Instrumental seismicity (PGA = 0.01 g) and historical seismicity (6 PGA-to-Intensity relationship + 4 ranges of Intensity) – Updated weight.
5. Incorporating instrumental seismicity
5.2. Implementation
5.1. Principle
The process followed to incorporate instrumental seismicity is the one that was previously developed and published (Humbert et al., 2008b; Humbert, 2015). The updating process is the same than for historical seismicity (based on Eqs. (7)–(9)). The PGA level selected for the Bayesian updating is equal to 0.01 g (Humbert et al., 2008b). The
The principle develop here is to incorporate instrumental seismicity in addition to historical seismicity in the PSHA Bayesian updating process, in order to account for all sources of available observation. 105
Nuclear Engineering and Design 350 (2019) 98–106
E. Viallet, et al.
same stations of observation are also used, cf. Fig. 2. The full logic tree considering instrumental and historical seismicity in described in Fig. 10.
ground motion in the central United States and California. BSSA 97 (2), 497–510. Colombo, A., 1985. Bayesian Methods in Reliability Analysis. 8th Conference on Structural Mechanics in Reactor Technology (SMiRT 8), Brussels, Belgium, M21-10. Faenza, L., Michelini, A., 2010. Regression analysis of MCS intensity and ground motion parameters in Italy and its application in shakemap. Geophys. J. Int. 180 (3), 1138–1152. Heinfling, G., Courtois, A., Viallet, E., 2005. Reliability Based Approach for the Ageing Management of Prestressed Concrete Containment Vessel. 18th Conference on Structural Mechanics in Reactor Technology (SMiRT 18), Beijing, China, H06-5. Humbert, N., Viallet, E., 2008a. An evaluation of epistemic and random uncertainties included in attenuation relationship parameters. 14th World Conference on Earthquake Engineering October 12–17, 2008, Beijing, China. Humbert, N., Viallet, E., 2008. A method for comparison of recent PSHA on the French territory with experimental feedback. 14th World Conference on Earthquake Engineering October 12–17, 2008, Beijing, China. Humbert, N., 2015. PSHA updating technique with a Bayesian framework: innovations. OECD/NEA/CSNI Workshop on Testing PSHA Results and Benefit of Bayesian Techniques for Seismic Hazard Assessment, 4-6 February 2015. Eucentre Foundation, Pavia, Italy. ISBA: International Society for Bayesian Analysis, https://bayesian.org/. Keller, M., Pasanisi, A., Marcilhac, M., Yalamas, T., Secanell, R., Senfaute, G., 2014. A Bayesian methodology applied to the estimation of earthquake recurrence parameters for seismic hazard assessment. Qual. Reliab. Eng. Int. 30 (7). Madsen, H.O., 1985. Model updating in First-Order Reliability Theory with application to fatigue crack growth. 2nd International Workshop on Stochastic Methods in Structural Mechanics, Pavia, Italy. Marin, S., Avouac, J.P., Nicolas, M., Schlupp, A., 2004. A probabilistic approach to seismic hazard in metropolitan France. BSSA 94 (6), 2137–2163. OECD/NEA, 2008. Recent findings and developments in probabilistic seismic hazards analysis methodologies and applications. OECD/NEA/CSNI Workshop Proceedings, Lyon, France 7–9 April 2008, https://www.oecd-nea.org/globalsearch/download.php?doc=8592. OECD/NEA, 2015. Workshop on Testing Probabilistic Seismic Hazard Analysis Results and the Benefits of Bayesian Techniques. OECD/NEA/CSNI Workshop proceedings, Pavia, Italy 4–6 February 2015, https://www.oecd-nea.org/nsd/docs/2015/csnir2015-15.pdf. Ordaz, M., Singh, S., Arciniega, A., 2007. Bayesian Attenuation regression: an application to Mexico City. Runge, A., Händel, A., Riggelsen, C., Scherbaum, F., Kühn, N., 2013. A smart elicitation technique for informative priors in ground-motion mixture modeling. EGU General Assembly 2013, held 7–12 April, 2013 in Vienna, Austria, id. EGU2013-990. Secanell, R., Martin, C., Viallet, E., Senfaute, G., 2015. Bayesian methodology to update the Probabilistic Seismic Hazard Assessment. OECD/NEA/CSNI Workshop on Testing PSHA Results and Benefit of Bayesian Techniques for Seismic Hazard Assessment, 4–6 February 2015. EUcentre Foundation, Pavia, Italy. Secanell, R., Martin, C., Viallet, E., Senfaute, G., 2017. A Bayesian methodology to update the probabilistic seismic hazard assessment. Bull. Earthq. Eng doi: 10.1007/s10518017-0137-3 SisFrance: http://www.sisfrance.net/. Stirling, M., Gerstenberger, M., 2010. Ground motion-based testing of seismic hazard models in New Zealand. BSSA 100 (4), 1407–1414. Tasan, H., Beauval, C., Helmstetter, A., Sandikkaya, A., Guéguen, Ph., 2014. Testing Probabilistic Seismic Hazard estimates against accelerometric data in two countries: France and Turkey. Geophys. J. Int. 198, 1554–1571. Viallet, E., Humbert, N., Martin, C., Secanell, R., 2008. A Bayesian methodology to update the Probabilistic Seismic Hazard Assessment. 14th World Conference on Earthquake Engineering October 12-17, 2008, Beijing, China. Viallet, E., Humbert, N., Mottier, P., 2017. Updating of a PSHA Based on Bayesian Inference with Historical Macroseismic Intensities. 17WCEE, Santiago, Chile. Wang, M., Takada, T., 2009. A Bayesian framework for prediction of seismic ground motion. BSSA 99 (4), 2348–2364. Wald, D.J., Quitoriano, V., Heaton, T.H., Kanamori, H., 1999. Relationships between peak ground acceleration, peak ground velocity and modified mercalli intensity in California. Earthquake Spectra 15, 557–564. Worden, C.B., Gerstenberger, M.C., Rhoades, D.A., Wald, D.J., 2011. Probabilistic relationships between ground-motion parameters and modified mercalli intensity in California. BSSA 102 (1), 204–221.
5.3. Results The final updated weights of the 42 hazard curves based accounting for instrumental and historical seismicity is presented in Fig. 11. 5.4. Discussion We can observe that incorporating instrumental seismicity in the full Bayesian updating logic tree has a slight impact on the results, but does no transform drastically the trends. This instrumental seismicity is then a useful complementary and independent source of observation. Through this process, the update weight is calculated by covering a large range of PGA (from 0.01 g to ∼0.5 g). In addition, the cumulated year of observation among the 19 sites (and through instrumental and historical seismicity) reaches around 5800 years, which is a significant period of time. 6. Conclusion and perspectives As recommended by the OECD/NEA/CSNI (OECD/NEA, 2015), a state-of-the-art PSHA should include a testing phase against any available local observation (including any kind of observation and any period of observation) and should include testing not only against its median results but also against its whole distribution (percentiles). This study, by integrating historical and instrumental observed seismicity through a Bayesian updating technique with the proper propagation of all sources of uncertainties, and performed on every branch of the prior PSHA, allows to fulfil this requirement. In our case, this process of developing the Bayesian inference by covering a large range of PGA (from 0.01 g to ∼0.5 g) and consequently a long cumulated time of observation (around 5800 years) confers a strong confidence to the updated PSHA. Such approaches are now fully matured and should be systematically used as a final step of any PSHA in order to reduce epistemic uncertainties that could not be addressed through experts’ judgments. In the next stages, we will implement this Bayesian updating process on a full scale PSHA develop on the French metropolitan territory, using both instrumental and historical observations. Finally, we strongly recommend that such PSHA testing and Bayesian updating should be included as a final step to SSHAC procedure. References Atkinson, G., Sonley, E., 2000. Empirical relationships between modified mercalli intensity and response spectra. BSSA 90 (2), 537–544. Atkinson, G., Kaka, S., 2006. Relationships between felt intensity and instrumental ground motion for the new madrid shakemaps. Atkinson, G., Kaka, S., 2007. Relationships between felt intensity and instrumental
106
Contents lists available at ScienceDirect
Nuclear Engineering and Design journal homepage: www.elsevier.com/locate/nucengdes
Updating a probabilistic seismic hazard assessment with instrumental and historical observations based on a Bayesian inference
T
E. Vialleta, , N. Humbertb, P. Mottierc ⁎
a
EDF, Nuclear New Build Engineering Division, France EDF, Hydropower Engineering Division, France c Polytechnique Montreal, Canada b
ARTICLE INFO
ABSTRACT
Keywords: Probabilistic seismic hazard assessment PSHA Bayesian inference Bayesian updating
Since the basic work of Cornell, many studies have been conducted in order to evaluate the probabilistic seismic hazard (PSHA) at a given site or at a regional scale. In general, results of such studies are used as inputs for regulatory hazard maps or for risk assessments. Such approaches are nowadays considered as well established and come more and more used worldwide, in addition to deterministic approaches. Nevertheless, some discrepancies have been observed recently in some PSHA, especially from studies conducted in areas with low to moderate seismicity. The lessons learnt from these results lead to conclude that, due to epistemic uncertainties inherent to such a domain, some deterministic choices have to be made and, depending on expert judgments, may lead to strong differences in terms of seismic motion evaluation. In this context, the objective of this paper is to present a methodology that can be used to take into consideration instrumental and historical observations in order to reduce epistemic uncertainties in a PSHA. The method developed here is based on a Bayesian inference technique used in order to quantify the likelihood of the prior estimation, and finally update the PSHA. The period of observation under consideration is the completeness period for each set of observation used. The updating process is developed at a regional scale, over a significant number of stations. The potential correlation between points of observation is also discussed and accounted for. Finally, a case of application is proposed on the French metropolitan territory to demonstrate the efficiency of this updating method and draw perspectives for further applications.
1. Introduction and background
consideration, cf. Eq. (1).
Bayes theory on conditional probabilities was developed in the second half of the 18th century by Reverend Thomas Bayes and was strongly promoted by Pierre-Simon de Laplace in the early 19th century. The use of Bayesian theorem of conditional probability (also known as Bayesian inference), which has drastically increased during the second half of the 20th century, is nowadays largely disseminated in a wide range of fields and applications such as biology, pharmaceutical statistics, economics, finance, environmental and earth science, industrial statistics, mechanics (ISBA, 2018) … The principle of Bayes theorem of conditional probabilities is simple to present (but not always simple to solve!). Based on a probabilistic assessment of a given parameter A (usually called as the prior estimation), Bayesian theorem consists in calculating the likelihood of the prior estimation using actual observations B of the parameter under
P (A |B ) =
⁎
P (B |A). P (A) P (B )
(1)
where P(A) is the prior estimation based on a probabilistic assessment P(B|A) is conditional probability of the observed event B according to the prior assessment P(B) is the total probability of the observed event B according to the prior assessment P(A|B) is the expected updated (posterior) estimation
Corresponding author at: EDF, Nuclear New Build Engineering Division, Direction Technique, 19 rue Pierre Bourdeix, 69007 Lyon, France. E-mail address: [email protected] (E. Viallet).
https://doi.org/10.1016/j.nucengdes.2019.04.034 Received 29 December 2018; Received in revised form 22 April 2019; Accepted 24 April 2019 Available online 20 May 2019 0029-5493/ © 2019 Elsevier B.V. All rights reserved.
Nuclear Engineering and Design 350 (2019) 98–106
E. Viallet, et al.
2. Application of Bayes theorem in the field of nuclear industry and earthquake engineering
3.2. Prior estimation: PSHA complemented by PGA-to-Intensity relationships
In the field of nuclear engineering, Bayesian approaches have been used for a long time, especially in the field of reliability analyses (Colombo, 1985) and mechanical and structural engineering (Madsen, 1985; Heinfling et al., 2005). In the field of earthquake engineering, Bayesian approaches are also largely used, in various aspects, such as seismology and structural dynamics. In seismology, the Bayesian inference is used in elicitation techniques (Runge et al., 2013), for comparison of PSHA hypothesis (Humbert et al., 2008b; Secanell et al., 2017) or for direct updating of GMPEs on local data (Ordaz et al., 2007; Wang et al., 2009; Stirling et al., 2010). In the last decade, many developments and applications were published in the field of PSHA testing and Bayesian updating, especially in the framework of OECD/NEA (OECD/NEA, 2008) and SIGMA Project (Keller et al., 2014; Tasan et al., 2014). PSHA testing and Bayesian updating are now considered as state of the art and should be systematically implemented. This process is not antagonist with expert judgments, which are of course necessary to guaranty the robustness of the prior PSHA and the correct propagation of random and epistemic uncertainties. PSHA testing and Bayesian updating have to be considered as an additional step that allows actual observations at a given site to be taken into consideration to reduce epistemic uncertainties (i.e. uncertainties that could not be reduced by experts), taking into account random and epistemic uncertainties they intrinsically incorporate (OECD/NEA, 2015).
Theoretically, any kind of observation and any period of observation, including instrumental seismicity, historical seismicity and paleoseismicity data can be used in such a process. In this section, we’ll focus on historical seismicity. The objective is to determine an exceedance rate of a given macroseismic Intensity based on the hazard curve under consideration. The principle of the method was described previously (Secanell et al., 2017) and extended in this paper as described below. In order to use historical observations, PGA-to-Intensity relationships are used to transform PGA to macroseismic intensities, based on Eq. (2).
I = a . Log (PGA) + b ( + /
)
(2)
Random and epistemic uncertainties are taken into consideration through the full propagation of the standard deviation of the PGA-toIntensity relationship and through the use of different PGA-to-Intensity relationships. The goal is to determine the rate of exceedance of a given intensity from a given hazard curve expressed in term of PGA. The calculation is performed based on Eq. (3). Amax
NIi =
(NA +
A
NA ). PIi | A
Amin
(3)
where
NIi is the annual rate of exceedance of a given Intensity Ii NA is the annual rate of exceedance of a given PGA A PIi | A is the probability of a given PGA to produce the given intensity Ii according to Eq. (2), including σ (i.e. percentage of the contribution of the given PGA level to the considered class of intensity) ΔA is a discretization step small enough in order to get stable results of NIi (in this study, the range of PGA from 1 to 1000 cm/s2 is discretized into 40 steps)
3. The innovating method developed in this paper The method developed in this paper is based on the Bayesian updating process described in previous papers (Viallet et al., 2008; Humbert et al., 2008b). In these previous references, instrumental seismicity is taken into consideration as observations to update the prior estimation. In addition, a previous work (Viallet et al., 2017) demonstrated that historical seismicity can also be taken into consideration through the same Bayesian process. In the current paper, previous works are extended and the Bayesian inference is applied based on both historical and instrumental seismicity, in order to update the prior PSHA. The objective is thus to demonstrate that such an approach can be applied combining multiple sources of observation (instrumental and historical seismicity in this case) and to show that the proposed process has an improving impact on the weighting process of the prior PSHA branches.
Finally, the annual exceedance rate is multiplied by the completeness period of the historical catalogue in order to determine the prior rate of exceedance of a given hazard curve in terms of macroseismic Intensity. The illustration of this process is given in Fig. 1. This allows to express the prior PSHA estimation in term of rate of exceedance of a given Macroseismic Intensity. Through this process, the full range of uncertainties is propagated in order to guaranty that none of the possibilities are disregarded.
3.1. Principle for incorporating historical observations
3.3. Historical observations consideration
First of all, the approach needs a PSHA following the state of the art approach as input data. This prior estimation has to take into consideration all sources of random and epistemic uncertainties. This PSHA has to be carried out though a well-documented procedure, which should explicitly describe the elicitation of experts’ judgments, such as the so-called SSHAC procedure. Hazard curves expressed in term of PGA (Peak Ground Acceleration) rate of exceedance will be considered as the prior estimation for the implementation of Bayes theorem of conditional probability. The hazard curve is used to determine an exceedance rate of a given macroseismic Intensity based on Intensity to PGA relationships. Random and epistemic uncertainties are fully propagated, as described in the following sections. Historical observations are expressed in term of macroseismic Intensity observed at the site under consideration. Finally, the Bayesian theorem of conditional probability is applied for each single hazard curve of the prior estimation in order to calculate its likelihood and though to update its weight accordingly. The likelihood function becomes the posterior weight of the hazard curves.
The historical observations have to be well documented. This means that a historical database has to be available and its quality has to be high enough to allow to define completeness period and to be sure that it is possible to quantify the number of event that were felt on a given location over this completeness period. Once the previous requirements are fulfilled, the historical observations are the number of historical events that produced a Macroseismic Intensity higher or equal to a given one. An application is presented in the next section. We would like to mention that the historical observation which are used here are also usually used in the PSHA as input data (for earthquake recurrence parameters’ determination for instance). This means that these data could be used twice in the process, which may be confusing. However, we are dealing here with a reduced number of historical observations compared with the huge number that are considered in the prior PSHA. Then, it can be easily shown that removing these observations from the set of data used in the prior PSHA would have no impact on the PSHA result. Then, there is no risk of giving them an arbitrary high weight in the overall process (by using them twice). 99
Nuclear Engineering and Design 350 (2019) 98–106
PGA to Intensity relationship : Illustration
1.E+01 1.E+00
Annual rate of exceedance
1.E-01 1.E-02 1.E-03 1.E-04 1.E-05 1.E-06
Rate of exceedance (over the completeness period : XXX y)
E. Viallet, et al.
35
30
N Ii =
25
Amax
∑ (N
A + ∆A
VII-VIII / 265y
VIII / 515y
Amin
− N A ).PI i
A
20
15
10
5
0 V / 135y
V-VI / 135y
VI / 165y
VI-VII / 165y
VII / 265y
VIII-IX / 515y
IX / 715y
Intensity range (above) + Illustration of PGA contribution (below)
1.E-07 10
100
1000
PGA (cm/s²)
1.0
1.2
1.4
1.6
1.9
2.3
2.6
3.1
3.6
4.3
5.0
5.8
6.9
8.0
9.4
11.0
12.9
15.2
17.8
20.8
24.4
28.6
33.5
39.3
46.1
54.0
63.3
74.2
87.0
101.9
Fig. 1. Process used to transform hazard curve in PGA into macroseismic Intensity. Left: Discretized hazard curve – Middle: PGA-to-Intensity relationship – Right: Intensity rate of exceedance.
3.4. The updating process: accounting for random occurrence of events
multiple sites of observation, the Poisson’s occurrence model Eq. (6) is generalized through a negative binomial distribution that allows to account for correlation between stations of observation, as described in detail in a previous paper (Humbert, 2015) and shown in Eq. (7).
The updating process (Eq. (1)) applied to our case is shown in Eq. (4).
P (Model|Observation) =
P (Observation|Model) . P (Model) P (Observation)
P (Observation|Ai ) = P (n, t ) = exp ln
(4)
where
ln (1 + n) .
P(Model) is the prior estimation, direct output of the PSHA based on a probabilistic assessment P(Observation|Model) is the probability of the observed event according to the prior assessment P(Observation) is the total probability of the observed event according to the prior assessment P(Model|Observation) is the expected updated (posterior) estimation
P (Observation|Ai ) . P (Ai ) P (Observation)
(Ai ). n!
(A ) 1
n
(7)
P (Observation|Ai )
(8)
3.4.3. Likelihood function and posterior estimation Finally, the likelihood function, Eq. (9), becomes the updated weight of the hazard curve Ai under consideration, which allows to recalculate the range of (reduced) epistemic uncertainties of the PSHA, and becomes the posterior estimation.
P (Observation|Ai ) P (Observation)
3.4.1. Calculation of P(Observation|Ai) Based on the rate of exceedance calculated as described in Eq. (3), P (Observation|Ai) is calculated assuming a Poisson’s occurrence model, Eq. (6).
e
1 K
N
P (Ai ) is the prior weight of the hazard curve under consideration, which results directly from the weighting process of the PSHA under consideration, including expert judgment and epistemic uncertainties propagation P (Observation|Ai ) is the probability of the observed event according to the prior assessment P (Observation) is the total probability of the observed event B according to the prior assessment
P (Observation|Ai ) = P (n, t ) =
1
K
3.4.2. Calculation of P(Observation) P(Observation) is simply calculated through Equation (8).
(5)
P (Observation) =
t )n
1
ln
lnΓ is the logarithmic gamma function λ(A) is the cumulative rate of exceedance of a given acceleration A for all the stations under consideration n is the number of events with an acceleration higher that A K is the average number of sites impacted by one earthquake (this parameter is indicative of the correlation among stations. If K tends to 1, the negative binomial distribution goes towards a Poisson distribution).
where
.t (
A
(A ) +n 1
where
This equation is applied to each single hazard curve of the prior PSHA (Ai) so that the updating process will allows to calculate a likelihood of each single hazard curve of that prior estimation, Eq. (5).
P (Ai |Observation) =
1K K
K
(9)
4. Application: French metropolitan territory 4.1. Description of the PSHA used for this application
(6)
In order to illustrate the efficiency of this method, the application is performed based on the French Metropolitan territory. The prior PSHA is based on a previous study (Humbert et al., 2008a). This study is not described in detail in this paper in order to concentrate on the Bayesian updating process. However, this PSHA includes the propagation of a large range of random and epistemic uncertainties, especially on the GMPEs and the random model as shown by Eq. (10).
where λ is the annual rate of exceedance of a given acceleration Ai t is the period of time under consideration n is the number of event with an acceleration higher that Ai In order to perform this updating process at a regional scale using 100
Nuclear Engineering and Design 350 (2019) 98–106
E. Viallet, et al.
Fig. 2. Location of the 19 EDF NPPs on the French Metropolitan territory (green diamonds). NB: Yellow triangles are French accelerometric networks stations not used in this study. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)
( ) (±a ) + (±a ). M . a a R + ( ± ) + (± ) + (± )
( )
a1 a2 Log (PGA) = ± .M+ ± . M2 + 1 2 Log
2
2
5
5
3
relationships from different authors were used. The only criterion that was used to select PGA-to-Intensity relationships is the domain of validity in term of PGA that should be large enough to allow to estimate Intensities without any bias (especially boundary effects). These selected PGA-to-Intensity relationships are presented in Table 1 and plotted in Fig. 4.
4
3
4
6
6
(10)
where PGA is the Peak Ground Acceleration, M is the magnitude, R is the hypocentral distance, Every parameter (a1 to a6) has its own uncertainty, even the random part σ.
4.3. Completeness periods for historical seismicity Based on SisFrance historical seismicity database (SisFrance, 2018), the period of completeness of historical seismicity is known for the French Metropolitan territory, as shown in Table 2. This assumption has to be consistent with the one that is made for the prior PSHA, in term of period and uncertainties, if any. This consistency was confirmed in our case. In addition, it should be pointed out that the completeness of observation data is crucial and should be checked and confirmed in order to guaranty that the updated process does not introduce any bias (if observation database is incomplete, the Bayesian approach will artificially increase the weight of branches that are consistent with this biased observation set).
This process allows to obtain hazard curves for a large number of locations. In our case, we are using the results obtained for 19 EDF Nuclear Power Plant locations, in order to use historical seismicity as observation, which is well documented for these sites. The location of EDF NPPs is shown in Fig. 2. For this paper, 42 hazard curves have been generated among the whole logic tree that contains much more than this number. This choice (only 42 hazard curves) was made in order to ease the interpretation and to better illustrate the updating process. In this purpose, the 42 selected hazard curves were randomly chosen among the whole logic tree of the prior PSHA. Fig. 3 shows the 42 hazard curves obtained at one site.
4.4. Observations: historical earthquakes felt on sites Based on SisFrance database and Seismic Hazard Assessment of the EDF’s NPP, it is possible to quantify the number of historical events felt on the 19 EDF sites. This number of observed events, quantified per range of magnitude, is presented in Table 3.
4.2. PGA-to-Intensity transformation In order to propagate a large range of uncertainties, different 101
Nuclear Engineering and Design 350 (2019) 98–106
E. Viallet, et al.
1.E+01
1 592 466 564 626 688 5 143 189 551 557 593 1551 1557 1593 2551 2557 2593 3551 3557 3592
1.E+00
Annual rate of exceedance
1.E-01 1.E-02 1.E-03 1.E-04 1.E-05 1.E-06
62 746 139 363 259 781 30 176 550 554 573 1550 1554 1573 2550 2554 2573 3550 3554 3573 3593
1.E-07 10
100
1000
PGA (cm/s²) Fig. 3. Illustration of the 42 hazard curves resulting from the logic tree under consideration at a given site.
4.5. Logic tree developed for Bayesian updating based on historical observation
is quite comparable for all hazard curves. This illustrates that for such high Intensities, the updating process is not able to rank the hazard curves, simply because at this level of Intensity (which corresponds to long return periods), all hazard curves do not predict any observation, as actually observed. When there is not enough observation, the updating process does not modify prior weights, as also observed in previous studies (Secanell et al., 2015). Finally, the updated weight of the hazard curves resulting from the historical seismicity logic tree is shown in Fig. 9. This figure clearly shows the efficiency of the updating process. The updating process, which does not modify any of the parameters of the prior PSHA, only modifies the weight of the hazard curves. One can finally observe that none of the prior hazard curves are completely eliminated. Their weights are reduced according to the likelihood that is calculated through Bayes theorem of conditional probabilities, according to observations.
The logic tree developed to update the weight of the 42 hazard curves is described in Fig. 5. This process allows to account for different PGA-to-Intensity relationships and to account for different Intensity ranges (i.e. different return periods). 4.6. Results The updated weight (Bayesian likelihood function) of the hazard curves are presented in the next figures. The likelihood function obtained for each PGA-to-Intensity relationship is shown in Fig. 6. In this figure, the height of the colored bar corresponds to the value of the likelihood of the hazard curve under consideration, for the PGA-to-Intensity relationship under consideration. The likelihood function obtained for each Intensity range included in the logic tree is shown in Fig. 7. In this figure, the height of the coloured bar corresponds to the value of the likelihood of the hazard curve under consideration, for the Intensity range under consideration. Fig. 8 “radar view” illustrates how the likelihood evolves with Intensity. This figure shows that for high Intensity (I = IX), the likelihood
4.7. Discussion First of all, we can conclude that the Bayesian updating procedure is successfully implemented and allows to update a PSHA based on actual historical seismicity, as observed on the sites under consideration. We also observe that if observations are not sufficient, the updating process
Table 1 PGA-to-Intensity relationships used in this study. Author
Domain of application
Relationship (PGA in cm.s−2)
σ
Atkinson and Kaka (2006)
Missouri et California I ∊ [II;IX] PGA ∊ [0.4;260] Centre US et California I ∊ [II;IX]
I = 2.315 + 1.319*log(PGA) + 0.372*log(PGA)2
0.93
If log(PGA) ≤ 1.69, I = 2.65 + 1.39*log(PGA)
1.01
If log(PGA) > 1.69, I = −1.91 + 4.09*log(PGA) I = 2.58*Log(PGA) + 1.68
0.35
I = 2.3*Log(PGA(g)) + 10 If log(PGA) ≤ 1.57, I = 1.78 + 1.55*log(PGA)
0.3 0.73
If log(PGA) > 1.57, I = −1.60 + 3.70*log(PGA) If I ≤ V, I = 7.58 + 2.2*log(PGA(g))
1.7
Atkinson and Kaka (2007) Faenza and Michelini (2010) Marin et al. (2004) Worden et al. (2011) Wald et al. (1999) and Atkinson and Sonley (2000)
Italy I ∊ [I;VIII] France California I ∊ [II;VIII-IX] PGA ∊ [0.07;800] California-France
If I > V, I = 10.18 + 4.35*log(PGA(g))
102
Nuclear Engineering and Design 350 (2019) 98–106
E. Viallet, et al.
10 9 Atkinson et Kaka (2007)
Macroseismic Intensity
8 7
Marin & All (2004)
6 5
Faenza & Michelini (2010)
4
Worden & All (2011)
3
Wald, Atkinson&Sonley (BRGM-2000)
2
Atkinson et Kaka (2006)
1 0 1
10
100
1000
PGA (cm/s²) Fig. 4. PGA-to-Intensity relationships used in this study.
Table 2 Completeness periods for historical seismicity. Zone
Intensities IV to IV-V
Intensities V to V-VI
Intensities VI to VI-VII
Intensities VII to VII-VIII
Intensities VIII to VIII-IX
Intensities IX to IX-X
Entire France Alpes Britany North Provence Pyrenees Rhine Other areas
1920 1950 1920 1920 1925 1930 1885 1960
1880 1880 1880 1880 1880 1910 1870 1850
1850 1830 1750 1850 1850 1910 1750 1750
1750 1800 1700 1750 1750 1750 1630 1750
1500 1500 1500 1500 1500 1500 1500 1500
1300 1300 1300 1300 1300 1300 1300 1300
Table 3 Total amount of historical events felt on all EDF NPP sites. Intensity range Completeness period
V 135y
V-VI 135y
VI 165y
VI-VII 165y
VII 265y
VII-VIII 265y
VIII 515y
VIII-IX 515y
IX 715y
IX-X 715y
X 715y
Total number of events (whole 19 EDF’s sites)
24
12
3
2
0
0
0
0
0
0
0
PGA to Intensity relationships Weight Atkinson&Kaka (2006) Weight = 1/6 Atkinson&Kaka (2007) Weight = 1/6
Bayesian updating Logic tree => Implemented for each hazard curve of the PSHA
Faenza & Michelini, 2010 Weight = 1/6 Marin, Avouac, Schlupp, & Nicolas, 2004 Weight = 1/6 Worden, Gerstenberger, Rhoades, & Wald, 2011 Weight = 1/6 Wald et Atkinson & Sonley (BRGM-2000) Weight = 1/6
Intensity Weight
I = VI Weight = 1/4 I = VII Weight = 1/4
Updated weight
I = VIII Weight = 1/4 I = IX Weight = 1/4
Fig. 5. The logic tree developed to update the weight of the 42 hazard curves based on historical seismicity.
103
Nuclear Engineering and Design 350 (2019) 98–106
E. Viallet, et al.
Baysian updating - Historical seismicity Effect of Intensity - PGA relationship Atkinson et Kaka (2006)
0.6
Updated weight
0.5
Wald, Atkinson&Sonley (BRGM-2000) Worden & All (2011)
0.4 0.3
Faenza & Michelini (2010) Marin & All (2004)
0.2 0.1 0 1
259 564 688
5
Atkinson et Kaka (2007)
176 551 573 1551 1573 2551 2573 3551 3573
Hazard curve id Fig. 6. Likelihood function depending on PGA-to-Intensity relationships.
Baysian updating - Historical seismicity Effect of Intensity range 0.4
Updated weight
0.35 0.3 0.25
Weights (I=IX)
0.2
Weights (I=VIII)
0.15
Weights (I=VII)
0.1
Weights (I=VI)
0.05 0 1
259
564
688
5
176
551
573 1551 1573 2551 2573 3551 3573
Hazard curve id Fig. 7. Likelihood function depending on range of Macroseismic Intensity used.
35923593 3573 3557 3554 3551 3550
1
Baysian updating - Historical seismicity
62 139 259 363 466 564 592
2593
Effect of Intensity range
626
2573
688
2557
746
2554
781
2551
5
2550 1593 1573 1557 1554 1551 1550593
30
573
143 176 189 550 551 557 554
Weights (I=VI) Weights (I=VII) Weights (I=VIII) Weights (I=IX)
Fig. 8. Likelihood function depending on range of Macroseismic Intensity used – Radar view.
does not affect the prior estimation (typically the case for Intensity IX in our study). This means that there is no risk of performing updating: the process is self-reliable. This also means that even with few (even no) observation, Bayesian updating can be implemented. In addition, it is important to point out that the updating process does not modify any of the assumptions or input data of the prior estimation. It only updates the weight of the hazard curves (output of the PSHA), based on observations. This process allows to reduce the
epistemic uncertainties that are contained in the prior model, based on their intrinsic (and objective) likelihood. Finally, the procedure developed in this paper, based on historical seismicity, takes into consideration all sources of uncertainties (random and epistemic) in the PGA-to-Intensity transformation process and in the random occurrence of events, which avoid any deterministic or inappropriate choices…
104
Nuclear Engineering and Design 350 (2019) 98–106
E. Viallet, et al.
Updated weight 0.1
Updated weight
0.075
Historical seismicity logic tree
0.05
0.025
0 1
259
564
688
5
176
551
573 1551 1573 2551 2573 3551 3573
Hazard curve id Fig. 9. Historical seismicity logic tree (6 PGA-to-Intensity relationship + 4 ranges of Intensity): Updated weight.
Instrumental seismicity - PGA threshold = 0.01g – Weight = 0.2 Bayesian updating Logic tree => Implemented for each hazard curve of the PSHA
Historical seismicity (see dedicated logic tree in fig. 5) Total weight = 0.8 - Intensity = VI Weight = 0.2 - Intensity = VII Weight = 0.2 - Intensity = VIII Weight = 0.2 - Intensity = IX Weight = 0.2
Updated weight
Fig. 10. The full logic tree developed to update the weight of the 42 hazard curves based on both instrumental and historical seismicity (see dedicated logic tree in Fig. 5 for historical seismicity).
Updated weight 0.100
Updated weight
0.075
Full logic tree updated weight (likelihood function)
0.050
0.025
0.000 1
259
564
688
5
176
551
573
1551 1573 2551 2573 3551 3573
Hazard curve id Fig. 11. Full logic tree: Instrumental seismicity (PGA = 0.01 g) and historical seismicity (6 PGA-to-Intensity relationship + 4 ranges of Intensity) – Updated weight.
5. Incorporating instrumental seismicity
5.2. Implementation
5.1. Principle
The process followed to incorporate instrumental seismicity is the one that was previously developed and published (Humbert et al., 2008b; Humbert, 2015). The updating process is the same than for historical seismicity (based on Eqs. (7)–(9)). The PGA level selected for the Bayesian updating is equal to 0.01 g (Humbert et al., 2008b). The
The principle develop here is to incorporate instrumental seismicity in addition to historical seismicity in the PSHA Bayesian updating process, in order to account for all sources of available observation. 105
Nuclear Engineering and Design 350 (2019) 98–106
E. Viallet, et al.
same stations of observation are also used, cf. Fig. 2. The full logic tree considering instrumental and historical seismicity in described in Fig. 10.
ground motion in the central United States and California. BSSA 97 (2), 497–510. Colombo, A., 1985. Bayesian Methods in Reliability Analysis. 8th Conference on Structural Mechanics in Reactor Technology (SMiRT 8), Brussels, Belgium, M21-10. Faenza, L., Michelini, A., 2010. Regression analysis of MCS intensity and ground motion parameters in Italy and its application in shakemap. Geophys. J. Int. 180 (3), 1138–1152. Heinfling, G., Courtois, A., Viallet, E., 2005. Reliability Based Approach for the Ageing Management of Prestressed Concrete Containment Vessel. 18th Conference on Structural Mechanics in Reactor Technology (SMiRT 18), Beijing, China, H06-5. Humbert, N., Viallet, E., 2008a. An evaluation of epistemic and random uncertainties included in attenuation relationship parameters. 14th World Conference on Earthquake Engineering October 12–17, 2008, Beijing, China. Humbert, N., Viallet, E., 2008. A method for comparison of recent PSHA on the French territory with experimental feedback. 14th World Conference on Earthquake Engineering October 12–17, 2008, Beijing, China. Humbert, N., 2015. PSHA updating technique with a Bayesian framework: innovations. OECD/NEA/CSNI Workshop on Testing PSHA Results and Benefit of Bayesian Techniques for Seismic Hazard Assessment, 4-6 February 2015. Eucentre Foundation, Pavia, Italy. ISBA: International Society for Bayesian Analysis, https://bayesian.org/. Keller, M., Pasanisi, A., Marcilhac, M., Yalamas, T., Secanell, R., Senfaute, G., 2014. A Bayesian methodology applied to the estimation of earthquake recurrence parameters for seismic hazard assessment. Qual. Reliab. Eng. Int. 30 (7). Madsen, H.O., 1985. Model updating in First-Order Reliability Theory with application to fatigue crack growth. 2nd International Workshop on Stochastic Methods in Structural Mechanics, Pavia, Italy. Marin, S., Avouac, J.P., Nicolas, M., Schlupp, A., 2004. A probabilistic approach to seismic hazard in metropolitan France. BSSA 94 (6), 2137–2163. OECD/NEA, 2008. Recent findings and developments in probabilistic seismic hazards analysis methodologies and applications. OECD/NEA/CSNI Workshop Proceedings, Lyon, France 7–9 April 2008, https://www.oecd-nea.org/globalsearch/download.php?doc=8592. OECD/NEA, 2015. Workshop on Testing Probabilistic Seismic Hazard Analysis Results and the Benefits of Bayesian Techniques. OECD/NEA/CSNI Workshop proceedings, Pavia, Italy 4–6 February 2015, https://www.oecd-nea.org/nsd/docs/2015/csnir2015-15.pdf. Ordaz, M., Singh, S., Arciniega, A., 2007. Bayesian Attenuation regression: an application to Mexico City. Runge, A., Händel, A., Riggelsen, C., Scherbaum, F., Kühn, N., 2013. A smart elicitation technique for informative priors in ground-motion mixture modeling. EGU General Assembly 2013, held 7–12 April, 2013 in Vienna, Austria, id. EGU2013-990. Secanell, R., Martin, C., Viallet, E., Senfaute, G., 2015. Bayesian methodology to update the Probabilistic Seismic Hazard Assessment. OECD/NEA/CSNI Workshop on Testing PSHA Results and Benefit of Bayesian Techniques for Seismic Hazard Assessment, 4–6 February 2015. EUcentre Foundation, Pavia, Italy. Secanell, R., Martin, C., Viallet, E., Senfaute, G., 2017. A Bayesian methodology to update the probabilistic seismic hazard assessment. Bull. Earthq. Eng doi: 10.1007/s10518017-0137-3 SisFrance: http://www.sisfrance.net/. Stirling, M., Gerstenberger, M., 2010. Ground motion-based testing of seismic hazard models in New Zealand. BSSA 100 (4), 1407–1414. Tasan, H., Beauval, C., Helmstetter, A., Sandikkaya, A., Guéguen, Ph., 2014. Testing Probabilistic Seismic Hazard estimates against accelerometric data in two countries: France and Turkey. Geophys. J. Int. 198, 1554–1571. Viallet, E., Humbert, N., Martin, C., Secanell, R., 2008. A Bayesian methodology to update the Probabilistic Seismic Hazard Assessment. 14th World Conference on Earthquake Engineering October 12-17, 2008, Beijing, China. Viallet, E., Humbert, N., Mottier, P., 2017. Updating of a PSHA Based on Bayesian Inference with Historical Macroseismic Intensities. 17WCEE, Santiago, Chile. Wang, M., Takada, T., 2009. A Bayesian framework for prediction of seismic ground motion. BSSA 99 (4), 2348–2364. Wald, D.J., Quitoriano, V., Heaton, T.H., Kanamori, H., 1999. Relationships between peak ground acceleration, peak ground velocity and modified mercalli intensity in California. Earthquake Spectra 15, 557–564. Worden, C.B., Gerstenberger, M.C., Rhoades, D.A., Wald, D.J., 2011. Probabilistic relationships between ground-motion parameters and modified mercalli intensity in California. BSSA 102 (1), 204–221.
5.3. Results The final updated weights of the 42 hazard curves based accounting for instrumental and historical seismicity is presented in Fig. 11. 5.4. Discussion We can observe that incorporating instrumental seismicity in the full Bayesian updating logic tree has a slight impact on the results, but does no transform drastically the trends. This instrumental seismicity is then a useful complementary and independent source of observation. Through this process, the update weight is calculated by covering a large range of PGA (from 0.01 g to ∼0.5 g). In addition, the cumulated year of observation among the 19 sites (and through instrumental and historical seismicity) reaches around 5800 years, which is a significant period of time. 6. Conclusion and perspectives As recommended by the OECD/NEA/CSNI (OECD/NEA, 2015), a state-of-the-art PSHA should include a testing phase against any available local observation (including any kind of observation and any period of observation) and should include testing not only against its median results but also against its whole distribution (percentiles). This study, by integrating historical and instrumental observed seismicity through a Bayesian updating technique with the proper propagation of all sources of uncertainties, and performed on every branch of the prior PSHA, allows to fulfil this requirement. In our case, this process of developing the Bayesian inference by covering a large range of PGA (from 0.01 g to ∼0.5 g) and consequently a long cumulated time of observation (around 5800 years) confers a strong confidence to the updated PSHA. Such approaches are now fully matured and should be systematically used as a final step of any PSHA in order to reduce epistemic uncertainties that could not be addressed through experts’ judgments. In the next stages, we will implement this Bayesian updating process on a full scale PSHA develop on the French metropolitan territory, using both instrumental and historical observations. Finally, we strongly recommend that such PSHA testing and Bayesian updating should be included as a final step to SSHAC procedure. References Atkinson, G., Sonley, E., 2000. Empirical relationships between modified mercalli intensity and response spectra. BSSA 90 (2), 537–544. Atkinson, G., Kaka, S., 2006. Relationships between felt intensity and instrumental ground motion for the new madrid shakemaps. Atkinson, G., Kaka, S., 2007. Relationships between felt intensity and instrumental
106