Probability distribution of the macroseismic intensity attenuation in the Italian volcanic districts

We present the probabilistic version of the analysis performed in Azzaro et al. (2006a) on the attenuation of the seismic intensity in Italian volcanic districts. The main results are the estimate of the probability distribution of the intensity at site IS, conditioned on the site-epicenter distance d and on I0, and then, assuming the mode of this distribution as estimator of IS, the forecasting of future macroseismic fields given I0. To this end we have modified the method presented in Rotondi and Zonno (2004) by inserting the following innovative elements: identification of possible different trends and exploitation of knowledge from prior experience or data. Data set. The intensity dataset considered in the present analysis is the same used in the study by Azzaro et al. (2006a), based on a deterministic approach. We consider a total of 38 earthquakes located in the Italian volcanic areas, so distributed: Etna region (24 events), Aeolian Islands (6 events), Vesuvius-Ischia (3 events) and Albani Hills (5 events). The CMTE local earthquake catalogue (Azzaro et al., 2000, 2002, 2006b) has been used for the Etna region while for the other Italian volcanic districts (Aeolian Islands, Ischia, Vesuvius and Albani Hills) the CPTI04 Italian seismic catalogue (Gruppo di lavoro CPTI, 2004) and the DBMI04 associated database (Stucchi et al., 2007) have been considered (Tab. 1). For the analysis, subsets of earthquakes with epicentral intensity I0 ≥ VII MCS and I0 ≥ VI MCS were used for the Etna region and for the other Italian volcanic districts, respectively. Probability model. We cite here the key-elements of the probabilistic method, referring to Rotondi and Zonno (2004) for a detailed description. Instead of adding a gaussian error to deterministic relationships which express the intensity decay as a function of some factors (epicentral intensity, site-epicenter distance, depth, site types, and styles of faulting), we treat the decay as an aleatory variable defined on the domain {0, I0}. Consequently, we assume that the intensity IS is a discrete binomial distributed variable Bin(I0 , p) where pI0 means the probability of null decay, and p belongs to [0,1]. According to the Bayesian approach, p is considered as a random variable following the beta distribution Beta(α, β). Since mean and variance of p are functions of the α, β hyperparameters, we can express our initial knowledge on the decay process through these parameters. To do this, we have divided each macroseismic field in bins of fixed width and the intensity data points in subsets according to this spatial subdivision. For each bin we have repeated the following procedure: a) assessing the prior values to α, β, that is a prior distribution for p; b) updating, through Bayes’ theorem, the hyperparameters on the basis of the current observations; c) estimating the p parameter through the mean of its posterior distribution. By substituting this estimate in the distribution Bin(I0 , p), we obtain an updated binomial distribution indicated as plug-in distribution. Its mode has been assumed as the expected value of the intensity at the sites within the corresponding bin. To predict the intensity at any distance we have smoothed the p’s estimated in the different bins through a monotonically decreasing function; the lowest mean squared error was given by the inverse power function . Hence, the mode of the plug-in distribution obtained by setting p=g(d) provides an expected value for IS at any distance. If, on the contrary, we assume that, from the attenuation viewpoint, the sites inside any bin behave in the same way, we can average over the domain [0,1] of p by integrating the product of the likelihood with respect to the posterior Beta distribution of p. In this way we have obtained the so-called predictive distribution for every bin and its mode is taken as expected value for IS at any site inside that bin. Trends in the intensity decay. We have analysed the macroseismic field of the 38 earthquakes constituting our dataset (Tab. 1) by drawing the decay versus the site-epicenter distance of each data point. A quick look at these graphical representations suggests that these earthquakes do not show an homogeneous decay. To identify different trends in the decay, we have synthetized the information contained in each field by collecting, in a matrix, median, mean, and quartile of each set of distances from the epicenter of the points with the same ΔI. Then we have applied to this matrix a clustering algorithm based on the evaluation of the distance between each pair of rows of the matrix. The dataset has been thus partitioned into two groups of events according to their attenuation trend: the first set mainly formed by the earthquakes of Mt. Etna and Vesuvius-Ischia areas, the second one including the events of the Aeolian Islands and Albani Hills. The set 1 shows an higher decay than the set 2, so two different spatial scales are required: bins of width 1 km for the set 1 and of width 25 km for the set 2. A similar classification analysis was performed in Zonno et al. (2008) on 55 earthquakes representative of the Italian territory; in that case three classes were identified. The probabilistic analysis above described has been separately applied to the two sets, discriminating the events of from those of , and using as a priori distributions for the parameters p’s those indicated in Zonno et al. (2008) for the class of earthquakes with the highest attenuation. The hyperparameters α’s and β’s have been then updated through the observed intensity data points according to the expressions α=α0 + ΣNn=1 IS (n) and β= β0 + ΣNn=1 (I0 - IS (n)). Some results. For each bin the values of the predictive probability function of for the Etna area and Aeolian Islands, are shown in Fig. 1; the squares indicate the values of the intensity decay computed through the logarithmic regressions (Tab. 2) obtained by Azzaro et al. (2006) with the same dataset. These values can be compared with the mode of the predictive function in each bin. The fit between the two methods is good but much more information is provided by the probabilistic approach. In addition to the estimate of the intensity at any site, the probability distribution of IS provides a measure of the uncertainty and its values can be directly used in the software “SASHA” (D’Amico and Albarello, 2007) to calculate the probabilistic seismic hazard at the site. Conclusions. The identification of different decay trends produced by the clustering algorithm matches well with that already presented in the literature (Azzaro et al. 2006), and this suggests that the method could be successfully applied to other cases. Only two earthquakes in Albani Hills - 1876/10/26, I0 VI-VII, 1927/12/26, I0 VII-VIII - are unexpectedly included in the set 1 together with the events of Mt. Etna and Vesuvio-Ischia areas; further, detailed analyses are required to explain such an anomaly. Some problems are still open: a) most of the earthquakes here considered have epicentral intensity I0 VII or VIII, so that we have evaluated the probability functions of IS conditioned on these two values of I0. Also other values of I0 must be used in the analysis; b) the method should be also validated on other earthquakes not included in the dataset of Tab. 1, on the basis of probabilistic measures of the degree to which the model predicts the decay in the data points of a macroseismic field (Rotondi and Zonno, 2004).PublishedTrieste4.2. TTC - Scenari e mappe di pericolosità sismicaope