Bayesian analysis of a probability distribution for local intensity attenuation

Intensity attenuation and its variation as a function of the distance and earthquake size is still a critical issue in evaluating seismic hazard. We present a method that allows us to incorporate additional information from the historical earthquake felt reports in the probability estimation of local intensity attenuation. The approach is based on two ideas: a) standard intensity versus epicentral distance relationships constitute an unnecessary filter between observations and estimates; and b) the intensity decay process is affected by many, scarcely known elements (the physical parameters of the source, propagation path effects, building vulnerability, the semi-qualitative character of macroseismic scales, etc.). Hence intensity decay should be treated as a random variable as is the macroseismic intensity. We assume here that decay, defined on the set {0,1, ..., I0}, follows a binomial distribution with parameters (I0, p); p depends on the distance from the epicenter and is related to the probability of null decay at that distance. According to the Bayesian approach this p parameter is, in turn, a Beta random variable. The observations related to earthquakes with their epicenter outside the area concerned, but belonging to homogeneous zones, are used as prior knowledge of the phenomenon, while the data points of events inside the area are used to update the estimates through the posterior means of the quantities involved. Our methodology is described in detail in an application to the Umbria-Marche area in Central Italy. The data sets examined are the macroseismic intensity database DOM4.1 and the zonation ZS.4, both compiled by the Italian Group for Defence against Earthquakes (GNDT). The method is validated by comparing the observed and the estimated intensity data points of the Camerino (28/07/1799) and of the Colfiorito (26/09/1997) earthquakes. Mailing address: Dr. Gaetano Zonno, Istituto Nazionale di Geofisica e Vulcanologia, Sezione di Milano, Via Bassini 15, 20133 Milano, Italy; e-mail: zonno@mi.ingv.it


Introduction
In the procedures for assessing seismic hazard at site, the hardest issues to tackle are the quality of the data, attenuation, and modelling.Where long historical catalogues are available, as in Italy, it is quite natural to take the macroseismic intensity as a measure of the size of earthquakes.This ordinal quantity, often measured on the 12 degrees of the MCS scale, is in practice treated, as best one can, as an integer variable on {1,12}.The uncertainty in the assessment of the degree of intensity of an event induces us to express the intensity as a range of degrees, VIII-IX, for example, or, correspondingly, by the number 8.5.A more accurate way would assign to each degree a different weight expressing the experts' belief in the membership of the event in that intensity class (Rotondi et al., 1993).This approach is employed in part in Appendix B to take into account half a degree in the counting of null decays.Like all collections of historical data, a parametric catalogue of earthquakes is affected by incompleteness, a factor that can seriously corrupt any statistical analysis performed on that data set.Because of the large amount of information needed few papers that give a direct evaluation of the quality of the data by historiographical analysis can be found in the literature (Lee and Brillinger, 1979).The complete part of a catalogue is generally identified by statistical techniques (e.g., see Rotondi and Garavaglia, 2002); Albarello et al. (2001) use advanced historiographical results to validate the results obtained by a statistical procedure.
As for modelling the definition of a stochastic model that can reproduce the time evolution of the seismic phenomenon in a region is still an open problem.Its solution will take greater knowledge of the physical mechanism generating earthquakes.It is generally accepted that the occurrence probability is best developed from point processes which enable us to estimate a time-varying risk; for a complete survey on this topic we refer to the paper of Vere-Jones (1995).
We address here the problem of attenuation, and local intensity attenuation in particular, since our objective is to assess seismic hazard on the basis of a historical parametric catalogue.Many studies on this topic have appeared in the literature; in the large majority of these the key role is played by the deterministic function which expresses the link between the ∆I intensity decay and factors such as epicentral intensity, site-epicenter distance, depth, site types, and styles of faulting (Cella et al., 1996;Peruzza, 1996;Dowrick and Rhoades, 1999).In some cases a normally distributed random error is added to take into account the scatter (expression of epistemic uncertainty) of the observations around the Is site intensity value predicted through the attenuation relationship.More emphasis is given to the uncertainty when the intensity decay is considered an aleatory variable modelled by a probability distribution.For the intensity decay normalized on I0, a Beta distribution with a mean proportional to an attenuation law and varying deviation was first proposed by Zonno et al. (1993), while a logistic model was used by Magri et al. (1994) to estimate the probability that the attenuation exceeds a threshold value, given a certain distance from the epicenter.Assuming that ∆I is a random variable, the probability of the intensity at site is no longer given by the probability of I0 translated by the attenuation ± the deviation of the error, but is the convolution of the probability distribution of I0 and that of the intensity decay (Magri et al., 1994;Zonno et al., 1995;Tsapanos et al., 2002).Albarello and Mucciarelli (2002) exploit Magri et al.'s probabilistic attenuation (1994) for hazard calculation, while Tsapanos et al. (2002) apply Zonno et al.'s method (1995) to estimate the attenuation for some Japanese cities.
We have made a complete probabilistic analysis of the attenuation issue, avoiding the use of any deterministic attenuation relationship.On the contrary, the emphasis here is on exploiting information from seismogenically homogeneous zones, that is, on assigning and updating prior parameters in the Bayesian framework (Sections 2.1, 2.2 and Appendix A).In Section 3 we give the predictive distribution of the intensity at site, conditioned on I0 and on the d distance from the epicenter, with both discrete and continuous d.We present some validation criteria for these two classes of distributions, and in Section 4 apply these criteria first to an historical event, the Camerino 1799 earthquake, to explore the goodness of fit of the model, and then to a recent event, the Colfiorito 1997 earthquake, to validate the procedure.Appendix B reports further details of the data analysis.In Section 5 we compare our approach with that of Magri et al. (1994), estimating the decay through the mode of the respective probability distributions.

