Rapid scenarios and observed intensities

After a destructive earthquake, national Governments need to know the approximate amount of damage, the number of casualties, and the financial losses as soon as possible. Rapid scenarios are also used to inform the general public; see the widely used Shakemap package [Wald et al. 1999, 2006] of the US Geological Survey (USGS) and the one modified by the Istituto Nazionale di Geofisica e Vulcanologia (INGV; National Institute of Geophysics and Volcanology), which is reproduced for Figure 1. The general matter of the use of intensities in damage scenarios was discussed in a special session at the 2008 Annual Meeting of the Seismological Society of America (http://www.seismosoc.org/meetings/2008/specialsessions.html), and was also discussed in the NIS-1 session of the European Congress in Moscow, in August 2012 (http://www.esc2012-moscow.org/esc_thematicareas.html). The purposes of the present report are to: (i) compare different types of intensities; (ii) check two rapid scenarios of intensity; and (iii) understand whether the KF formula [Sirovich 1996, Sirovich et al. 2009] can be used as a new 'attenuation' relationship to improve rapid scenarios. […]


Introduction -applicative interest
After a destructive earthquake, national Governments need to know the approximate amount of damage, the number of casualties, and the financial losses as soon as possible.Rapid scenarios are also used to inform the general public; see the widely used Shakemap package [Wald et al. 1999[Wald et al. , 2006] ] of the US Geological Survey (USGS) and the one modified by the Istituto Nazionale di Geofisica e Vulcanologia (INGV; National Institute of Geophysics and Volcanology), which is reproduced for Figure 1.The general matter of the use of intensities in damage scenarios was discussed in a special session at the 2008 Annual Meeting of the Seismological Society of America (http://www.seismosoc.org/meetings/2008/specialsessions.html),and was also discussed in the NIS-1 session of the European Congress in Moscow, in August 2012 (http://www.esc2012-moscow.org/esc_thematicareas.html).
The purposes of the present report are to: (i) compare different types of intensities; (ii) check two rapid scenarios of intensity; and (iii) understand whether the KF formula [Sirovich 1996, Sirovich et al. 2009] can be used as a new 'attenuation' relationship to improve rapid scenarios.
Here, we focus on the meizoseismal area of the May 20, 2012, Emilia (Po Valley, north-eastern Italy) earthquake, which includes all of the damage and was the immediate target of the Italian Department of Civil Protection (Dipartimento Protezione Civile, DPC).In this area, the DPC needed to estimate the amount of damage, the number of casualties, and the financial losses in the first tens of minutes after the event.Strictly speaking, none of the pictures we are going to present can fulfill all of the requirements of the DPC, but we hope that this analysis is useful to improve the available algorithms for rapid scenarios.
Our reference for comparisons are the intensities observed in the field by seismologists and engineers.
The surveyed sites are shown as black dots in Figures 1  and 2. We used the natural neighbor approach (n-n; based on the Voronoi polygons) to draw isoseismals.Let us recall that the n-n isolines [Sirovich et al. 2002] are obtained without any explicit or implicit assumption about the observed phenome-non, and thus they strictly honor the data.Therefore, in Figures 1C-F the intensity observed at each site can be read simply by looking at the color of the area where the black dot is, because each site is inside its correct isoseismal area.With the same data, everybody obtains the same n-n isoseismals; the use of these isolines helps in the understanding of the phenomenon better than the Italian tradition of plotting dots in colors.
Let us remember that any collection of field data (no matter if from field intensities or web answers) samples the regional seismic response.The sampling must obey the spatial version of the Nyquist principle along the coordinate axes [Press et al. 1992, pp. 494-495].We will see that this is crucial to avoid ill-posed comparisons of contoured isoseismals.
We want to stress that we take as reference the field intensities produced only by the main shock of May 20, 2012 (not the cumulative damage due to the subsequent shocks).
Finally, three warnings about the Figures.First, Figure 1A-F shows an area of the Po Valley between Ferrara (right side, middle) and Mantova (top left corner; both not shown); this can be checked from the longitudes and latitudes on the X and Y axes.Second, we apologize for the quality of some Figures, which were prepared in a hurry, by different Institutes, with different tools, and at different scales, and they were greatly enlarged.Third, it also needs to be considered that different Institutions adopt different color conventions for intensities.Some of them use soft-shade colors; the Istituto Nazionale di Oceanografia e di Geofisica Sperimentale (OGS; National Institute of Oceanography and Experimental Geophysics) uses full colors, in agreement with the traditional definition of intensities as degrees of a scale.When reviewing the new interpolation scheme of Shakemap [Worden et al. 2010], we suggested adding thin black isoseismal lines that correspond to the full intensity degrees.Indeed, soft-shade colors call to mind more continuous variations of real numbers than degrees of a scale, and render the understanding of the shape of the regional intensity pattern difficult.Galli et al. [2012a, b]   central shades of the colors adopted by Shakemap and other packages (see the figures for clarification of the situation).

Source parameters used
We used the source parameters available on May 20, 2012, soon after the earthquake to calculate our scenario in quasi-real-time.Seven source determinations were available, from the INGV (first report Relazione sul terremoto e la sequenza in corso agg.20/5 ore 8 and the Time Domain Moment Tensor, USGS Centroid Moment Solution, body-waves solution and W-phase, EMSC-CSEM agency, OGS Udine).We discarded two of these (USGS body-waves and OGS) for seismological reasons and also because they were discordant.Table 1 lists the ranges of values explored parametrically and the sampling steps used to build our scenario; a barycentric position for the nucleation was used.The limits of the ranges coincide with the extreme parameters available.In this way, 108 sources were obtained.

Rapid forecasts, and web intensities
Figure 1A shows a strong enlargement of the MCS intensity map produced by the INGV using its improved version [Michelini et al. 2008] of the Shakemap software package [Wald et al. 1999[Wald et al. , 2006]]; only the meizoseismal area is shown.See http://shakemap.rm.ingv.it/shake/archivefor the map we used, which is version #7 of June 2, 2012.Note that the INGV can consider the influence of the extended source (Figure 1A, small rectangle).A slight tendency towards W-E elongation is seen; the early versions of the INGV Shakemap of the May 20, 2012, event had a more rounded concentric aspect.
The MCS intensity map produced shortly after the earthquake from web answers (see www.haisentitoilterremoto.it) by P. Tosi, V. De Rubeis and P. Sbarra of the INGV is shown in Figure 1B.Note that the complete map available at the time used 10,294 answers (sites not shown in Figure 1B).It is worth commenting that the INGV procedure resembles the well-known "Did you feel it?"reports of the USGS, but it adopts its own forms [see details in Sbarra et al. 2010].Note that in Figure 1B the red signifies VII; the area of degree VI is also marked by a black line; regarding this, Figure 1B is selfexplanatory, however.Both isoseismals show E-W orientations, with a slight tendency toward WNW-ESE.
Figure 1C repeats the tentative rapid scenario of the OGS that we sent to some colleagues of the DPC and the INGV just after the May 20, 2012, earthquake, for a totally blind test; in other words, we offered them the opportunity to check in the future on the quality of our results by comparing them with the intensities that were still to be collected via the web, or observed in the field.At that time, we sent the median scenario, which was calculated on a regular grid.In Figure 1C, the scenario already sent to the DPC was resampled in the sites used to build Figure 1B.
Usually, we compare field results and forecasts by calculating the summed squares of the differences in the same sites, and their distribution.Here, we can propose only qualitative comparisons.As said, geographical sampling is crucial to draw contour lines of any geophysical more or less regular field [see also Sirovich et al. 2002].Thus, the comparison between field intensities (observed in sparse sites) and intensities calculated on a regular grid would be ill-posed.For this reason, in view of an advanced blind test, with the purpose of a well-posed graphical comparison between Figure 1B and 1C, we plotted again our rapid scenario using the coordinates of the web answers collected by the www.haisentitoilterremoto.it Group [De Rubeis, written communication, 2012]; obviously, the parameters used on May 21, 2012, remained fixed.The median values minus one absolute deviation are shown in Figure 1C (the deviation is a few tenths of a degree); only the coordinates with at least three web answers were used.
The OGS scenario was obtained with the KF formula model, which simulates the horizontal component of the body-wave radiation from an extended source.Let us recall that KF uses 11 source parameters, which identify a fault of unit width, with the three coordinates of its nucleation, the fault-plane solution, the seismic moment, the rupture velocities, and the shear-wave velocity in the half space.These parameters are available from agencies in the first hours after a destructive earthquake; the four further parameters used (see these and the KF in Sirovich and Pettenati [2009]) are not rapidly available, and thus they must be fixed.In the study case, the rupture was kept symmetric, and the fixed parameters are in Table 1.For the description of the direct, parametric, use of KF to build scenarios, we refer to Sirovich and Pettenati [2009].The scenario of Figure 1C used the data available on May 21, 2012, but could have been produced in some tens of minutes after the event.Similar KF scenarios were tested for the Los Angeles 2008 M 5.4 earthquake [Sirovich et al. 2009] and for the L'Aquila event of 2009 [Pettenati et al. 2011].
Figure 1D repeats the scenario of Figure 1C, with the only difference being in the geographical sampling.Indeed, we used the sites surveyed by the DPC with the purpose of a well-posed graphical comparison with Figure 1E.Note the importance of the sampling, which in this case was inadequate (see Section 6).
Figure 1E shows the DPC results (52 sites, MCS scale).We stress that the intensities shown were produced by the main earthquake of May 20, 2012 [Galli et al. 2012a, b].The n-n contour was used, with the same color convention as in Figure 1C, D and F. Most likely, the uniform green area (V) in Figure 1E is due to the absence of peripheral sampling points (sites surveyed by the DPC).
The intensities from Quest [Arcoraci et al. 2012] are shown in Figure 1F (in the EMS98 scale).We found the geographical coordinates of 38 intensity values and two 'D' classifications (for hamlets with a few buildings; we put D = VI for Buonacompra and Reno Centese).The same comment of the green area applies as for Figure 1E, but the absence of sampling points is also responsible for the yellow (VI) area that opens eastwards.It can be reasonably assumed that all of the damaged sites eastwards were visited (and thus the VI isoseismal should close), but the V degree was not sampled there.

Preliminary comparisons
Figure 2A compares the EMS98 intensities by Quest with the MCS values by the DPC in the 26 localities visited by both.There is a slight discrepancy around the damage/ no damage limit (twelve VI degrees for Quest, four for the DPC), but no more systematic differences are seen.
Figure 2B shows the same comparison as Figure 2A, but in the Emilia territory (see longitude and latitude in Figure 2B).The black dots are the 26 sites visited by both the INGV and the DPC. Figure 2B shows the DPC intensity minus the INGV intensity.The Voronoi polygons help to understand the local differences.The residuals can be seen to be scattered.There are only three dark gray polygons (DPC -INGV = +0.5),which confirms that the EMS98 estimates were slightly higher than the MCS estimates.
Figure 2C compares the MCS intensities of www.haisen titoilterremoto.itwith the MCS field values of the DPC in the 24 localities they have in common.The numbers added to the dots indicate the number of the localities in common; consider that each locality of www.haisentitoilterremoto.it corresponds to one or more web answers.For example, in Figure 2C there is one locality (the inhabited center of Bondeno) where the DPC obtained degree 5.5, while www.hai sentitoilterremoto.itreceived six answers, which produced degree VIII.This strong discrepancy remains to be understood.However, what is more important is that there is some dispersion around the V and VI degrees, without systematic effects (the misfit of high intensities with the web data is not surprising as people are frightened and have many difficulties in the most damaged areas).So, Figure 2C confirms the acceptable performance of www.haisentitoilterremoto.it for degrees V and VI, as already commented on in Figure 1B.

Discussion and comments
We have verified beyond any doubt that the black dots in Figure 1C are more than sufficient to describe the intensity pattern in this Figure .Instead, spurious effects are seen in Figure 1D, E (sampling by the DPC), and in Figure 1F (sampling by Quest), due to scarce peripheral sampling, and thus peripheral divergence of the Voronoi polygons.This is the reason why degree V remains constant with distance in these Figures, which is not realistic; the same holds for the open shape of VI in Figure 1F.Insufficient sampling in Figure 1F might have applicative fallout in the definition of the area that was most stricken by the earthquake.Still on sampling, the differences between Figure 1C and D (same model, same source parameters) are due only to the denser sampling in the central part of Figure 1D than in Figure 1C.On the contrary, the almost total disappearance of degree IV in Figure 1D is the spurious effect of the lack of surveyed sites in peripheral areas.The isolated point of degree V to the north of Figure 1D is one more piece of evidence that the n-n approach honors the data.
Then we assume that the reference intensity data for validations are the intensities observed by seismologists in the field.From this point of view, the existence of two surveys in two different scales complicates our attempt at validation of rapid scenarios and web answers.First of all, for the use of different scales, we believe that the uncertainties brought into our problem by the MCS-EMS98 scale differences are less important than the many uncertainties in the present comparisons.As a very first approximation, we are satisfied by the authoritative opinion by three renowned experts [Musson et al. 2006]: "Ideally, direct conversion between intensity scales should never be made", but, "these values are likely to vary more between two seismologists using the same scale than between two scales used by the same seismologist", and thus the "relationship between the major 12-degree scales, such as MSK, MMI and MCS, and EMS-98 is more or less 1:1".We are aware that the matter is questionable, however.Then, Figures 1E, F and 2A, B show that 26 localities were visited by both the DPC and the INGV, and that different groups often gave different interpretations of damage.
Both references in Figure 1E, F show that the area of VI is elongated in the WNW-ESE direction, which is the orientation of the fault source.The DPC gives two sites with VII, Quest gives four (one is in common).The five sites are also approximately aligned in the WNW-ESE direction.We conclude that the heaviest damage is there, even if the two groups of experts diverged a little bit with their final interpretations.We conclude that Figure 1E, F can be used together as reference.Thus, our comments are as follows.
Notwithstanding some difficulty in the interpretation of the shaded colors in Figure 1A, its area of VI appears to be overestimated; the influence of the recording site (Figure 1A, triangle) needs to be better understood.The web intensities of Figure 1B behave well; also degree VII is approximately in the correct position (this is not trivial).The OGS rapid scenario overestimates the DPC and Quest references, very roughly by half a degree in the central area.The KF rapid scenario was able to match the WNW-ESE orientation of the meizoseismal area of the reference, which follows the orientation of the known causative fault at depth.In Figure 1D, the limit between VII and VI is not well constrained northwards by sampling.In the KF model, this WNW-ESE orientation is the combined effect of the orientation of the extended fault, its directivity, and the radiation pattern of the horizontal component of shear waves produced by the parameters of Table 1.As said, KF uses many sources parametrically (e.g., different positions, mechanisms), and some overestimation is intrinsic.
Both references are not reliable close to the borders, although it can be said that Figure 1C (OGS, with good sampling) and Figure 1B are in reciprocal agreement on degrees V and IV, at least northwards.Finally, let us stress that the rapid KF scenario is the only forecasting procedure mentioned in this report, which starts only from source information.In other words, our scenarios cannot take advantage of any rapid information from the field.Thus, they cannot compete with ShakeMap.More appropriately, the KF formula can be used to improve ShakeMap or other packages by using KF as a new 'attenuation' relationship for intensity.

Figure 1 .
Figure 1.(A) Close-up of the INGV Shakemap (MCS scale) [Michelini et al. 2008] (version #7 June 2, 2012).(B) Close-up of the INGV map www.haisentitoilterremoto.it(MCS) by Tosi, De Rubeis and Sbarra, from 10,294 web answers.Note that the color scale is different from (A): red represents VII.(C) The rapid scenario sent from the OGS to the DPC and the INGV one day after the event; resampled in the web answers of www.haisentitoilterremoto.it [De Rubeis, written communication 2012].(D) As for (C), but resampled at the sites surveyed by the DPC.(E) MCS intensity map obtained by the DPC (52 sites) [Galli et al. 2012a, b]; n-n contours, same color convention as in (C, D, F) (see text).(F) EMS98 intensities by Quest [Arcoraci et al. 2012] for the main shock of May 20, 2012; same contour and color as in (C, D, E).

Figure 2 .
Figure 2. (A) Scattergram of the EMS98 intensities by Quest versus the MCS values by the DPC in the 26 localities visited by both groups.The numbers added to the dots indicate the number of the localities in common.(B) Same comparison as (A), but in the Emilia territory.Black dots, 26 sites visited by both Quest-INGV and the DPC.The gray scale shows the difference (DPC intensity minus the Quest INGV intensity).(C) Scattergram of the MCS intensities of www.haisentitoilterremoto.it versus the MCS field values by the DPC in the 24 localities in common.The numbers added to the dots indicate the number of the localities in common (see text for the web answers).
(see later) and the 'www.haisentitoilterremoto.it' Group of the INGV added these isoseismals, but use their own color scales.The full colors of the OGS try to follow the

Table 1 .
Source parameters explored parametrically to calculate the scenarios of Figures1C and 1D.