Earthquake forecasting: Statistics and Information

We present an axiomatic approach to earthquake forecasting in terms of multi-component random fields on a lattice. This approach provides a method for constructing point estimates and confidence intervals for conditional probabilities of strong earthquakes under conditions on the levels of precursors. Also, it provides an approach for setting multilevel alarm system and hypothesis testing for binary alarms. We use a method of comparison for different earthquake forecasts in terms of the increase of Shannon information. 'Forecasting' and 'prediction' of earthquakes are equivalent in this approach.


Introduction
The methodology of selecting and processing of relevant information about the future occurrence of potentially damaging earthquakes has reached a reasonable level of maturity over the recent years.However, the problem as a whole still lacks a comprehensive and generally accepted solution.Further efforts for optimization of the methodology of forecasting would be productive and well justified.
One of our objectives here is to address some drawbacks of the approach to earthquake prediction based on pattern recognition theory.It is developed by "Keilis-Borok team" and presented in Keilis-Borok and Soloviev [2003].The pattern recognition is used in forecast with the hope that 'training sets' in which strong earthquakes have higher values of precursors may help to decide whether or not a strong earthquake will happen in the future.As it is demonstrated by below published Krichevets's theorem, this hope has no basis.
A comprehensive review of the modern earthquake forecasting state of knowledge and guidelines for utilization can be found in Jordan et al. [2011].Note that all methods of evaluating the probabilities of earthquakes are based on a combination of geophysical, geological and probabilistic models and considerations.Even the best and very detailed models used in practice are in fact only 'caricatures' of immensely complicated real processes.(For example, the equations of the theory of elasticity are valid only for infinitesimal strains, but used even at the edges of the cracks where strains tend to infinity.It is also clear that we do not know about the inhomogeneity of elastic properties in the real medium, and we have no way of knowing about it with sufficient accuracy.) A mathematical toolkit for earthquake forecasting is well presented by Harte and Vere-Jones [2005].This work is based on the modeling of earthquake sequences in terms of the marked point processes.However, the mathematical technique used is quite sophisticated and does not provide direct practical tools to investigate the relations of the structure of temporal-spatial random fields of precursors to the appearance of strong earthquakes.
The use of the multicomponent lattice models (instead of marked point processes) gives a different novel way of investigating these relations by a more elementary way.Discretization of space and time allows us to separate the problem in question into two separate tasks.The first task is the selection of relevant precursors, i.e., observable and theoretically explained physical and geological facts which are causally related to a high probability of strong earthquakes.Particularly, this task involves the development of models of seismic events and computing probabilities of strong earthquakes in the framework of these models.Such probabilities are used as precursors in the second task.
The second task is the development of methodol-
ogy of working with these precursors in order to extract the maximum information about the probabilities of strong earthquakes.This is the main topic of this paper.
Our approach allows us to obtain the following results: -Estimates of probabilities of strong earthquakes for given values of precursors are calculated in terms of the frequencies of historic data.
-Confidence intervals are also constructed to provide reasonable bounds of precision for point estimates.
-Methods of predictions (i.e., binary alarm announcement [Keilis-Borok and Kossobokov 1990, Keilis-Borok 1996, Holliday et al. 2007]) and forecasting (i.e., calculating probabilities of earthquakes [Kagan and Jackson 2000, Harte and Vere-Jones 2005, WGNCEP 2007, Jordan et al. 2011]) are equivalent in the following sense: the setting of some threshold for the probability of earthquakes allows to update the alarm level.On the other hand, the knowledge of the alarm domain based on historical data allows us to evaluate the probabilities of earthquakes.In a sense, the prediction is equivalent to hypothesis testing as well (see Section 11).
-In our scheme we propose a scalar statistic which is the ratio of the actual increment of information to the maximal possible increment of information.This statistic allows us to linearly order all possible forecasting algorithms.Nowadays the final judgement about the quality of earthquake forecasting algorithms is left to experts.This arrangement puts the problem outside the scope of natural sciences which are trying to avoid subjective judgements.
The foundation of our proposed scheme is the assumption that the seismic process is random and cannot be described by a purely deterministic model.Indeed, if the seismic process is deterministic then the inaccuracy of the forecast could be explained by the non-completeness of our knowledge about the seismic events and non-precision of the available information.This may explain, at least in principle, attacks from the authorities addressed to geophysicists who failed to predict a damaging earthquake.However, these attacks have no grounds if one accepts that the seismic process is random.At the end of the last century (February-April 1999) a group of leading seismologists organized a debate via the web to form a collective opinion of the scientific community on the topic: 'Is the reliable prediction of individual earthquakes a realistic scientific goal?' (see http://www.nature.com/nature/debates/earthquake/).
Despite a considerable divergence in peripheral issues all experts taking part in the debate agreed on the following main principles: -the deterministic prediction of an individual earthquake, within sufficiently narrow limits to allow a planned evacuation programme, is an unrealistic goal; -forecasting of at least some forms of time-dependent seismic hazard can be justified on both physical and observational grounds.The following facts form the basis of our agreement with this point of view.The stringblock Burridge-Knopov model, generally accepted as a mathematic tool to demonstrate the power-like Gutenberg-Richter relationship between the magnitude and the number of earthquakes, involves the generators of chaotic behaviour or dynamic stochasticity.In fact, the nonlinearity makes the seismic processes stochastic: a small change in the shift force may lead to completely different consequences.If the force is below the threshold of static friction the block is immovable, if the force exceeds this threshold it starts moving, producing an avalanche of unpredictable size.
This mechanism is widespread in the Earth.Suppose that the front propagation of the earthquake approaches a region of enhanced strength of the rocks.The earthquake magnitude depends on whether this region will be destroyed or remains intact.In the first case the front moves further on, in the second case the earthquake remains localized.So, if the strength of the rocks is below the threshold the first scenario prevails, if it is above the threshold the second scenario is adapted.The whole situation is usually labelled as a butterfly effect: infinitesimally small changes of strength and stress lead to macroscopic consequences which cannot be predicted because this infinitesimal change is below any precision of the measurement.For these reasons determinism of seismic processes looks more doubtful than stochasticity.
The only comment we would like to contribute to this discussion is that the forecasting algorithms based exclusively on the empirical data without consistent physical models could hardly be effective in practice (see Section 13 for more details).
Finally, note that our approach may be well-applicable for the space-time forecasting of different extremal events outside the scope of earthquake prediction.

Critical analisis of the prediction based on pattern recognition
The main objective of this paper is to present a mathematical formalism of the earthquake forecasting alternative to the approach based on the pattern recognition theory.So, it is appropriate to summarize the latter approach.
An advantage of the formalism of the Keilis-Borok team [Keilis-Borok and Kossobokov 1990, Keilis-Borok 1996, Keilis-Borok and Solovoev 2003] is a separation of the forecast into two separate tasks: the selection of relevant precursors and the technique of their process-ing.On the other hand, some drawbacks of the scheme should be pointed out.

Pattern recognition
The method is relaying on the pattern recognition technique that can provide only empirical but not theoretical foundation.The possibility to base a reliable forcast on the idea of supervised pattern recognition is highly questionable in principle.Indeed, the features selection based on the observations only and not related to the physics of earthquakes is not helpful and any hopes for valid 'training sets' are not grounded.A successful supervised recognition is possible only if the features has proved causal relation with pattern.This principle is illustrated by a simple but important theorem by A.N. Krichevets.
Theorem.Let A be a finite set, a half classifies X as an object of sample B 1 and a half classifies X as an object from B 2 .
Proof.It is easy to define a one-to-one correspondence between classifications.Indeed, if {A 1 , Corollary.A supervised pattern recognition is impossible.After the training procedure the probability to classify correctly a new object is the same as before training.
This result implies that a clairvoyance is involved in a successful forecast: one should know in advance (or guess correctly) a "relevant" set of precursors.

Voting procedure
As a decision rule for a forecast a voting procedure is adopted: an alarm is announced when the number of precursors taking abnormally large values exceeds certain thershold.This procedure is well-justified only if the pattern in question is unique.Generically, this is not true: if two or more possible scenarios are present, different precursors may be strongly affected in different scenarious.In this case a combination of abnormal values for precursors related to different scenarios may produce a false alarm.

Alarm domain
The values of precursors are computed at the point of some lattice, and are assigned a domain surrounding the point.The above scheme employs a circle with a center at a lattice site and a radius proportional to the minimal magnitude of the forecasted earth-quake.However, this choice leads to a contradiction.
Indeed, suppose we announce an alarm for earthquakes with minimal magnitude 6 in a domain A 6 .Note that outside the alarm domain we expected that the earthquake will not occur.But the announcement of the alarm with minimal magnitude 7 is also correct for this point and the alarm should be announced in the domain A 7 as well.By definition A 6 ⊂ A 7 and we expect an earthquake with a magnitude of at least 7 and do not expect an earthquake with à magnitude of at least 6 in the domain A 7 \A 6 .This is absurd.
Futhermore, the alarm area is too large with this approach.

Molchan's error diagram
Below we demonstrate that the procedure of alarm announcement is equivalent to the classical problem of hypothesis testing.So, the standard error diagram provides the dependence of the second type error (probability of missing a target) from the first type error (probability of false alarm).
In line with the general case the well-known Molchan's error diagram [Molchan 1990] displays the dependence of the first kind error probability versus the share of alarms.
H 0 : the distribution of the precursor is equal to its unconditional distribution (earthquake "could either happen or not happen").
H 1 : the distribution of the precursor is equal to its conditional distribution under the condition that the earthquake will happen.
Note that the rejection of hypothesis H 0 leads to absurd results: "earthquake will not happen and will happen".
Of course, Molchan's error diagram can be interpreted not in terms of hypothesis testing.Then there is no contradiction, but the standard error diagram is different from the Molchan's error diagram by a small value proportional to the share of cells occupied by the strong earthquakes (in our approach the spacetime is divided into disjoint cells).

The loss function, risk function and the economic criterion for the quality of the forecast
The Keilis-Borok team refers to the economic loss function for the quantitative evaluation of the quality of forecast.First, we calculate the expected loss correctly.A natural economic measure for a quality of binary forecast is the economic risk function (it should not be confused with seismic risk, there is the concept of the theory of statistical decisions) or damage r related to the earthquakes and the necessary protective measures.In mathematical statistics the risk function is defined as the expectation of the loss function, which is an expected losses for different outcomes.In our case there are two types of losses: damage and expenses related to protection.For each cell of our grid (in our approach the space-time is divided into disjoint cells) the risk may be specied by the formula r = bPr {no earthquake; alarm} + aPr {earthquake; no alarm} + cPr {earthquake; alarm} here the symbol Pr{•} denotes a probability, a stands for the average damage from a seismic event, b stands for the average expenses for protection after a seismic alarm is announced, c stands for the average damage after the alarm, c = a +b − d, where d is the damage prevented by the alarm.The coefficient in front of Pr {no earthquake; no alarm}, obviously, equals 0, because in the absence both of a seismic event and an alarm there is no loss of any kind.Clearly, only the case when d > b is economically justified, i.e., the gain from the prevention measures is positive.Obviously, d should be less than a +b, i.e., an earthquake cannot be profitable.Taking into account that a, b and c depend on the geographical position of the cell, we write the total risk as the summation over all cells in the region of a given forecast.In the simplest case of the absence of the spacial component, when a single cell represents a region of forecast, the expression for the risk is simplied as follows r = amo +bx+cm(1−o), where x is a spatiotemporal share of alarm, m is a share of cells occupied by the earthquakes, o is the share of missed targets.In Molchan's theory related to the one-dimensional case, the third component is missing.
However, the risk function r, which is very useful for economical considerations and as a basis for an administrative decision, could hardly be used as a criteria for the quality of seismic prediction.First of all, it cannot be computed in a consistent way because the coefficients a, b and c are not known in practice, and hence no effective way of its numerical evaluation is known.The computation of these coefficients is a difficult economic problem and goes far beyond of the competence of scientists.On the other hand, the readiness of the authority to commit resources to solving the problem depends on the quality of the geophysical forecast.This situation leads to a vicious circle.
The second drawback of the economic risk as a criterion for the quality of prediction is related to the fact that it depends on many factors which have no relation to geophysics or earthquake prediction.These factors include the density of population, the number and size of industrial enterprises, infrastructure, etc.It also depends on subjective factors such as the willingness of authorities to use resources for prevention of the damage from earthquakes.The natural sciences could hardly accept the criteria for the forecast quality which depend on the type of state organization, priorities of ruling parties, results of the recent elections, etc.
In this paper, we introduce a scalar efficiency (information efficiency) forecast that only depends on information about future earthquakes contained in the precursors.

The measure of an area
In Molchan and Keilis-Borok [2008] the area of the alarm domain is defined in terms of a non-homogeneous measure depending on spacial coordinates.Specifically, the measure of an area is proportional to the number of earthquake epicenters from a sample catalog.This approach is used to eliminate the decrease of the share of alarmed sites x with the extension of the domain when a purely safe and aseismic territory is included into consideration.It would be well-justied if the quantity x could be accepted as an adequate criterion of the quality of forecast in its own right.On the other hand, it can be demonstrated that the our information efficiency converges to a non-zero value 1−o when the number of cells with an alarm is fixed but the total number of cells tends to infinity.An inhomogeneous area of the territory under forecast which is proportional to "activity" does not enable us to calculate the informational efficiency.Moreover, it possesses a number of unnatural features from the view of evaluating economical damage.A seismic event in the territory of low seismicity is more costly because no precautions are taken to prevent the damage of infrastructure.However, in this inhomogeneous area an alarm announced in an aseismic territory will have a smaller contribution than an alarm in a seismically active territory where the losses would be in fact smaller.We conclude that this approach 'hides' the most costly events and does not provide a reasonable estimate of economic damage.
Next, we build a forecasting method free from this kind of flaws.

Events and precursors on the lattice
In order to define explicitly estimates of probabilities of strong earthquakes we discretize the two-dimensional physical space and time, i.e., introduce a partition of three-dimensional space-time into rectangular cells with the space partition in the shape of squares and time partition in the shape of intervals.These cells should not intersect to avoid an ambiguity in computing the frequencies for each cell.In fact, the space cells should not be perfect squares because of the curvature of the Earth's surface, but this may be neg-lected if the region of forecasting is not too large.
So, we obtain a discrete set X K with N = I × J × K points which is defined as follows.Let us select a rectangular domain A of the two-dimensional lattice with I × J points x = (x i , y j ), x i = a × i, i = 1,...,I and y j = a × j, j = 1,...,J, a is the step of the lattice.A cell in X K takes the shape of parallelepiped of height Dt with a square base.Clearly, any point in X K has coordinates (x i , y j , t k ); We say that a seismic event happens if an earthquake with magnitude greater than some pre-selected threshold M 0 is registered.For any cell in our spacetime grid we define an indicator of an event, i.e., a binary function h (i, j, k).This function takes the value 1 if at least one seismic event is registered in a given cell and 0 otherwise.Suppose that for all points (x i , y j , t k ) the value of a vector precursor f (i, j, k) = { f q (i, j, k), q =1,...,Q} is given.The components of the precursor f q (i, j, k), q =1,...,Q are the scalar statistics constructed on the base of our understanding of the phenomena that precede a seismic event.

Mathematical assumptions
A number of basic assumptions form the foundation of the mathematical technique of earthquake forecasting.In the framework of mathematical theory they can be treated as axioms but are, in fact, an idealization and simplication with respect to the description of the real phenomena.Below we summarize the basic assumptions which are routinely used in existing studies of seismicity and algorithms of earthquake forecasting even the authors do not always formulate them explicitly.We accept the following assumptions or axioms of the mathematical theory: (i) The multicomponent random process {h(i, j, k), f (i, j, k)}, describing the joint evolution of the vector precursors and the indicator of seismic events, is stationary.
This assumption provides an opportunity to investigate the intrinsic relations between the precursors and the seismic events based on the historical data.In other words, the experience obtained by analysing the series of events in the past, is applicable to the future as the properties of the process do not depend on time.In reality, this assumption holds only approximately and for a restricted time period.Indeed, plate tectonics destroys the stationarity for a number of reasons including some purely geometrical considerations.For instance, the movements of the plates leads to their collisions, their partial destruction and also changes their shapes.Nevertheless, the seismic process can be treated as a quasi-stationary one for considerable periods of time.At the time when the system changes from one quasi-stationary regime to another (say, nowadays, many researchers speak about the abrupt climate change) the reliability of any prediction including the forecast of seismic events is severely restricted.
(ii) The multicomponent random process {h(i,j,k), Any quantitative characteristic of seismicity more representative than a registration of an individual event is, in fact, the result of averaging over time.For instance, the Gutenberg-Richter law, applied to a given region relates the magnitude with the average number of earthquakes where the averaging is taken over a specific time interval.In order to associate with the time averaging a proper probabilistic characteristic of the process and make a forecast about the future one naturally needs the assumption of ergodicity.This exactly means that any averaging over a time interval [0, T] will converge to the stochastic average when T → ∞.In view of ergodicity one can also construct unbiased and consistent estimates of conditional probabilities of strong earthquakes under the condition that the precursors take values in some intervals.Naturally, these estimates are the frequencies of observed earthquakes, i.e., ratios of the number of cells with seismic events and prescribed values of precursors to the total number of cells with the prescribed values of precursors.(Recall that an unbiased point estimate î of parameter i satisfies the condition E [i] = î and a consistent estimate converges to the true value i when the sample size tends to infinity).
(iii) Any statement about the value of the indicator of a seismic event h(i, j, k) in the cell (i, j, k) or its probability should be based on the values of the precursor f (i, j, k) only.
This assumption means that the precursor in the given cell accumulates all the relevant information about the past and the information about the local properties of the area that may be used for the forecast of the seismic event in this cell.In other words, the best possible precursor is used (which is not always the case in practice).As in the other cases, this assumption is only an approximation to reality, and the quality of a forecast depends on the quality of the selection and accumulation of relevant information in the precursors.
Below we present some corollaries and further specifications.
In practice this assumption means that the forecast for the time t k = t 0 + kDt cannot be affected by the values related to the future time intervals (t k , t k +Dt].In reality all of these events may be dependent, but our forecast does not use the information from the future after t k . (iii-b) The conditional distribution of the random variable h at a given cell depends on the values of the precursors f at this cell and is independent of all other variables.
(iii-c) The conditional probabilities Pr ij {h|u( f )} of the indicator of seismic events h in the cell (i, j, k), under condition u( f ) in this cell do not depend on the position of the cell in space (the time index k related to this probability may be dropped due to the stationarity of the process).
In other words, the rule for computing the conditional probability Pr based on the values of precursors is the same for all cells, and the space indices of probability Pr may be dropped.This condition is widely accepted in constructions of the forecasting algoritms but rarely formulated explicitly.However, the probability of a seismic event depends to a large extent on the local properties of the area.Hence, the quality of the forecasting depends on how adequately these properties are summarized in the precursors.This formalism properly describes the space inhomogenuity of the physical space because the stationary joint distribution of Pr ij {h, f ≤ x} for an arbitrary precursor f depends, in general, on the position of the cell in the domain A. Below we will use the distributions of precursors and indicators of seismic events in domain A that do not depend on the spatial coordinates and have the following form (iii-d) Note that assumption (iii) implies that the conditional probabilities Pr ij {h|u( f )} are computed via the probabilities Pr A {h, f ≤ x} only.
The properties listed above are sufficient to obtain the point estimates for the conditional probabilities of seismic events under conditions formulated in terms of the values of precursors.However, additional assumptions are required for a testing of the forecasting algorithm: (iv) The random variables f (i, j, k) are conditionally independent under condition that h(i, j, k) = 1.
Again, these conditions are not exactly true, however they may be treated as a reasonable approximation to reality.Indeed, if the threshold M 0 is sufficiently high then the strong earthquakes may be treated as rare events, and the cells where they are observed are far apart with a high probability.Any two events re-lated to cells separated by the time intervals t are asymptotically independent as t → ∞ because the seismic process has decaying correlations (the mixing property in the language of random processes).The loss of dependence (or decaying memory) is related to the physical phemonema such as healing of the defects in the rocks, relaxation of strength due to viscosity, etc.As usual in physical theories, we accept an idealized model of the real phenomena applying this asymptotic property for large but finite intervals between localizations of seismic events.The independence of strong earthquakes is not a new assumption, in the case of continuous space-time it is equivalent to the assumption that the locations of these events form a Poisson random field.(Note that the distribution of strong earthquake should be homogeneous in space a priori, because by condition there is no information about the heterogeneity outside precursors.The information about the location of epicenters should be contained in the precursors) The Poisson hypothesis is used in many papers, see, e.g.Harte and Vere-Jones [2005].It is very natural for the analysis of the tails of the Gutenberg-Richter law for large magnitudes [Pisarenko et al. 2008].Summing up, the development of the strict mathematical theory of earthquake forecasting does not require any additional assumption except those routinely accepted in the existing algorithms but usually not formulated explicitly.

The standard form for precursors
The correct solution of the forecasting problem given the values of precursors f (i, j, k) is provided by the estimate of conditional probability Pr ij {h|f (i, j, k)} of the indicator of a seismic event in the cell (i, j, k).In practice this solution may be difficult to obtain because the number of events in a catalog might be not sufficient.
Indeed, the range of value of a scalar precursor is usually divided into a number M of intervals, and only a few events are registered for any such interval.For a Q-dimensional precursor the number of Q-dimensional rectangles, covering the range, is already M Q , and majority of them contains zero events.Only a small number of such rectangles contains one or more events, that is the precision of such an estimate of conditional probability is usually too low to have a practical value.
For this reason one constructs a new scalar precursor in the form of the scalar function of component of the vector precursor, and optimizes its predictive power.This approach leads to additional complication as the units of measurement and the physical sense of different components of a precursor are substantially different.In order to overcome this problem we use some transformation to reduce all the components of the precursor to a standard form with the same sense and range of values.

, , Pr
Let us transform all the precursors f q (i, j, k), q =1, ...,Q to variables with the values in [0, 1] providing estimates of conditional probabilities.So, after some transformation F we obtain an estimate of Pr{h|u( f (i, j, k)= 1}, where u is a characteristic function of some interval B, i.e., the probability of event h(i, j, k) = 1 under condition that this precursor takes the value f (i, j, k)∈ B.
The transformation F of a scalar precursor f (i, j, k) is defined as follows.Fix an arbitrary small number f.Let L be a number of cells (i, j, k) such that h(i, j, k) = 1, and Z l , l = 1, ...,L, be the order statistics, i.e., the values f (i, j, k) in these cells listed in non-decreasing order.Define a new sequence z m , m = 1, ...,M, from the ordered statistics Z l by the following recursion: z 0 =−∞, z m is defined as the first point in the sequence Z l , such that z m − z m−1 ≥ f.Next, construct the sequence z m * = z m + (z m+1 − z m )/2, m = 1, ...,M − 1, and add the auxiliary elements z 0 =−∞, z M = 1.Define also a sequence n m , m = 1, ...,M, where n m equals to the number of values in the sequence Z l , such that The transformation F is defined as follows This definition implies that transformation F replaces the value of a precursor for the frequency, i.e., the ratio of a number of cells containing a seismic event and the values of a precursor from [z m * −1 , z m * ) to the number of cells with the value of precursor in this range.These frequencies are the natural estimates of conditional probabilities Pr {h(i, j, k) = 1|z m * −1 ≤ f (i, j,k) < z m *}, m = 1, ...,M, computed with respect to stationary distribution P A (x): (This conditional probability can be written as Pr{h=1|u( f )}, where u is the characteristic function of interval [z m * −1 , z m * ).The function g has a stepwise shape, and the length of the step is bounded from below by f.It can be checked that there exists the limit g = lim -- The estimates of conditional probabilities in terms of the function g are quite rough because typically the numbers n m are of the order 1.As a final result we will present below more sharp but less detailed estimates of conditional probabilities and confidence intervals for them.

Combinations of precursors
There are many ways to construct a single scalar precursor based on the vector precursor (F f q , q =1,..., Q).Each such construction inevitably contains a number of parameters or degrees of freedom.These parameters (including the parameters used for construction of the precursors themselves) should be selected in a way to optimize the predictive power of the forecasting algorithm.The optimization procedure will be presented below, its goal is to adapt the parameters of precursors to a given catalog of earthquakes, that is to obtain the best possible retrospective forecast.However, this adaptation procedure creates a "ghost" information related with the specic features of the given catalog but not present in physical propertities of real seismicity.This "ghost" information will not be reproduced if the algorithm is applied to another catalog of earthquakes.It is necessary to increase the volume of the catalog and to reduce the number of free parameters to get rid of this "ghost" information.Clearly, the first goal requires the considerable increase of the observation period and may be achieved in the remote future only.So, one concentrates on the reduction of numbers of degrees of freedom.The simplest variant including Q −1 parameters is the linear combination As a strictly monotonic function of a precursor is a precursor itself the log-linear combination is an equally suitable choice Here c q , q = 1, ..., Q −1 are free parameters.The result of the procedure has the form g = F f * .

Alarm levels, point and interval estimations
In view of Equation ( 5.3) the precursor g is the set of estimates for probabilities (5.1) (5.2) (6.1) (6.2) (5.3) (5.4) Its serious drawback is that typically {z * l−1 ≤ f (i; j; k) < z * l } correspond to single events, and therefore the precision of these estimates is very low (the confidence intervals discussed below may be taken as a convenient measure of precision).In order to increase the precision it is recommended to use the larger groups of cells containing a larger number of events, that is a more coarse covering of the space where the precursor takes its values.In a sense, the precision of the estimation and the resolving power are related by a kind of "uncertainty principle": the more precise estimate one wants to get the more coarse is the time-space range of their values and vice versa.
We adapt the following approach in order to achieve a reasonable compromise.
(1) For fixed thresholds a s , s = 1, ..., S+1, a 1 = 1, a s < a s +1; a S+1 = 0, we define S possible alarm levels a s+1 ≤ g(i; j; k) < a s and subsets X s , s = 1, ..., S, of the set X K corresponding to alarm levels, i.e., X s is a set of cells of X K , such that a s+1 ≤ g(i; j; k) < a s .
There are different ways to choose the number S of alarm levels and the thresholds a s , s = 2, ..., S. Say, fix S = 5, and select a s = 10 −a(s−1) .This is a natural choice of the alarm level because at a = 1 it corresponds to decimal places of the estimate of the conditional probability given by the precursor.The problem with S = 2, i.e. two-level alarm, may be reduced to the hypothesis testing and is discussed in more details below.
(2) Compute the point estimates of probabilities obtained via the distribution P X (x) of a precursor g in the same way as in Equation (4.4).The property (iv) implies that for any domain X s the binary random variables h(i, j, k) are independent and identically distributed, i.e. and the unbiased estimate of p s takes the form where n s stands for the number of cells in domain X s , and by m s we denote the number of cells in X s containing seismic events.
(3) The random variable m s takes integer values between 0 and n s .The probabilities of these values are computed via the well-known Bernoulli formula Let us specify the confidence interval covering the unknown parameter p s with the confidence level c.In view of the integral Mouvre-Laplace theorem for n s large enough the statistics is approximately Gaussian N(0, 1) with zero mean and unit variance.Note that the values n s increase with time.Omitting straightforward calculations and replacing the parameter p s by its estimate i s we obtain that i s − < i s < +i s + ; where and t c is the solution of equation Here U stands for the standard Gaussian distribution function.
(4) As a result of these considerations we introduce 'the precursor of alarms' which indicates the alarm level: R( f (i, j, k)) = s (i, j, k).It will be used for calculations of point estimates and the confidence interval in the form This result will be used for a prospective forecasting procedure.

The information gain and the precursor quality
The construction of a 'combined' precursor R involves parameters from formula (6.1) or (6.2) as well as parameters which appear in the definition of each individual precursors f q .It is natural to optimize the forecasting algorithm in such a way that the information gain related to the seismic events is maximal.In a onedimensional case the information gain as a measure of the forecast efficiency was first intoduced by Vere-Jones [1998].Here we exploit a modification of his ideas in the case of multi-dimensional space-time process.Remind the notions of the entropy and information.Putting aside the mathematical subtlety (see Kelbert and Suhov [2013] for details) we follow below an intuitive approach by Prohorov and Rozanov [1969].The information containing in a given text is, basically, the length of the shortest compression of this text without the loss of its content.The smallest length S of the sequence of digits 0 and 1 (in a binary code) for counting N different objects satises the relations 0 , , s i j k s i j the quantity S ~log 2 N characterizes the shortest length of coding the numbers of N objects.
Consider an experiment that can produce one of N non-intersecting events A 1 , ..., A N with probabilities q 1 , ..., q N , respectively, q 1 + ...+ q N = 1.A message informing about the outcomes of n such independent identical experiments may look as a sequence (A i 1 , ..., A i n ), where A i k is the outcome of the experiment k.But for long enough series of observations the frequency n i /n of event A i is very close to its probability q i .It means that in our message (A i 1 , ..., A i n ) the event A i appears n i times.The number of such messages is By the Stirling formula the length of the shortest coding of these messages The quantity S n measures the uncertainty of the given experiment before its start, in our case we are looking for one of the possible outcomes of n independent trials.The specific measure of uncertainty for one trial is known as Shannon's entropy of distribution q 1 , ..., q N (in physical literature it is also known as a measure of chaos or disorder).After one trial the uncertainty about the future outcomes decreases by the value S =S n −S n−1 , this decrement equals to the information gain I = S, obtained as a result of single trial.The quantity is the (unconditional) entropy of distribution for the indicator of the seismic event h in a space-time cell in the absence of any precursors.The conditional entropy S(h|a s+1 ≤ g < a s ) under condition that in the cell (i, j, k) the alarm level s is set up equals The expected conditional entropy S R (h) of indicator of seismic events where the averaging in performed by the distribution of precursors R takes the form We conclude that the knowledge of the precursor values helps to reduce the uncertainty about the future experiment by S(h)−S R (h) which is precisely information I(R,h) obtained from the precursor.Taking into account Equations (8.1), (8.2) and the fact that we specify the information gain as By analogy with the one-dimensional case [Kolmogorov 1965] the quantity I(R,h) may be called the mutual information about the random field h that may be obtained from observations of random field R. It is known that the information I(R,h) is non-negative and equals to 0 if and only if the random fields h and R are independent.This mutual information I(R,h) takes its maximal value S(h) in an idealized case of absolutely exact forecast.The mutual information is the information that the distributions of precursors contribute to that of the indicator of a seismic event.For this reason it may be considered as an adequate scalar estimate for the quality of the forecast.
The quantity I(R,h) depends on the cell size, i.e., on the space discretization length a and time interval Dt.We need a formal test to compare precursors defined for different size of the discretization cells.For this aim let us introduce the so-called 'efficiency' of precursors as the ratio of information gains This efficiency varies between 0 and 1 and serves as a natural estimate of information quality of precursors.It allows comparing different forecasting algorithms and select the best one.
A natural estimate of S(h) based on Equations (5.1) and (5.2) is defined as follows Taking into account Equation (8.1) and using an estimate of P A (a s+1 ≤ g < a s ) in the form of ratio x s = ns -N , we construct an estimate of I(R,h) as follows Remark.It seems reasonable to introduce a penalty related to the number of superficious parameters in evaluating the quality of forecast pointing to the natural analogy with the Akaike test [Akaike 1974] and similar methods in information theory.In our context the main parameter of importance is r(R; h) and its limit as t →∞.This quantity does not involve the number of parameters directly.Probably, the rate of convergence depends on the number of parameters but this dependence is not studied yet.

The forecasting procedure
The number of time intervals, i.e., the number of observations N used in the construction of estimates increases with the growth of observation time.So, the computation procedure requires constant innovations.On the other hand some computation time is required to 'adapt' the model parameters to the updated information about seismic events via an iterative procedure.For these reasons we propose the following forecasting algorithm.
(1) Given initial parameter values at the moment t K−1 = t 0 +(K−2)Dt we optimize them to obtain the maximum of efficiencys r(R; h) of a precursor in domain X K−1 .For this aim the Monte-Carlo methods is helpful: one perturbs the current values of parameters randomly and adapts the new values if the efficiency increases.The process continues before the value of efficiency stabilized, this may give a local maximum, so the precedure is repeated for a sufficient number of times.The choice of the initial value on the first step of the optimization procedure is somewhat arbitrary but a reasonable iteration procedure usually leads to consistent results.The optimization procedure takes the period of time t K−1 < t ≤ t K .
(2) Next, we construct the forecast in the following way.At the moment t K the values of precursor ĝ in each cell (i, j, K + 1) is computed with optimized parameters.Based on these parameter values the alarm levels, the point estimates and confidence intervals are computed in each cell as well as the values of efficiency of precursors.
(3) The estimates of stationary probabilities of seismic events in the cell î(i, j) are defined as follows: they can be used for creation of maps of seismic danger (but not "hazard", this term is reserved) in the region.

Retrospective and prospective informativities
The efficiency of precursor which is achieved as a result of parameters optimization could be considered as retrospective as it is constructed by the precursors adaptation to the historical catalogs of seismic events.The prospective efficiency for the space-time domain containing the cell in the 'future' is based on the forecast.It is computed via formulas (8.3), (8.4), (8.5) with the only difference to the retrospective efficiency, however, approaches this value with time.In principle, the prospective efficiency is an ultimate criterion of the precursors quality and the retrospective efficiency could serve only for the preliminary selection of precursors and their adaptation to the past history of seismic events.

Testing of the forecasting algorithm
The efficiency of precursor could be computed exactly only in an idealized case of infinite observation time.However, its estimate may be obtained based on the observation over a finite time interval.So, if an estimate produces a non-zero value not necessarily the real effect is present.It may be simply a random fluctuation even if the precursor provides no information about the future earthquake.For this reason we would like to check the hypothesis H 0 about the independence of a precursor and an event indicator with a reasonable level of confidence.In case the hypothesis is rejected one has additional assurance that the forecasting is real, not just a "ghost".So, consider the distributions and The function P A (P A −1 (y)) = y of variable y = P A (x) provides an uniform distribution F * (y) = Pr {p ≤ y} of some random variable p on [0,1].Next, consider a distribution function G(y) = P' A (P A −1 (y)) on [0,1], and use a parametric representation for abscissa P A (x) and ordinate P' A (x).If random fields g and h are independent the distribution functions P A (x) = P' A (x) and G(y) are uniform.So, the hypothesis about the absence of forecasting, i.e., about the independence of g and h is equivalent to the hypothesis H 0 that the distribution G(y) is uniform.
The empirical distribution G L (y) related to G(y) is defined as follows.Denote by u l , l = 1, ..., L the values of the function g (i, j, k) sorted in the nondecreasing order and belonging to the cells where h (i, j, k) = 1.Let n l be the numbers of cells such that h (i, j, k) = 1, g (i, j, k) = u l .Denote by m(u l ) the numbers of cells such that g (i, j, k) < u l , and define the empirical distribution G L (y) as a stepwise function with G L (0) = 0 and positive jumps of the size The well known methods of hypothesis testing requires that the function G L (y) has the same shape as for independent trials, i.e., random variables u l , l = 1, ..., L are independent in view of axiom (iv).Naturally, we accept the precursors such that the hypothesis H 0 is rejected with the reasonable level of confidence.(Remind, that the hypothesis is accepted if and only if its logical negation could be rejected based on the available observations.The fact that the hypothesis cannot be rejected does not mean at all that it should be accepted, it only means that the available observations do not contradict this hypothesis.Say, the well-known fact that "The Sun rises in the East" does not contradict to our hypothesis, however it may not be considered as a ground for its acceptance.)For large values of L the Kolmogorov statistics [Kolmogorov 1933a] is helpful for this aim with an asymptotic distribution or Smirnov's statistics [Smirnov 1939] with asymptotic distribution The asymptotic expressions for these statistics can be used for L > 20 [Bolshev and Smirnov 1965].

The binary alarm and the hypothesis testing
The prediction is the form of forecast when an alarm is announced in a given cell without a preliminary evaluation of the probability of a seismic event.
In this case we can estimate the probabilities of events, too.(If the alarm is announced in an arbitrary domain X we set up an alarm if at least half of the cell of our model is occupied by alarm.) Let M be the number of cells in X which are in the state of alarm, M 0 be the number of cells where the seismic event is present but no alarm was announced (the number of 'missed targets').Denote by x = M -N the share of the cells with alarm announced, m = L -N the share of the cells with seismic events, and o = M 0 -M the share of missed targets.Let a random variable h (i, j, k) equal 1 if an alarm is announced in the cell (i, j, k) and 0 otherwise.Obviously, the estimate of conditional probability Pr{h (i, j, k) = 1|h (i; j; k) = 1} of the seismic event under the condition of alarm is m ; and the esti-mate of conditional probability Pr{h (i, j, k) = 1|h (i; j; k) = 0} of the seismic event under the condition of no alarm is mo ---1−x .If the alarm is announced according to the procedure described in Section 5 the threshold a 1 specifying the acceptable domain of values for g (i, j, k) should be treated as a free parameter and selected by maximizing the information efficiency r (h; h).The estimate of information increase for given values of x and o equals The value of h (i; j; k) characterizes the results of checking two mutually exclusive simple hypothesis: H 0 : the distribution of f * (i, j, k) has the form implying 'no seismic events', or H 1 : the distribution of f * (i, j, k)) has the form implying the presence of seismic event.
Statistics for checking of these hypothesis is the precursor g (i, j, k), and the critical domain for H 0 has the form {g (i, j, k) ≥ a 1 }.(If usual method of alarm announcement is used the relevant precursor plays the role of statistics and the critical domain is defined by the rule of the alarm announcement).The probability of a first type error it is estimated as . The probability of second type error is estimated as o.(Note that due to condition (iii) any test used for checking these hypothesis should not depend on the coordinates of the cell).
The Neyman-Pearson theory allows to define the domain of images of all possible criteria: in coordinates (a, b) it is a convex domain with a boundary C which corresponds to the set of uniformly most powerful tests.This family may be defined in terms of the likelihood ratio K (x) =

EARTHQUAKE FORECASTING: STATISTICS AND INFORMATION
tions P 1 A (x) and P 0 A (x) have the densities p 1 A (x) and p 0 A (x).
where ~denotes the threshold.In Gercsik and Kelbert [2004] we demonstrated that among all the tests with the images on the boundary C there exists three different best tests.Here the term "best" may be understood in three different senses, i.e., maximizing the variational, correlational and informational efficiency.The most relevant criteria is the informational efficiency r(h; h).

The choice of precursors
We use the term 'empiric precursor of an earthquake' for any observable characterisric derived from the catalog only, which provides for this catalog a reasonable retrospective forecast of seismic events, and which is not derived from basic physical concepts of seismicity (say, the periods of relatice calm, deviation of some basic characteristic from a long-time average , etc).In contrast, the physical precursors are derived from some physical processes and characterize physical quantities (stress fields, strength, concentration of cracks, etc.) or well-defined physical processes (i.e., phase transitions, cracks propagations, etc.)In the meteorological forecast the danger of using empirical precursors was highlighted by A. Kolmogorov in 1933[Kolmogorov 1933b].From that time the meteorological forecast relies on the physical precursors which are theoretically justified by the models of atmospheric dynamics.Below we will present Kolmogorov's argument adapted to the case of seismic forecast.This demonstrates that the purely empirical precursors work well only for the given catalog from which they are derived.However, their efficiency deteriorates drastically when they are applied to any other independent catalog.
Consider a group of k empirical precursors used for a forecast and selected from a set of n such groups.According to Kolmogorov's remark the number k is typically rather small.This is related to the fact that a number of strong earthquakes in catalog is unlikely to exceed a few dozen.As the values of precursors are random there exists a small probability p that the efficiency of the forecast exceeds the given threshold C. Then the probability of event r(R, h) ≤ C equals 1 − p, and the probability of event r(R, h) > C for at least one collection of precursors equals P = 1 − (1 − p) n and tends to 1 as n → ∞. (According to Kolmogorov some arbitrariness of the assumption of independence is compensated by the large number of collections.) Summing up, if the number of groups is large enough with probability close to 1 it is possible to find a group giving an effective retrospective forecast for a given catalog.In practice this is always the case as the number of empirical precursors could be increased indefinitely by variation of real parameters used in their construction.It is important to note that for such a group, which is highly efficient for the initial catalog, the probability that the efficiency is greater than C is still equal to p for any other catalog.In other words: the larger the number of the groups of empirical precursors the less reliable the forecast is.So, the collection of a large list of the empirical precursors is counter-productive.Much more reliable are the physical precursors intrinsicly connected with the physical processes which preserve their values with the change of sample.The probability to find such a set of precursors by pure empirical choice is negligible because they are very rare in the immense collection of all possible precursors.

Demonstration of algorithm
A preliminary version of the forecast algorithm described above was used in the paper [Ghertzik 2008] for Califorrernia and the Sumatra-Andaman earthquake region.These computations serve as a demonstration of the efficiency of the method but their actual results should be taken with a pinch of salt because the selection of precursors does not appear well-justied from the modern point of view: the number of free parameters (34) to be adapted in the precursor "stress indicator" is too large.For this reason the method was not fully exploited.Nowadays, the progress is resumed with a particular emphasis on the selection of a new set of precursors strongly related to the critical phenomena in geophysics.
The following selection of precursors was used: (1) The stress indicator is a scalar predictor characterizing the local level of the stress field.It is conceivable as a stress tensor component.The predictor increases due to tectonic movements, decreases with viscosity and with increasing fracturing, and is affected by earthquake sources considered as cracks.
(2) The damage indicator is a predictor characterizing the stress-induced microfracturing level.The predictor increases with the value of applied stresses and decreases as cracks are healed.
(3) The third predictor is the concentration Zhurkov-Sobolev criterion of seismically active faults based on the calculation of the parameter concentration (density) seismogenic ruptures and modified in such a way that it can be calculated for a local cell.
(4) The far-order predictor is the integral of the radial correlation function for seismic energy released in a spatiotemporal cell.Its construction is based on the idea of the strong earthquake as a phase transition be-, , , , tween consolidation and fracture states of a medium.The phase transition is preceded by the appearance of a far order, an abrupt increase in the correlation radius for the energy of shocks released in a cell.
For detailed description of precursors we refer to original, for it is only an illustration but not the topic of this work.It is useful, however, to summarize the main results of this paper.
California region.The catalog Global Hypocenter Data Base CD-ROM NEIC/USGS, Denver, CO, 1989, together with data from the site NEIC/USGS PDE (ftp://hazard.cr.usgs.gov)for earthquakes with magnitudes M ≥ 4.0 with epicenters between 113°-129° of western longitude and 31°-43° of northern latitude was used for parameter adaptation.The initial time t 0 was selected by subtracting from the time of actual computation, 08.03.2006, an integer number of halfyear intervals such that t 0 fits the first half of the year 1936.(The final time 08.03.2006 could be considered as an initial moment for constructing half year forecast forward up to the date of the latest earthquake available in the catalog).During the computation the time interval from the first half of 1936 to the first half of 1976 was used for relaxation of the zero initial data used for precursors.After this date the catalog for the earthquakes with magnitude M ≥ 4.0 was used to estimate the probabilities of strong earthquakes with magnitude M ≥ 6.0 up to the moment of actual forecast.Note that the adapted restriction to include into consideration only earthquakes with magnitudes M ≥ 4.0 is a severe restriction.It decreases the precision of precursor computation and therefore, if a prediction is successful, increases the degree of condence in the predictor choice.We choose a = 150 km as a size of the spacial lattice, and Dt = 6 months as a time-step.Retrospective forecast was performed with 5 alarm levels defined by the thresholds a s = 10 −a(5−s) , s = 1, ..., 4 and a = 0.75.(Due to too short time step no alarms was registered on the lowest level when parameter a = 1 was selected).In order to reduce the influence of the boundary conditions the large square covering all the seismic events in the catalog used in the computations was reduced by two layers of elementary cells from each boundary.As a result of optimization the forecast information efficiency of 0.526 was achieved, i.e., the forecasting algorithm applied to the given catalog extracted from it about 53% of all available information about seismic events.It seems that this result could be only partially explained by a lucky selection of precursors: another contributor to the high efficiency of the algorithm is the adaptation of the parameters to the features of the specific catalog.The influence of this artificial information may be reduced only with the increase of the observation interval.Accepting the rule of binary alarm announcement in the cells from group 1 and 2 from 5 levels possible one obtains that the space-time share of alarmed cells is 3.4% and the share of missed targets is 18.2%.This result is comparable with the best forecasts available in the literature and obtained by other methods (in the cases when the quantitative parameters of algorithms are presented in the publications).When the forecast was constructed in the future we obtained that the estimate of probability of a strong earthquake anywhere in the area under study is 0.174, and the maximal point estimate of an event in any individual cell is 0.071.As a whole the seismic situation in California did not look too alarming.Indeed, there were no strong earthquakes in California in the next half-year.Table 1.Number of events and alarm cells for each level of the retrospective analysis for California.
Table 1 presents the number of events and the total number of alarm cells for each level of the retrospective analysis, as well as their ratios, which are wellgrounded point estimates of the probability that at least one earthquake with a magnitude M ≥ 6.0 will occur in a cell at a given level of alarm.The 5% condence intervals are also given in the table for each of these probabilities.The true values of event probabilities lie within these intervals with a 95% probability (condence coefficient).In all, 5856 cells participated in the prediction, and 44 of these cells (0.75% of all cells) contained events.
SAE region.We have conducted a retrospective forecast of strong earthquakes with magnitudes M ≥ 7.0 for the whole region where the Sumatra-Andaman earthquake (SAE) happened on 26.12.2004 with magnitude M = 9.3.The catalog Global Hypocenter Data Base CD-ROM NEIC/USGS, Denver, CO, 1989, together with the data from the site NEIC/USGS PDE (ftp://hazard.cr.usgs.gov) for earthquakes with magnitudes M ≥ 5.5 with epicenters between 84.3°-128° of eastern longitude and 20°-26° of northern latitude was used for the parameters adaptation.The initial moment of time t 0 was selected by subtracting from the time of actual computation, 10.11.2004, an integer number of half-year intervals such that t 0 fits the first half of the year 1936.(The final time was selected in such a way that the next half-year period covers SAE and its powerful aftershock).During the computation the time interval from the first half of 1936 to the first half of 1976 was used for relaxation of the zero initial data used for precursors.After this date the catalog was used to estimate the probabilities of strong earthquakes with magnitude M ≥ 7.0 up to the moment of actual forecast.(For magnitude M ≥ 7.5 the number of seismic events was not sufficient for reliable forecast because the 5%-confi-dence intervals strongly overlapped).In this case the restriction to include into consideration only earthquakes with magnitudes M ≥ 5.5 was adapted.As before, it decreases the precision of precursor computation and therefore, if the prediction is successful, increases the degree of confidence to the predictor choice.We selected the size of the spacial grid as a = 400 km and the size of time-step Dt = 6 months.Retrospective analysis was conducted following the same scheme as in the previous case.In order to reduce the influence of the boundary conditions the large square covering all the seismic events in the catalog used in the computations was reduced by two layers of elementary cells from each boundary.As a result of optimization the forecast information efficiency was 0.549, i.e., the forecasting algorithm extracted around 55% of all available information about seismic events when applied to the given catalog.In the case of binary alarm announcement the space-time share of alarmed cells was 3.1% and the share of missed targets was 8.3%.This result is comparable with the best forecasts available in the literature and obtained by the other methods (in the cases when the quantitative parameters of algorithms are presented in the publications).In the case of forward forecast the two most powerful earthquakes, i.e., SAE and its major aftershock, happened in the alarm zone of the second level, and two other events with smaller magnitudes in the zone of alarm level 4.
The results of the half-year forward prediction beginning from November 10, 2004, and M ≥ 7.0 earthquakes that occurred in this period are plotted in Figure 2. As seen from the figure, the two strongest events (the SAE and its strongest aftershock) occurred in the level-2 alarm zone, while two weaker events occurred in the level-4 alarm zone.It is noteworthy that, if monitoring in this region were conducted according to the proposed scheme, a nine-cell square cluster containing one cell of the level-1 alarm would be regarded as hazardous, giving every reason for anxiety in relation to the impending SAE.
Table 2 presents the number of events and the total number of alarm cells for each level of the retrospective analysis, as well as their ratios, which are wellgrounded point estimates of the probability that at least one earthquake with a magnitude M 7.0 will occur in a cell at a given level of alarm.The 5% condence intervals are also given in the table for each of these probabilities.In all, 8352 cells participated in the forecasting and 36 of these cells (0.43% of all cells) contained events.
Figure 3 plots the dependence of the fraction of target misses on the spatiotemporal fraction of alarms for points at which the boundaries between alarm levels divide the set of all cells into two classes: union of more anxiety cells and union of less anxiety cells.It is an estimate for Molchan's error diagram.

Conclusions
This article demonstrates that the prediction method based on the pattern recognition does not have sufficient foundation and contains a number of errors.We propoused here the alternative approach that is free of this errors and enhances the computing capability of the forecasting.Discretization of spacetime, using elements of the theory of multicomponent random processes and axiomatization of the problem allows to drastically simplify the procedure of earthquake forecasting and to obtain the following results: (1) The method of announcement to declare a multi-level alarm was found.
(2) The method of calculating the point estimate of probability of an event for each alarm level was found.
(3) The method of calculating interval estimates of event probability for each level of anxiety was found.
(4) The equivalence of the announcement of a two-phase alarm and solving the problem of testing hypotheses was shown.
(5) The method of testing an algorthim of forecast was found.
(6) The scalar criterion for comparison of different algorithms of forecast was found.
(7) The equivalence of the method of the announcement of alarms and method assessing the probability of event was shown.
It appears to the authors, that a satisfactory mathematical formalism of earthquake forecasting based on arbitrary set of precursors is fully fledged.The most important result we believe is the discovery of the scalar quality criterion of forecast algorithms.The criterion allows to find out which one is the best from existing algorithms and to focus on finding the most effective precursors.
and use the ratios as estimate of unconditional probability of a seismic event in a given cell

Figure 1 .
Figure 1.Mapped areas of five levels of half-year forward forecasting for California on March 8, 2006.

Figure 2 .
Figure 2. Mapped areas of five levels of half-year forward forecasting for the SAE region on November 10, 2004.