Description of the model
The key element of our approach is the assumption that probabilistic tools are particularly well suited to investigate intensity decay, a process affected by several, not always known and quantifiable features.We have accordingly treated ∆I intensity decay as a random variable that depends fundamentally on the d distance from the epicenter and on the I0 epicentral intensity.Since the variable ∆I is discrete and belongs to the domain {0, I 0}, it is reasonable to choose for Is= I0− ∆I, at a fixed distance, the binomial distribution Bin(I0, p) conditioned on I 0 and p i i - , , , where p∈ [0,1].Consequently, intensities at site can be considered as independent identically distributed random variables conditioned on I 0.Moreover, it is well-known that ground shaking may differ even among sites located at the same distance.To account for this variability the parameter p has been considered a random variable in the Bayesian paradigm.The great flexibility as well as the tractability, of the Beta distribution motivates its choice as conjugate prior distribution for p (2.2) with mean and variance given respectively by 3) The value of the prior hyperparameters α, β expresses our initial state of information on the decay process we might obtain, for instance, from examination of seismogenic areas similar to the zone under study for seismotectonic characteristics, kinematic context, expected rupture mechanism and for seismic wave attenuation.Observation of the intensity data points relative to earthquakes that occurred within the area we are analysing provides additional information D which modifies our previous beliefs concerning the process.Formally this means updating the prior distribution (2.2) of p by Bayes' theorem and obtaining its posterior distribution which, in the present case, belongs to the same family of Beta distributions, since this distribution is conjugate for the binomial likelihood.The learning process is continuously in progress: Be(p|D) represents our beliefs posterior to conditioning on D, but prior to conditioning on any further data provided by future events, which, as we shall see, can be easily included in the process.

