Earthquake forecasting in Italy , before and after the Umbria-Marche seismic sequence 1997 . A review of earthquake occurrence modeling at different spatio-temporal-magnitude scales

The main goal of this work is to review the scientific research carried out before and after the Umbria-Marche sequence related to the earthquake forecasting/prediction in Italy. In particular, I focus the attention on models that aim addressing three main practical questions: was (is) Umbria-Marche a region with high probability of occurrence of a destructive earthquake? Was a precursory activity recorded before the mainshock(s)? What was our capability to model the spatio-temporal-magnitude evolution of that seismic sequence? The models are reviewed pointing out what we have learned after the Umbria-Marche earthquakes, in terms of physical understanding of earthquake occurrence process, and of improving our capability to forecast earthquakes and to track in real-time seismic sequences. Mailing address: Dr. Warner Mazzocchi, Istituto Nazionale di Geofisica e Vulcanologia, Via di Vigna Murata 605, 00143 Roma, Italy; e-mail: warner.marzocchi@ingv.it


Introduction
The italian territory is characterized by a generally high seismic risk.Although the strongest earthquakes range only from medium to relatively large magnitudes (up to about 7), the high density of inhabitants and the age and quality of buildings make the vulnerability and exposure rather high.
The reduction of seismic risk is one prominent (and obvious) objective of National Au-thorities in Italy as well as in many other countries.Such mitigation requires a prioritization of actions devoted to this aim, such as a proper land use planning, the reinforcement of building and infrastructures to diminish vulnerability, and a sound management during seismic crises.All of these issues entail, amongst others, the use of reliable earthquake forecasting/prediction models.
The paramount importance of these topics was highlighted abruptly by the occurrence of Umbria-Marche earthquakes in 1997 with the consequent damages and death toll.In this review, I focus the attention on defining the state of the art on this matter just before the Umbria-Marche earthquake and on the improvements made in the following ten years.In particular, I point to the description of three different topics: 1) Large earthquake forecasting models at National scale; 2) Earthquake prediction of large events in Italian territory; 3) Modeling the space-time-magnitude evolution of the seismic sequence.
The distinction between points 1 and 2 recalls the distinction between forecast and prediction that is only heuristically defined.Basically, a forecast is the probability of earthquake occurrence for specific time-space-magnitude windows; we have a prediction when the probability of occurrence is high and the space-time window small enough to develop rational strategies for earthquake risk mitigation, such as evacuations or similar actions.In such a case, the output of the model is not a probability, but an «alarm» (or not) for a specific space-time window.The third point involves forecasting models in a time scale of days that is often called «aftershock forecasting».
Notwithstanding the past and recent efforts devoted to these issues, the proposed models still follow different approaches, philosophies, datasets, and physical theories that often led, not surprisingly, to almost opposite results.Besides discussing the main basic ingredients of such contradictory results, that range from the choice of the spatial features (regular grids, seismotectonic zonations, faults dataset) to the choice of the time-dependent model, it is also emphasized that one of the major problems in achieving a general consensus on this important issue is due to the difficulty of testing and comparing the forecasting capability of different models.In practice, several models applied to Italy are not testable at all.It is believed that the possibility to test a model is the only way to maintain this branch of seismology inside a scientific domain (avoiding a pseudo-science dominated by personal beliefs), and consequently to reduce the still large epistemic uncertainty on this field.For this reason, many Italian researchers are actively involved in national and international initiatives, whose final goal is to produce a new generation of earthquake forecasting codes that can run in (almost) real-time and that can be therefore quantitatively evaluated in a forward perspective.Particularly important is the CSEP project (Collaboratory for the Study of Earthquake Predictability) that exports the positive experience made by RELM («Regional Earth-quake Likelihood Models») in California.More details on that can be found in a dedicated special issue of Seismological Research Letters (February 2007), and in the CSEP webpage (http://www.cseptesting.org/).
The main goal of this review is to give a complete picture of the topics reported above.In order to achieve this objective fairly, I do not report the whole range of efforts made in modeling space-time earthquake distribution, but only those addressing the earthquake forecasting issue for Italian territory.

Forecasting large earthquakes in Italy:
was Umbria-Marche a region with high probability of occurrence for destructive earthquakes?
In retrospect, when the Umbria-Marche earthquakes occurred, we cannot say that it they were «unexpected».The most recent hazard map set up before 1997 (Slejko et al., 1998) showed that the epicentral area overlaps a region where significant ground shacing was expected.
As regards the earthquake forecasting, just few years before the Umbria-Marche events, Boschi et al. (1995) published the first probability map on a national scale.They used the «official» (at that time) seismogenic zonation derived on the basis of geological/tectonic information (Scandone, 1992), and applied two different mechanisms of earthquake occurrence at some selected zones.In particular, they considered only zones where sufficient and complete historical catalog information were available (i.e., the zones with highest historical seismic rates), and they used a Poisson process or a renewal process with a lognormal distribution for the interevent times depending on how the past historical earthquake occurrence looks like in that regions.Specifically, the authors analyzed the coefficient of variation α, defined as α = σ/µ, where σ and µ are, respectively, the standard deviation and the mean of the interevent times.A α value less than one implies a tendency towards some regularity in the interevent time, and this was the basis to use a lognormal distribution that mimics a sort of recurrence.A value of α close to 1 was taken as evidence of a pure time-independent stochastic process modeled through a Poisson process.
Here, I do not deepen the physical implications that stand behind the choice to use different processes for different zones, but I just focus the attention on the forecasting made by Boschi et al. (1995).In particular, the model aimed forecasting the occurrence of the M 5.9+ next earthquakes in Italy for different time windows.The Umbria-Marche region has the sixth highest probability of occurrence in the interval 1995-2000 out of 20 regions (see fig. 1).In this respect, the occurrence of the 1997-1998 earth-  space distribution of past events reported in a seismic catalog (Faenza et al., 2003;Akinci et al., 2004).
The most remarkable differences among these models is related to the spatial distribution adopted and the relative assumptions that stand behind it, as well as to the assumptions of the kind of time-dependency.Models belonging to group 1 assume the availability of a complete dataset of seismogenic faults, and the knowledge about the time-magnitude behavior of each of them.In particular, Romeo (2005) assumes that all faults have a characteristic magnitude (Schwartz and Coppersmith, 1984) and a time distribution described by a Brownian Passage Time (Matthews et al., 2002); Pace et al. (2006) also assume a complete fault dataset, but they conjecture faults have different behavior, ranging from characteristic magnitudes and BPT time distribution to magnitudes that follow a Gutenberg-Richter law and a Poisson time distribution.These issues are still under discussion and a general consensus has not been reached yet (see, for instance, Marzocchi, 2007;Peruzza et al., 2007).
Models of groups 2 and 3 relax the strong quakes in this region cannot be considered a «surprise».Nevertheless, the model also assigns a very high probability of occurrence (larger than 0.87 and close to 1) in the time interval 1995-2000 to two areas, Southeastern Sicily and Appennino Abruzzese, where no events occurred.Here, I do not attempt a rigorous evaluation of such a model.Yet, I underline that, regardless of its scientific reliability, this model is one of the few cases where an independent check is possible.As I will discuss later, the possibility to validate a model is a mandatory requirement that is still today rarely fulfilled.
After the Umbria-Marche seismic sequence many more works were devoted to improving our forecasting capability and some works supplied earthquake forecasts for this area.Basically, the attempts can be divided into three groups: 1) estimating the probability on distinct seismic faults (Romeo, 2005;Pace et al., 2006); 2) estimating the probability on seismotectonic areas (Cinti et al., 2004;Gruppo di Lavoro 2004); 3) estimating the probability using the time- assumptions mentioned above.The plausible incompleteness of the fault dataset is accommodated by means of a broader spatial distribution.In other words, the incompleteness of the seismogenic faults dataset is taken into consideration increasing the epistemic uncertainty, i.e., widening the spatial distribution of potential future epicenters.In particular, models of class 2 assume that seismotectonic zonation is a reliable and objective entity that is homogeneous in terms of earthquake occurrence.Despite the geological/tectonic rationales that stand behind them are often quite reasonable, I argue that it is very difficult to verify this assumption; this means that different researchers can suggest different zonations without the possibility to compare them in a quantitative way.Models of group 3 assume that the past seismic catalog is a reliable spatial representation of the whole seismicity in Italy; this may be questionable, when we consider seismic cat-  4).The rationales behind this discrepancy are linked to the assumptions of the models described above; Romeo's (2005) and Pace et al. (2006) models assume that all faults capable of earthquakes with M 5.5+ for this region have been identified and produced recent earthquakes; therefore, the renewal models used lead to low probability for future occurrences.On the other hand, Faenza et al. (2003) and Cinti et al. (2004) model assumes that the high recent seismic activity may have perturbed significantly other faults that have not been recognized yet.For this reason, the probability of future occurrence is high.
The models just described have other important differences, mostly due to the time and magnitude statistical distributions.A full discussion of these important points is beyond the scope of this paper and can be found in the already cited articles.Here, I just only mention to the fact that models of group 1 do not explain the empirical evidence of time clustering observed at different scales (Console et al., 2003a, Faenza et al., 2003;Cinti et al., 2004) and they have been never tested against real data, even retrospectively.On the other hands, some models belonging to classes 2 and 3 have shown to work reasonably well if applied retrospectively to the last 50 years (see Faenza et al., 2003;Cinti et al., 2004).
I argue that one of the major limitations in achieving a general consensus on this important issue is due to the lack of efforts devoted to test and to compare the forecasting capability of different models; in practice, some models are not testable at all.So far, different strategies for validation and comparison among models have been proposed.The most promising are based on probability gain, entropy score (e.g., Harte and Vere-Jones, 2005), and different flavors of likelihood tests (e.g., Schorlemmer et al., 2007 and references therein).A formal validation and comparison among competitive models require the use of pure prospective tests and the definition of a commonly accepted framework in terms of selected area, authorized seismic datasets, etc.In this respect, the Italian territory has been very recently added as a «natural testing region», together with California, New Zealand, and western Pacific for CSEP experiments that share the principles just reported.In these regions, any testable forecasting model can be independently judged and compared to other competitive models in a pure prospective testing.In principle, this allows a drastic reduction of the epistemic uncertainty evaluating what is the best performing forecasting model on future seismicity.

