What is the impact of the August 24 , 2016 Amatrice earthquake on the seismic hazard assessment in central Italy ?

The recent Amatrice strong event (Mw6.0) occurred on August 24, 2016 in Central Apennines (Italy) in a seismic gap zone, motivated us to study and provide better understanding of the seismic hazard assessment in the macro area defined as “Central Italy”. The area affected by the sequence is placed between the Mw6.


I. INTRODUCTION
Central Apennines are one of the most seismically active areas in Italy, where a long history of earthquakes has strongly influenced the development of earthquake-resistant structural design.Seismicity is mostly characterized by the occurrence of earthquakes which are distributed within a 50-60 km wide mountain belt extending, in NNW-SSE Apennines direction, from Tuscany to Abruzzi.During the last century, the region has frequently experienced strong shaking caused by earthquakes.The updated historical/parametric catalog [CPTI15, Rovida et al., 2016] contains in the macro area of Central Italy 69 events with magnitude 5.0 or greater, and 7 with magnitude exceeding 6.0, considering a depth 30km.The most relevant of them are those of the 1703 Val Nerina seismic sequence (M w 6.92,I 0 XI MCS) that struck the whole central Italy killing roughly ten thousand people, and the 1915 Marsica earthquake (M w 7.08, I 0 XI MCS) about 60 km south of the latest site that killed around 32000 people.More recently, in 1997 the 6.0 magnitude Colfiorito earthquake, 50 km to the north, killed 11 people in the Umbria and Marche regions.The recent August 24, 2016 Amatrice earthquake (M w 6.0, 01:36 UTC), which had an epicenter roughly 10 km southeast of Norcia, occurred just over seven years after the April 6, 2009 L'Aquila earthquake (M w 6.3), killing more than 300 people only about 50 km away.The Amatrice area was severely damaged also in 1639 (M w 6.2, I 0 IX-X MCS) and successively in 1672 (M w 5.3, I 0 VII-VIII MCS).This part of the central Apennines has been classified as a first category seismic zone in the current Italian probabilistic seismic hazard map (MPS04) [Gruppo di lavoro MPS, 2004].Because the Apennines earthquakes are relatively shallow and structurally complex and many of the local towns and cities contain vulnerable buildings, the strong shaking from these earthquakes has the potential to cause major damage and loss of life in urban areas.The implementation of building codes mandating the use of earthquake resistant buildings has been highly useful and successful in mitigating the impact of earthquakes.In fact, the probabilistic seismic hazard (PSHA) maps have been widely used in building codes to develop design criteria for buildings and for policy makers to implement priorities for earthquake risk reduction.These maps are derived by estimating a variety of parameters for selected models that are used to forecast future seismicity and the resulting shaking.Probabilistic seismic hazard maps are revised regularly adding newly observed and studied material to improve the hazard assessments in order to reduce earthquake damages in United States, New Zealand and Japan, every 6, 7 and 10 years respectively.In Italy, following the issuing of the Prime Minister Ordinance n.3519 in 2006, the current probabilistic seismic hazard map (MPS04) is released as the official reference for seismic hazard values to be used in engineering applications [Gruppo di Lavoro, MPS04, 2004;Stucchi et al., 2011].The hazard assessment was conducted using the Cornell [1968] probabilistic approach and assuming a Poisson process, where the predicted occurrence of future earthquakes is independent of previous occurrence on the same source.Recently there has been an effort to update the current seismic hazard maps requested in 2015 by the Italian Civil Protection Agency to the Center for Seismic Hazard (CPS) of the Istituto Nazionale di Geofisica e Vulcanologia (INGV).In this study for the seismic hazard calculation in Italy, we utilized the smoothed-historical-seismicity method [Frankel, 1995] using two different catalog periods from 1871 to 2003 and from 1871 up to August 31, 2016 in order to investigate the effect of the recent seismicity on the seismic hazard maps.Results are presented both in terms of the earthquake rates as well as the hazard curves at Amatrice, interested by the 2016 sequence.
(http://iside.rm.ingv.it/iside/standard/index.jsp) Detailed information regarding to data and macro area, considered for the seismic hazard assessment, are reported in the Electronic Supplementary material as files.We used the declustered CPTI catalog and the period of completeness considered for this area as given by the Working Group, MPS16, 2016: the events from 1871 to 2014 with completeness magnitude M w ≥ 4.55, and depth ≤ 30 km.
We updated this catalog adding nine recent events provided from the ISIDE catalog for the period between 01/01/2015 and-08/31/2016, with magnitude M w ≥ 4.55 and depth ≤ 30 km within the study area.The local magnitude M L of the ISIDE catalog was converted into moment magnitude M w , using the equation suggested by Gasperini et al. (2013): In order to remove large fluctuations of seismicity rates in space and time due to aftershock sequences we declustered the ISIDE catalogue using the algorithm developed by Gardner and Knopoff [1974].The total number of declustered events was 198.

III. METHODS
To compute the Probabilistic Seismic Hazard Analysis (PSHA) in Central Italy we used the classical Cornell approach [Cornell, 1968], which provides an estimate of ground shaking at a given site for each earthquake magnitude and distance, through an earthquake rate model and a ground motion prediction equation (GMPE).The earthquake rate model considers both the spatial and the magnitude distributions of the seismicity.
For the spatial distribution we use the classical time independent smoothed seismicity approach [Frankel, 1995], with a fine spatial grid of 0.025° both in latitude and longitude; the spatial rate in each cell of the grid is proportional to: (2) where is the total number of events in the catalog, is the distance between the i-th cell and the j-th event, is the smoothing (or correlation) distance.
To obtain the final rate for each cell we normalized the values to preserve the total number of events by using the equation: (3) where is the total number of cells.
The larger is the value of , the more smoothed will be the seismicity, but also the larger will be the contamination between different areas.The result is given as number of events per cell size, per total time spanned by the catalog.For each cell i, the smoothed value can be obtained from Equation (3), and normalized to preserve the total number of events.
The correlation distance,  = 30 km, was determined for "Central Italy" considering the CPTI15 catalog from 1871 to 2014 (M w  4.55 and depth  30km) by maximizing the likelihood of the seismicity contained in half of the catalog, under the time-independent model obtained from the other half of the catalog and vice-versa.
For the magnitude distribution we use the tapered Gutenber-Richter law [Kagan and Jackson, 2000], that has a cumulative distribution function equal to: (4) where is the seismic moment, is the seismic moment of completeness, is the parameter that controls the distribution in the upper ranges of M ('upper corner moment'), and represents the slope of the distribution.Since this distribution holds for seismic moment, we convert all M w magnitude values with the Kanamori [1977] formula for estimating the seismic moment: To estimate , we used the fact that if << [Bird and Kagan, 2004], as in our case, and is the coefficient of the classical Gutenberg-Richter law [Gutenberg and Richter, 1954]; the bvalue was determined by the maximum likelihood method of Aki (1965): where Mw mean is the mean magnitude, and Mw compl is the completeness magnitude of the catalog.
For the corner magnitude (that once converted with the Kanamori formula becomes the corner moment ) we use the maximum observed earthquake in the "Central Italy" zone in the CPTI15 catalog (the Marsica earthquake of January 13, 1915 M w 7.08).
By adopting this simple approach we avoided the estimation of M max , the maximum possible magnitude, that is a difficult parameter to estimate also with a catalog that lasts centuries [Zöller and Holschneider, 2015].
We also attempted to include aftershocks and foreshocks into the probabilistic seismic hazard assessment as the non-inclusion may lead to underestimate the seismic risk.In order to incorporate clustered events into PSHA we used the procedure proposed by Marzocchi and Taroni [2014] for the new Italian seismic hazard map (Internal report for MPS 2016) by applying a corrective factor for the missing rate due to declustering in the catalogs.First, we calculated the rate using the declustered catalog, then we obtained the new clustered seismicity rates multiplying by a factor in a way that the total seismicity rate matches the whole observed seismicity rate.In this way we assumed that each cell of the grid had similar clustering properties.
In the Electronic Supplementary material we put additional information.
We calculated the annual rate of exceeding a specified ground motion at a site from a double summation over distance and magnitude, using both the time-independent model and Ground Motion Prediction Equations (GMPEs).We use two recent GMPEs among those selected during the SHARE EU Project [Delavaud et al., 2012] for active shallow crustal regions to assess ground shaking hazard, without considering the prevalent fault mechanism.The first derives from the Italian GMPE model of Bindi et al. [2011] (hereafter ITA10) and the second one was developed by Akkar and Bommer [2010] (hereafter AB10) for Europe using the Mediterranean, and the Middle East data.In the Electronic Supplementary files, we also put the matrix used following the Marzocchi and Taroni [2014] approach for the rate model and the seismic hazard.

IV. RESULTS
We considered two different time periods of the declustered seismicity to calculate the earthquake rate model in the macro area "Central Italy": 1) 1871-2003 (named hereafter "Model 2003") 2) 1871-August 31, 2016 (named hereafter "Model 2016").So that we verified if the latest strong earthquakes (Table 1) caused any change to the cumulative annual rate in this study area, with respect to the rate of the first period, when the official map of the seismic hazard for the Italian territory (MPS04) was released."Model 2003" contains 189 declustered events with Mw4.55 and depth30km (b-value=0.97±0.07)while "Model 2016" includes the 9 events (Table 1) occurred in the period 2004-2016 (198 total declustered events with a b-value=0.94±0.07)(Mw4.55 and depth30km).The total cumulative annual .rateand b-values obtained for "Model 2003" and "Model 2016", respectively are summarized in Table 2.A comparison of Model 2016 with Model 2003 shows a smaller decrease of the rate, which contributes to a hazard reduction, together with a smaller decrease of the b-value, which, on the other hand, increases the hazard (Table 3).    .The rate in the Amatrice zone is slightly decreased (5%), while in the L'Aquila area the rate is unchanged.We also verified if the seismic hazard in terms of ground acceleration in the Amatrice changed between Model 2003 and Model 2016.To do so, we computed the hazard curves using seismicity models described with the GMPE derived from the Italian strong motion database [Bindi et al., 2011] (ITA10).For the Model 2016 we also computed the hazard curve using the GMPE AB10 [Akkar and Bommer, 2010], applicable in the European, Mediterranean, and Middle Eastern regions ("Model 2016 AB10" in Figure 3).Figure 3 also reports the hazard curve for Amatrice site obtained using the non declustered catalog for Model 2016 adopting the procedure provided by Marzocchi and Taroni (2014) (named as "Model 2016 ITA10 not declustered" in Figure 3).We outline that the Marzocchi and Taroni [2014] approach works with return periods bigger or equal than 475 years, i.e. annual rate of exceedance lesser than 1/475.Figure 3 shows the four hazard curves related to the adopted models, having marked the lines referred to 1/475 and 1/2500 annual rate of exceedance.The following recurrence times, 475 and 2500 years are equivalent to an exceedance probability of 10% and about 2% in 50 years, respectively.Table 3 highlights the PGA (g) values obtained from the four models, considering the return times of 475 and 2500 years.We explored both the earthquake rate and the ground shaking changes over these two different time periods (1871-2003 and 1871-2016).The choice of these temporal periods was motivated by emulating the periods of time considered in the old ( 2003) and the ongoing new (2016) Italian seismic hazard map.
The purpose was to check whether the use of different temporal catalogs could lead to a difference in terms of earthquake rate and PGA.
As it is seen in Figure 1 the 2016 Amatrice earthquake occurred in an area where the seismic rate is the highest of central Italy.Moreover, the percentage change in the earthquake seismicity rate between two models (Model 2003 andModel 2016) in Amatrice is relatively small being around 5% (Figure 2).This means that in this area the occurrence of new earthquakes does not significantly change the earthquake rate.
We  [Akinci et al., 2009;Akinci et al. 2016;Pace et al. 2006].These studies showed that the contribution of the recently activated source faults vanishes in the overall seismic hazard; the time-dependent PGA values are 20% lower than the Poissonian ones.
On the other hand, some source faults with long elapsed time since the last characteristic event become the most hazardous sites, where the timedependent PGA values are about 50% higher than those of the Poissonian model.Therefore, it is quite important to perceive the impact of earthquake occurrence model on PSHA estimates.
The results achieved in this study brought to the following conclusions: 1) The latest activity occurred over the last 13 years in the central Apennines does not influence the existing seismic hazard model under the Poisson assumption using a smoothed seismicity approach.Two maps obtained in two different time periods are quite coherent in terms of earthquake rate activity.The 2016 Amatrice mainshock occurred in the most seismically active part of central Italy.However, renewal type earthquake rate models need to be taken into consideration for the future seismic hazard assessments in the region.
2) Adopting a proper GMPE is a fundamental step for reliable PSHA calculations, as it is also supported by several previous studies [e.g., Visini et al., 2014;Akinci et al., 2004].The uncertainties in ground motion estimation are very large and need to be reduced in the complex GMPEs.
3) The seismic hazard calculated using the clustered catalogs leads to increased rates of ground motion at both 10% and 2% probability of PGA exceedance in 50 years at Amatrice.The declustered seismic hazard maps may lead to underestimate the seismic hazard in the study region. .