Assigning prior parameters
The procedure described here provides initial knowledge of the phenomenon, based on the intensity data points of earthquakes occurring in zones that are homogeneous from the viewpoint of seismotectonics and seismic-wave propagation.Graphic representations of the spatial distribution of intensity report data show that in some cases the circular pattern (points at the same distance dc from the epi/ipocenter have the same attenuation) is altered by the extension of the source and in the direction along the fault strike.In fact the first cause of distortion is statistically significant only up to a distance of half the source length -about 20 km at most in the case of the largest shocks - (Gasperini et al., 1999), while the second cause is partially compensated by the wave radiation pattern, prevailingly perpendicular to the fault strike (Gasperini, 2001).However, if we want to take both of these features into account, we may assume that the pattern is ellipsoidal with the focuses at the extremes of the fault.This means that we assume that the attenuation is uniform at points such that the sum d e of their distances from the focuses varies within a 10 km range.Here the circular or ellipsoidal pattern is obtained by simply taking dc or de as the distance.
Let us consider an earthquake of I0 intensity and draw L distance bins {R1, R2, ..., RL} around the epicenter of width ∆r.The resulting area is divided into bands identified by their maximum distance rj, j=1, 2, ..., L. Taking a sufficiently small step we may assume that all the sites within each Rj band have the same distance dj = =(rj −∆r / 2) from the center of symmetry (epicenter or fault) and that their intensities have the same binomial distribution (2.1) with parameter pj.That is, the decay process is assumed to behave in the same way within each band.Moreover, we denote by " , the set of Nj felt intensities in the j-th band.As the probability that the decay is equal to 0, according to (2.1), is given by if we estimate this probability using simply the frequency of null decay, N j (I 0 )/N j (N j (I 0 ) being the number of any sites in which no intensity decay has occurred), the initial mean value for pj will be (2.5) Once we have assigned a value to the variance σ 2 j of pj, as we do here below, we invert (2.3) to obtain the values for the prior hyperparameters αj,0 and βj,0.This procedure, however, is not applicable in bands where there is no report of null decay.Consequently we have smoothed the valuable p j's through the function f(d)=(c1/d) c 2 and estimated the coefficients c1, c2 by the method of least squares.In this way both an initial expected value and a variance can be assigned to each variable p j, j = 1, ..., L: at each distance d j, j = 1, ..., L, the expected value for pj is given by f(dj); the variance σ 2 1 of p1 in the first band is set as equal to the mean-square error, while to the other variances σ 2 j monotonically increasing values are assigned in the range of variability of the absolute differences between frequency and smoothing function.Some conditions must be satisfied to ensure that the values thus chosen for the prior mean and variance of pj provide feasible values for the prior hyperparameters αj,0, βj,0 when substituted in (2.3).Appendix A contains a list of these conditions and additional remarks regarding this issue.Moreover, in order to have variances σ 2 j that are strictly positive, the pj parameters must vary within the range (0, 1); we have set p j =0.98 when f(dj)>0.98. Figure 2 shows some examples of smoothing related to the case study described in detail in Section 4.

Updating parameters
Let us consider all the earthquakes of I0 intensity with an epicenter within the area under study.We denote by " , the set of their intensity data points, subdivided into L subsets, each containing the intensities felt within the corresponding bands.
On the basis of this new information we update our knowledge on the attenuation process.We use Bayes' theorem to compute the posteri-or distribution Be (pj| Dj), and estimate pj through its posterior mean where i (n)  s is the intensity felt at the n-th site inside the Rj band, Nj is the total number of data points assumed at dj distance from the epicenter, and αj,0, βj,0 are the parameters of the prior Beta distribution for pj.
In order to let the p parameter of the binomial distribution for the intensity Is at site vary with continuity, we smooth the estimates p ˆj, j=1, ..., L with the method of least squares, again using an inverse power function g(d)=(γ1/d) γ 2 .In this way it is possible to assign the probability of the intensity decay Pr{∆I |I0,g(d)} at any distance from the epicenter, as required in the validation criteria presented in the next section.Figure 3 shows the comparison between prior and posterior smoothing curves obtained for different epicentral intensities in an Italian seismogenic zone.