Predicting large earthquakes:
Was a precursory activity recorded before the mainshock(s)?
Earthquake prediction has always been a very hotly debated issue.The lack of convincing evidence of consistent precursory phenomena, together with theoretical considerations, led to widespread skepticism about the real possibility to predict earthquakes (cf.Geller et al., 1997).Some researchers, for instance, described earthquake prediction as the alchemy of our century (Geller, Nature Debate, 1999).Others disagree, believing that we will be able to improve significantly our ability to forecast great seismic events (e.g., Wyss, 1997).Apart from these interesting philosophical speculations, it is unquestionable that a careful reading of the large body of literature on this issue might leave the reader with a hazy picture, typical of an undefined problem.For instance, the terms prediction and forecasting are usually, but not always, taken as synonymous.The meaning of the term precursor appears strongly subjective, embracing a vast number of diverse definitions (e.g., Kagan, 1997).Moreover, unambiguous and generally accepted procedures to set up and test prediction models are not yet available.Also, randomness and unpredictability (of earthquakes) are usually (and erroneously) considered synonymous concepts.The appraisal of such issues is beyond the scope of this review, even though an unbiased outlook of any earthquake prediction attempt deems a deep thought on all of them.Retrospectively, after the Umbria-Marche earthquakes, Di Giovanbattista and Tyupkin (2000) found some evidence of a «precursory» activity before the first shock of September 26.In particular, by applying the RTL prognostic parameter (Sobolev and Tyupkin, 1996), a quiescence stage followed by a period of foreshock activation is observed before the event.The main shock occurred soon after the recovery of the RTL parameter to its normal background level, and an acceleration of seismic energy release, as described by the log-periodic time-tofailure equation (see fig. 5).In comparison to the average background of the previous period, the increase of the area of rupture activated during the twelve months leading up to the Umbria-Marche concentrates in the vicinity of the mainshock.In recent years, the precursory acceleration of seismic energy before a mainshock has been deeply investigated (e.g., King and Bowman, 2003;Mignan et al., 2006).Despite some progress in this field, I note that such models require many empirical parameters that are quite difficult to establish before a mainshock.This limitation has prevented, until now,   1995.5 1996.0 1996.5 1997.0 1997.5 1998.0 the set up of models that are able to work in a forward perspective.After Umbria-Marche earthquakes, the most interesting prediction models reported in the literature for Italian territory are based on the M8 and CN algorithms originally proposed by Keilis-Borok and Kossobokov (1990), andKeilis-Borok et al. (1988).The prediction capability at a worldwide scale for M 7.5+ has been matter of different evaluations (Kossobokov et al., 1999;Marzocchi et al., 2003).The application of the models M8 and CN to the Italian territory have been described by Romashkova et al. (1998), and by Peresan et al. (1999), respectively (see also references therein); since 1 st January 1998, the CN model is applied to a forward earthquake prediction experiment.In this model, the Italian territory is divided in three tectonic macrozones partially overlapped, called NORD, CENTRO and SUD.The area of the zones is ~ 83,000 km 2 for zone NORD, 80,000 km 2 for zone CENTRO, and ~ 62,000 km 2 for zone SUD (see fig. 6).The target earthquakes have M=5.4+ for NORD, and M=5.6+ for CENTRO and SUD.Before discussing the prediction ability of this model, I note that the spatial-temporal window for the predictions made is very large.For instance, in recent years it happened that the code raised an alert for one year in the NORD region that spans from Lazio, to Liguria up to Slovenia (see fig. 6); it is not straightforward to understand what kind of practical actions could be taken differently for what is usually done to mitigate the long-term seismic risk.
That said, I analyzed the predictions made in six years, and I conclude that the prediction capability of the model is rather poor (Marzocchi, 2004).In summary, in the period January 1998 -June 2004, 9 target earthquakes occurred in Italy; 5 of them occurred outside the macrozones considered, therefore «not predicted».The other 4 target earthquakes were «predicted», keeping in alarm 46 out of 76 months an area (NORD) more than twice of Switzerland, and 34 out of 76 an area (CENTRO) twice of Switzerland.It is concluded that the model is not able to make a «better job» than what made by a model based on past seismicity of magnitude comparable to the target earthquakes (see, i.e., Marzocchi et al., 2003).Despite the negative results, I stress again that such a procedure has the undoubted merit of allowing a quantitative validation of the prediction capability.