Figure 1 :
Figure 1: Cumulative annual earthquake rate for Model 2003 (1871-2003), considering the shallow events with M w  4.55 for each spatial (0.025° x 0.025°) cell in a log 10 scale.Black dots are the locations of nine earthquakes occurred during the period 01/01/2004-08/31/2016 (M w  4.55).The values close to black dots represent the identification number of events as reported in Table1.With the arrow is indicated the Amatrice earthquake of August 24, 2016 (M w 6.0).

Figure 2
Figure 2 shows the difference in percentage between the two models: "Model 2003" and "Model 2016", given as: (rate "Model 2016" -rate "Model 2003")/ (rate "Model 2003") for each spatial cell.The earthquakes from 1871 to 2003 are indicated with light grey circles and the earthquake from 2004 to 08/31/2016 with black dots.The annual occurrence rate and b-value in the analyzed macro area "Central Italy" during the period 1871-2016 lowered by 4% and 3%, respectively, compared to the early part of the catalog.The rate in the Amatrice zone is slightly decreased (5%), while in the L'Aquila area the rate is unchanged.We also verified if the seismic hazard in terms of ground acceleration in the Amatrice changed between Model 2003 and Model 2016.To do so, we computed the hazard curves using seismicity models described with the GMPE derived from the Italian strong motion database[Bindi et al., 2011] (ITA10).For the Model 2016 we also computed the hazard curve using the GMPE AB10[Akkar and Bommer, 2010],

Figure 2 :
Figure 2: Percentage difference between the Model 2003 and the Model 2016; light grey circles are the events from 1871 to 2003, the black dots indicate the 9 events from 2004 to 08/31/2016 (M w  4.55).The values close to black dots represent the identification number of events as reported in Table1.

Figure 3 :
Figure 3: Hazard curves for Amatrice: Model 2003 -green curve, Model 2016 -red curve, Model 2016 not declusteredblue curve, Model 2016 AB10 -purple curve.Black dotted line represents the 1/475 annual rate of exceedance in 50 years, while the grey line represents the 1/2500 annual rate of exceedance probability.Both X axis (PGA in g) and Y axis (annual rate of exceedance) are in logarithmic scale.

Table 3
PGA values at Amatrice for each computed model