Validation
To predict the intensity points an earthquake will generate, on the basis of the knowledge accumulated before its occurrence, we have two possibilities.We can apply either a predictive probability function for all the points within every Rj band, or use a different binomial Bin(I 0, g(d)) probability function (2.1) for the points at distance d from the epicenter, where g(⋅) denotes the smoothing inverse power function given in Section 2.2.
In general, given a sample distribution f(x|θ) and the posterior density function p(θ |D) for the θ parameter, conditional on data D, the density of the predictive distribution is given by In our case we have where αj and βj, parameters of the Be(pj |Dj) posterior distribution, are updated on the basis of the past information as follows: One of the most frequently used probabilistic measures of the degree to which a model predicts the data is the so-called logarithmic scoring rule, based on the logarithm of a posterior probability (Lindley, 1987;Winkler, 1996).We use both the marginal likelihood and the binomial distribution to evaluate this measure, obtaining respectively the following expressions: where N′ s is the number of sites, within the L bands Rj, j=1, ..., L, at which the future event will be felt, D′= ∪ L j=1 D′ j is the set of its i (n′)  s intensity points, and d n′ being the distance of the n′-th site from the epicenter.
Another probabilistic measure of the fit is given by the p(A)/p(B) ratio between the probability that the fitted model assesses to the realization A and the probability of the predicted value B. The idea behind this measure is borrowed from the concept of deviance (Read and Cressie, 1988) and is based on a consideration of how much is gained from having predicted B when A occurs.Of course, the gain is maximum, that is equal to 1, when we have predicted what then really occurs.Let us choose the mode of a posterior probability as the predicted value, and denote the mode of the marginal likelihood, and that of the «smoothing» binomial distribution respectively by The errors, expressed in probabilistic terms, are given by the geometric means of the corresponding odds in logarithmic scale (3.5) The deterministic absolute discrepancy between observed and estimated intensities at site -applying both the marginal likelihood and the binomial distribution -is given instead by (3.6) We should like to make it clear that, in the case studies presented in Section 4, when the Is intensity observed at site is not expressed by an exact degree of the MCS scale, the integer part of its value has been used in evaluating the above-mentioned posterior probabilities which are defined on the discrete set {0, 1, ..., I0}.Instead the values actually recorded in the macroseismic intensity database DOM4.1 are used in the evaluation of the discrepancies.