Tracking the evolution of a seismic sequence: What was and what is now our capability to model the spatio-temporal-magnitude evolution of that sequence?
Forecasting aftershock activity does not have a great relevance for land use planning, but it may have a prominent role for managing at best seismic emergencies after large earthquakes; for instance, linked to an appropriate cost/benefit analysis, it may provide clues for Civil Protection to set practical recipes for inhabitants and workers during seismic crises.Compared to the large earthquake forecasting described in previous sections, here the spacetime windows are much smaller.In forecasting aftershocks, in fact, we focus the attention on small region and forecasting time window of 1 day (see, e.g., Gerstenberger et al., 2005).
Before 1997, the state of the art related to the description of aftershock sequences was mostly based on two empirical laws: the Bath (1965) and the Omori (1894) law.The latter describes the power-law decay in time of the rate of events, and, with only some modification, it is still considered reliable after more than one century of its discover (at the very beginning of seismology).The Bath law states that the differences in magnitudes between mainshocks and their largest detected aftershocks are approximately constant, independent of the magnitudes of mainshocks, and it is about 1.2.This empirical law was particularly important for Umbria-Marche earthquakes.As a matter of fact, it was the main basis to declare that the «worst has passed» after the first large shock of Sept. 26.Nevertheless, the largest event occurred about 12 hours after, killing technicians working inside the damaged Assisi cathedral.Studies carried out after this tragic event showed that Bath law could be due to the choice of different magnitude thresholds for mainshocks and aftershocks (Vere-Jones, 1969;Lombardi, 2002;Console et al., 2003b).Moreover, a real-time application of the Bath law requires a real-time identification of the mainshock; for the time being, a rigorous definition of mainshock is still lacking, while the heuristic recognition of a mainshock as the largest event of a cluster can be made only retrospectively.To summarize, as clearly demonstrated by the Umbria-Marche earthquakes, the Bath law cannot be interpreted in a strict deterministic sense, and, overall, it cannot be applied during emergencies.
Since 1997, a wide variety of models have been proposed to track aftershock sequences, mostly based on statistical and physical rationales.Despite I do not like such a distinction because I do not see a sharp distinction between statistical and physical models, in this review I will use this terminology to classify in a rough way the models so far proposed.It is worth saying that some models, overall the ones based on physical laws, were prominently made to «understand» the physics of the process (no matter what the term «understand» means!) and not directly for forecasting purposes.Both issues are clearly related, but not necessarily completely overlapping.In fact, even though a significant increase of our knowledge of physical processes behind the seismic sequence occurrence may potentially lead to a better real-time forecasting, the limited quality of the real-time data (i.e., we can have in real time hypocentral location, but hardly a reliable focal mechanism solution) may prevent some model to be profitably applied in real-time.
The Umbria-Marche 1997 seismic sequence occurred in a period of time dominated by the scientific concept of fault interaction via Coulomb Failure Function (CFF; i.e., Stein et al., 1992).Remarkably, this sequence represents a clear failure of this model to explain the whole or even the most important part of the spatial distribution of the sequence.In particular, Cocco et al. (2000) showed that CFF cannot explain the occurrence of some M 5.0+ events and, in general, of the high seismicity on the hanging wall of the fault where the «mainshock» occurred at the beginning of the sequence.This discrepancy was explained by introducing a significant influence of fluids flowing across the unclumped faults that would have favored the occurrence of earthquakes (Miller et al., 2004;Antonioli et al., 2005).Specifically, Miller et al. (2004) show that this sequence may have been driven by a fluid pressure pulse generated from the coseismic release of deep source of trapped high-pressure CO 2 .The au-thors find a strong correlation between the high-pressure front and the aftershock hypocentres over a two-week period, using a simple model of nonlinear diffusion (see fig.Today, the hypothesis of the role of fluids is still a matter of discussion and study, mostly because of the scarce direct evidence of fluid flowing across the earth.The only clear evidence is that Umbria-Marche region, as well as a large part of the Italian territory, is characterized by a very high CO 2 flux (Chiodini et al., 2004) that has been suggested to be related to the spatial distribution of the Italian seismicity.As regards the practical usefulness for earthquake forecasting purposes, all of the above mentioned models present two main limitations.First, they require very detailed input that is not yet available in real-time, such as the true fault geometry, a reliable complete description of the focal mechanism, and the heterogeneous slip distribution.Second, these models do not provide any explicit forecast of the size of the triggered earthquakes.In practice, such models could, at best, explain the spatial distribution of the earthquakes and, retrospectively, the succession of the largest earthquakes of the sequence.For the time being, they cannot be used in a real-time forward perspective.
The statistical models, instead, are more flexible and require input data of lower quality.In practice, this means that they can be already used to track in real time a future seismic sequence.The first outstanding application of statistical models to forecast seismic sequences in the Italian territory was proposed by Console et al. (2003a), using the so-called ETAS modeling (Epidemic-Type Aftershock Sequence; Ogata, 1988Ogata, , 1998)).The results show an impressive similarity between the time evolution of the real and expected seismicity.The only remarkable difference is that the expected number of events was a factor of three less than what was really observed.This discrepancy can be due mostly to two different factors: 1) the use of a minimum magnitude for the ETAS modeling equal to the completeness magnitude can induce a similar bias (Sornette and Werner, 2005); 2) a classical stationary ETAS model is not appropriate to explain the time-dependent effects of the fluids migration (see, i.e.Lombardi et al., 2006).One important feature of all ETAS models is that they forecast the magnitude of events using a Gutenberg-Richter law.This means that there is no explicit forecast of the largest aftershocks that are expected mostly in the highest probability regions.