Case studies
The sequence of the 1997-1998 Umbria-Marche (Italy) earthquakes began on 26 September 1997 with two main shocks.The first event, ML= 5.6, occurred at 00:33 GMT, the second, ML= 5.8, several hours later, at 09:40 GMT.A third quake, ML= 5.5, occurred on 14 October at 15:33 GMT (Amato et al., 1998).According to Cattaneo et al. (2000) the sequence, which lasted until April 1998, was characterized by shallow earthquakes (less than 9 km deep, with few exceptions) which show coherent focal mechanisms indicating normal faults with NW-SE striking and SW dipping planes.The damage ascertained in the days immediately following the main shock has been associated with the IX degree of the MCS scale.The epicenter was located at Colfiorito which, according to the ZS4 zonation (Scandone et al., 1992;Meletti et al., 2000), belongs to the seismogenic zone number 47.The only other shock of the same size recorded in this zone is the Camerino earthquake of 28 July 1799.
The ZS4 zonation (available at http://emidius.mi.ingv.it/GNDT/P511/home.html)describes seismogenic zones judged homogeneous from the viewpoint of kinematic context and expected rupture mechanism.Zone 47 belongs to the set Z = {28, 29,32,33,34,36,37,44,45,46,47,50,51, 52} of intermediate zones characterized by mixed expected rupture mechanisms, for the most part dip-slip.These zones, shaded in fig. 1, are linked with the passive sinking of the Adria lithosphere beneath the mountain chain in the Northern Apenninic Arc.As most of the area is located south of latitude 43°30′, it is, according to Gasperini (2001) and Mele et al. (1997), also homogeneous with respect to wave attenuation.The NT4.1 catalogue (Camassi and Stucchi, 1996) records 501 earthquakes within this set, 89 in zone 47 and 412 in the other zones.Macroseismic observations in the DOM4.1 database (Monachesi and Stucchi, 1997) cover only 169 of these: 32 of the 89 shocks in zone 47 and 137 of the 412 distributed over the remaining zones (here below denoted by D47 and DZ−47 respectively).The DOM4.1 data set dates from 1000 to 1992, and does not include the intensity points of the Colfiorito earthquake of 26 September 1997.It does contain some 37 000 macroseismic observations regarding more than 900 earthquakes at more than 10 000 sites, as well as 935 maps.This is just a portion of the general database of macroseismic observations for the Italian area, based on the studies of the Italian Group for Defence against Earthquakes (GNDT) and of other scientific organizations, used to compile the NT4.1 catalogue (Camassi and Stucchi, 1996).These macroseismic data would allow us to assess the physical dimensions and orientation of the source of the large historical earthquakes constituting the data set under study (Gasperini et al., 1999), but that is beyond the aim of this paper.For our purpose here we have adopted the circular pattern on the basis of the available information.
Let us suppose that the homogeneity of the zones belonging to set Z can be also extended to the decay process.On the basis of this assumption the data set DZ−47 constitutes the prior information that allows us to assign the parameters αj,0 and βj,0 of the prior density for the pj parameter characterizing the binomial distribution (2.1) of Is in the Rj band.Applying the smoothing procedure described in Section 2.2, we obtain the estimate of p at any distance for epicentral intensities in the range of V-XI MCS degrees.The expression of the corresponding smoothing curves and their graphic representation are given in fig. 2. Macroseismic observations that were not deemed sufficiently reliable have been ignored (see the comment on this point in Appendix B).The events of sets DZ−47 and D 47 , which have been analyzed, are listed in   ber of points is indicated when there are only 1 or 2 occurrences.The agreement is quite good up to around 40 km; beyond that the attenuation of Camerino appears greater, but the data are very scarse.The fit has been measured quantitatively, applying the validation criteria presented in Section 3: the logarithmic scoring rule (3.2), (3.3), the odds criterion (3.4), (3.5), and the discrepancy (3.6) between estimated and observed intensities, using both the predictive and the binomial distribution for Is .The values of these criteria are given in table IV.The graphs on the right of fig. 5 show the discrepancy between the intensity recorded at site and the posterior mode of the predictive (top) and that of the binomial distribution (middle).The graphs on the left display, besides the estimated (blue asterisks) and observed (red crosses) intensities at site, the shortest intervals (green bars) in which the Is intensity falls with a probability of at least 70% according to the above-mentioned distributions.The set of intervals provided by the predictive distributions covers 93% of the observations, while 67% are included in the intervals of the binomial distributions.Hence it is not surprising that the best results (table IV) are provided by the predictive distribution in a backward check.Now let us look at the forward validation.By plugging the updated parameters α j , β j in (3.1)   that the Colfiorito earthquake has an attenuation pattern which contrasts in part with our expectation based on prior information: in particular, the decay is quicker within the 40 km limit and slower at a greater distance.
We might add that the Colfiorito earthquake is characterized by the wide spread of the in-tensities reported within 20 km of the epicenter.In many earthquakes these intensities at settlements located a few kilometers from one another vary up to 2 or 3 degrees, due to peculiar geotechnical and topographic characteristics of the sites (Gasperini, 2001).In the Colfiorito earthquake variations of more than 2 degrees of intensity have been observed within 1 km, moving from the inner part of the basin and from bedrock sites towards the edge of the basin (Tertulliani, 2000).But further analysis of this anomaly is beyond to the scope of this paper.
Let us now evaluate the validation criteria (3.2) to (3.6) for the set of intensity data points summarized in table V.The results are reported in table VI.The attenuation trend, not fully consistent with either prior information or Camerino, leads to results less clear than those for Camerino.According to both the scoring rule and the odds criterion, the predictive distribution (3.1) should be preferred to the binomial with the p parameter given by the «smoothing» function.In fact the former assigns a greater probability to the future observations than the latter, and the error in probability that the predictive distribution produces is smaller than that of the binomial, as the smaller value of the odds of the predictive distribution shows.According to the discrepancy criterion (3.6), instead, if we use   as predictor the posterior mode of the binomial distributions, rather than that of the predictive distributions, the error in assessing the intensity at site is smaller -less than one degree on the average.The sensitivity of the criteria to the two distributions is represented graphically in the top and center graphs of fig. 9.The fact that the predictive distribution produces (at least) 70% highest probability intervals (the green bars in fig.9) longer than those of the binomial explains the results provided by the scoring rule and by the odds criterion: the greater number of intensity points (red crosses) covered by those intervals (84% for the predictive against 80% for the binomial) implies a greater probability assigned to the observations by the predictive than by the binomial distribution.In the same graphs the blue asterisks indicate the intensity at site estimated through the posterior mode.The graphs on the top and center, right in fig. 9 display the differences between the recorded intensity at site and the posterior mode of the predictive (top right) and of the binomial distribution (middle right).We point out that, since I s is a discrete variable, at least 70% highest probability intervals may actually have a probability of even 80%-90%.

Comparison with another probabilistic method of attenuation
The literature records two main ways of dealing with local attenuation relationships: a) the intensity decay is considered a function of I 0 and of the site-epicenter distance, and the unknown parameters are estimated on a representative earthquake of the area under study; or b) the decay depends only on the distance from the epicenter, and the estimation takes into account all the sites in the database.Examples of these deterministic relationships in Italy are given by Grandori et al. (1987) and Berardi et al. (1994).Magri et al. (1994), instead, presented a probabilistic approach to the study of attenuation.They consider the probability that an earthquake has been felt at site with an intensity less than or equal to a is value, conditioned on the epicentral intensity I0 and given a function of the epicentral distance d, e.g., the natural logarithm logd.On the basis of an exploratory analysis of the observed decays for earthquakes in Italy, they state that the logistic distribution describes this probability well   Now let us go back and repeat our validation of the Camerino 1799 earthquake applying the method proposed in Magri et al. (1994).The intensities felt and those estimated through the mode of the discrete distributions obtained at site from the logistic model (5.2) are represented in the graph on the bottom left in fig.5, while the graph on the right shows the respective set of discrepancies.The results of the validation criteria are given in table IV.In the Camerino case we can conclude that the predictive distribution provides the best results for all the criteria.We also see that 70% of the observations fall into the (at least) 70% probability intervals, represented by the green bars in the graph on the bottom left of fig. 5, against 93% for the predictive distribution and 67% for the binomial.
In the case of the Colfiorito 1997 earthquake, the divergence between prior information and observations discussed in Section 4 and shown in fig.4, is also reflected in the heterogeneous results provided by the validation criteria (see table VI).Each distribution gives the best result for a different criterion; in particular the value of the scoring rule is high for both the predictive and the binomial distribution because at sites within 20 km of the epicenter these distributions concentrate most of the probability on an intensity of the VII to IX degree, whereas in many of those sites the felt intensity was lower.The logistic model offers a good description of the attenuation trend within the 20 km interval, but predicts excessively rapid intensity decay beyond the 40 km limit.Only 73% of the observations fall in the (at least) 70% probability intervals (the green bars in the left-hand graphs of fig.9) compared with 84% and 80% of the observations included in the probability intervals of the predictive and of the binomial distribution respectively.
If an inadequate learning set were to blame, then adding the data for the Colfiorito 1997 earthquake to update the model parameters would be expected to improve the performance of the different methods.Repeating the inferential procedure with the addition of these data, we obtain the relationships a ∆I =3.418 + 0.21∆I and b∆I =−0.4457−1.453∆Iinstead of (5.3), and f(d)=(6.017/d) 0.223 as the smoothing function.The results for the validation criteria, given table VII, show that all the methods have good learning ability; the predictive distribution in particular scores the minimum total given by the sum of the values for the three criteria.We also add that when we used the parameters estimated at national scale by Magri et al. (1994) in the logistic model, these produced the worse results in all the cases; e.g., in the forward validation of the Colfiorito earthquake we had scoringlog=1.697,odds log =0.705, difflog=0.896,all values larger than those in table VI obtained using the regional estimates (5.3) of the parameters.

Conclusions
This paper presents a complete analysis of the problem of the assessment of intensity attenuation from the probabilistic viewpoint.Bypassing the traditional attenuation relationships and respecting the semi-qualitative character of the macroseismic intensity, it treats decay as a random variable.The probability distribution of the I s intensity -a binomial -is conditioned on the epicentral intensity and on the distance between site and epicenter.In the Bayesian context the p parameter of this distribution is also a random variable, and we have adopted a conjugate Beta distribution for it.This approach allows us to exploit various types of information in the inference procedure, and permits easy updating of the estimates when new data is available.Moreover, the uncertanty on the p parameter can be expressed explicitly, giving, together with the estimate, the shortest interval in which p varies with some probability.We have also described in detail how we elicit the prior hyperparameters.Some new validation criteria, both probabilistic and deterministic, have been proposed.In particular, backward analysis has been performed on the Camerino earthquake, to check the goodness of fit of the model, and both backward and forward analysis, on the Colfiorito earthquake, to test the forecasting capability of our method.The results have been compared with another probabilistic method for the analysis of seismic attenuation based on a logistic model.On the basis of the validation criteria proposed, we believe that the predictive distribution can be indicated as the all-around best of the probabilistic models for seismic attenuation examined, also considering that a conservative attitude is preferable in hazard studies.