Final remarks
The occurrence of the Umbria-Marche earthquakes in 1997 marked an important change point for Italian research on earthquake forecasting.After that date, many important contributions have been reported, resulting in a wide range of prediction/forecasting models for different space-time-magnitude windows.Nevertheless, up to now, these models have led to contradictory results, and no consensus models have yet been defined.The reduction of such a large epistemic uncertainty is the main objective of very recent international initiatives and projects like CSEP/RELM, NERIES («Network of Research Infrastructures for European Seismology»), and SAFER («Seismic eArly warning For EuoRope»); all of them point towards the formulation of earthquake forecasting codes that can be applied in a forward perspective.This is a pivotal issue, both for practical and scientific purposes, because it allows the codes to be applied in real-time and a comparison among competitive models to be performed.In this respect, I foresee significant steps ahead in next few years that will bring new and robust tools useful for seismic hazard assessment ad to manage at best seismic emergencies.

Fig. 1 .
Fig.1.Boschi et al.'s (1995)  forecasts.A) Map of the regionalization used.B) The probabilities estimated for M 5.9+ at different forecasting time windows (5, 20, and 100 years) are reported; the intervals represent the 68% confidence intervals.The red boxes mark the region where Umbria-Marche earthquakes occurred.

Fig. 2 .
Fig. 2. Pace et al. (2006) forecasts of M 5.5+ earthquakes in the next 50 years.A) Probabilities, recurrence times, and elapsed times for each seismogenic structure are reported.B) Map of the seismogenic structures; stars represent the location of Umbria-Marche earthquakes.In both panels, the red boxes mark the Umbria-Marche region.