Fig. 1 .Fig. 2 .
Fig. 1.The ZS4 zonation of Central Italy (courtesy of the Italian Group for Defence against Earthquakes).

Fig. 3 .
Fig. 3.The inverse power function approximating the pj's assigned by prior knowledge (dotted curve, solid dots) and updated by observations (solid curve, open squares), for different epicentral intensities.

Fig. 4 .
Fig. 4. Left: intensity data points for the earthquakes of intensity IX in DZ−47 (black) and for the Camerino (top) and Colfiorito (bottom) earthquakes (red); right: the arithmetic averages of the previous intensities at intervals of 10 km from the epicenter.The cardinality of the smallest -1 or 2 occurrences -data sets is indicated.

Fig. 5 .
Fig. 5. Camerino 1799 earthquake.Left: scatter diagrams of observed MCS intensities (red crosses) with respect to distance, of intensities (blue asterisks) estimated at the same sites through the (posterior) mode of the predictive distributions, of the «smoothed» probability densities, and of the logistic distributions; the green bars denote the shortest intervals with at least 70% (posterior) probability of Is.Right: discrepancies between estimated and observed intensities; the red open circles indicate the sites in which the discrepancy is null.

Fig. 6 .
Fig. 6.Estimate of the predictive probability function for the intensity Is.
values of ∆I, estimates of the parameters a and b are obtained by maximizing the likelihood.They have then modelled the dependence of a's and b's from ∆I by the linear approximations a∆I=q + t∆I and b∆I =u + v∆I.this method to the set formed by the DZ−47 earthquakes of I0= IX and by the Camerino 1799 earthquake, obtaining the following relationships: in the different dj, j=1, ..., L, we have, by difference, the discrete probability functions for I s represented in fig.8.

Fig. 8 .
Fig. 8. Estimate of the probability function for the intensity Is drawn from the logistic model (5.2).

Fig. 9 .
Fig. 9. Colfiorito 1997 earthquake.Left: scatter diagrams of observed MCS intensities (red crosses) with respect to distance, of intensities (blue asterisks) estimated at the same sites through the (posterior) mode of the predictive distributions, of the «smoothed» probability densities, and of the logistic distributions; the green bars denote the shortest intervals with at least 70% posterior probability of Is.Right: discrepancies between estimated and observed intensities; the red open circles indicate the sites in which the discrepancy is null.

Table I .
DZ−47 data set: intensity data points, drawn from the DOM4.1 database, for earthquakes in the seismogenic zones of the Z set, excluding zone 47.

Table II .
Summary of the DZ−47 data set employed as prior information in the inferential procedure.

Table III .
D47 data set: intensity data points, of earthquakes in the seismogenic zone 47, by means of which the prior parameters are updated.

Table IV .
Criteria of backward validation applied to the Camerino 1799 earthquake; * denotes the best result.

Table VI .
Criteria of forward validation applied to the Colfiorito 1997 earthquake; * denotes the best result.

Table VII .
Criteria of backward validation applied to the Colfiorito 1997 earthquake; * denotes the best results.