Fig. 3 .
Fig. 3.Romeo's (2005) forecasts of M 6.0+ earthquakes in the next 30 years.A) Probabilities for each seismogenic structure.B) Map of seismogenic structures; different colors (gray and black) stand for different reliability (seeRomeo (2005) for further information).In both panels, red boxes mark the Umbria-Marche region.

Fig. 4 .
Fig. 4. Cinti et al. (2004; left panel) and Faenza et al. (2003; right panel) forecasts for M 5.5+ events for the next 10 years.The left panel (A) reports the probability to have at least one event in next 10 years for each seismotectonic zonation.The right panel (B) reports the same quantity per unit of surface (in Km 2 ).The blue box represents the Colfiorito area.For more details see also http://www.bo.ingv.it/~earthquake/ITALY/forecasting/M5.5+/

Fig. 5 .
Fig. 5. Best fit solution and cumulative sum of the square root of the seismic energy released by events before the two main shocks of Sept. 26 2007 marked by black solid circles (Di Giovambattista and Tyupkin, 2000).

Fig. 6 .
Fig. 6.Target areas for earthquake predictions of M 5.4+ (A -North) and 5.6+ (B -Center and C -South Italy) made with CN algorithm in Italy(Peresan et al., 1999).The red box indicates the region where Umbria-Marche earthquakes occurred.
7).The triggering amplitude (10-20 MPa) of the pressure pulse overwhelms the typical (0.1-0.2 MPa) range from stress changes in the usual stress triggering models.Finally, they propose that aftershocks of large earthquakes in such geologic environments may be driven by the coseismic release of trapped, high-pressure fluids propagating through damaged zones created by the mainshock.Along the same line, Antonioli et al. (2005) model the spatial and temporal evolution of seismicity during the 1997 Umbria-Marche seismic sequence in terms of subsequent failures promoted by fluid flow.

Fig. 7 .
Fig. 7. Comparison of aftershock data to CFF changes and pore pressure variations.a) there is no correlation between aftershock occurrence and positive CFF variations; b) aftershocks occur where the pore pressure is higher.