Over 40 years ago (starting in 1962) I became involved in the registration and analysis of seismoacoustic impulses (microearthquakes with the size of the focal zone less than one meter) when conducting rockburst prediction experiments in a Ukrainian coal mine. I realized that the event occurrence was random and needed to be analyzed using statistical methods. Because standard statistical techniques were not powerful enough to analyze the spatio-temporal patterns of an earthquake occurrence, I began to study the theory of stochastic processes. Continuous processes cannot model the occurrence of almost instantaneous events properly, thus stochastic point processes should serve as the basic model for description of seismicity. After studying the problem extensively, I wrote a paper on the application of multidimensional point processes to analyze seismic occurrences. As I was finishing the manuscript (Kagan, 1973), I learned that David Vere-Jones, a mathematician at Victoria University in New Zealand, had published a paper (1970) in which he detailed an earthquake occurrence model based on the theory of one-dimensional point processes. I wrote to him and soon we met in Moscow, Russia and started our collaboration which continues to this day (Kagan and Vere-Jones, 1996).

In 1974 I came to UCLA and, working with Leon Knopoff, started to apply my model to seismicity analysis. Very soon it became clear to me that the size-time-space statistical distributions describing earthquake occurrences follow power-laws. At that time I saw a short article in Science News magazine about a new book by Benoit Mandelbrot (1977), and I immediately understood the importance of fractal concepts for the seismicity analysis. In 1980-82, with Knopoff and also by myself, I wrote several papers which presented a scale-invariant, fractal model of the earthquake process (Kagan and Knopoff, 1980, 1981; Kagan, 1981a,b; 1982). It appears now that these papers were well ahead of their time, since only recently they have attracted the attention of geophysicists and statisticians. These publications were also of interest to Mandelbrot and later to Per Bak (who passed away last year), both of whom I also corresponded with and discussed my results since the 1980s.

Since the early 1980s I have been trying to develop a theoretical model of earthquake occurrences and to improve statistical methods for their analysis. These efforts evolved in several directions; some of the topics are theoretical and some have a clear practical application. In most of these papers I use mathematical methods not employed previously in seismology. These attempts caught the attention of seismologists and statisticians which, I think, explains the relatively high level of citation. Among these problems I would mention the following:

1. Earthquake occurrence exhibits scale-invariant properties (Kagan, 1994):

1. Earthquake size distribution is a power-law (the Gutenberg-Richter relation for magnitudes or the Pareto distribution for seismic moment) with an index of about 0.6. Preservation of energy principle requires that the distribution should be limited on the high side; thus I use the generalized gamma or modified Gutenberg-Richter distribution. Both index and the “corner” (maximum) moment have universal values for shallow earthquakes (depth interval 0-70 km), occurring in continents and continental boundaries (Kagan, 1999b; 2002a,b).

2. Earthquake occurrence has a power-law temporal decay of the rate of the aftershock and foreshock occurrence (Omori's law), with the index 0.5 for shallow earthquakes, and more than 0.8 for deeper events.

3. Spatial distribution of earthquakes is fractal: the correlation dimension of earthquake hypocenters is equal to 2.2 for shallow earthquakes. The value of the fractal dimension declines to 1.8-1.9 for intermediate events (depth interval 71-300 km) and to about 1.5-1.6 for deeper events.

4. The stochastic 3-D disorientation of earthquake focal mechanisms is approximated by the rotational Cauchy distribution.

2. Because the fractal dimensions of the earthquake process are significantly smaller than the dimension of the embedding space, this suggests possibilities for temporal and spatial forecasting of earthquakes. The fractal spatial dimension of shallow earthquake hypocenters is approximately 2.2 (less than 3). This suggests that seismic energy release is localized on a small subset of the Earth's crust. Moreover, the temporal clustering of shallow earthquakes has a fractal dimension of approximately 0.5 (less than 1.0), i.e., the seismic energy release is concentrated in the foreshock-mainshock-aftershock sequences occupying vanishingly small time intervals (Kagan and Knopoff, 1981; 1987). As an experiment, David Jackson and I (2000) proposed long-term and short-term forecasts based on smoothed seismicity for magnitude 5.8 and larger earthquakes. We have been forecasting future seismic activity in the western Pacific region where more than 50% of global earthquakes concentrate, daily since 1999 and we annually test the performance and predictive power of our algorithm (see also http://scec.ess.ucla.edu/~ykagan/predictions_index.html ).

3. If the conjecture on the universality of the earthquake size distribution is correct, it opens a new approach for the estimation of the seismic hazard. Due to the development of quantitative plate tectonics and spatial geodesy techniques, the tectonic deformation rate can be evaluated for all the plate boundaries and continental regions of extended plate deformation. Knowing the earthquake size distribution we can integrate it to obtain an estimate of seismic moment rate. By comparing this rate with the rate of tectonic deformation, one can evaluate the upper bound for earthquake size, and hence translate tectonic deformation rate into an estimate of seismic activity (Kagan, 1999b; 2002a,b; Bird et al., 2002; Kreemer et al., 2002). Such estimates can be compared with the estimates based on the seismic history (see the paragraph above) to improve our knowledge of earthquake hazard.

4. An important issue in the study of earthquake occurrence and the seismic hazard is verification of seismicity models. Until recently seismic event models and predictions were based exclusively on case histories and it was widely believed that the long-term earthquake occurrence at least for large earthquakes is quasi-periodic (seismic gap and characteristic earthquake models). The Parkfield earthquake prediction experiment and recent estimates of earthquake probabilities in the San Francisco area are based on these models. However, when David Jackson and I (1991; 1995; see also Rong et al., 2003) tested the seismic gap models against the earthquake record, it turned out that the performance of the gap hypothesis is worse than a similar earthquake forecast (null hypothesis) based on random choice (Poisson model). It is most likely that instead of being quasi-periodic, large earthquakes are still clustered in time and space (Kagan and Jackson, 1999). The importance of rigorous testing of earthquake prediction hypotheses was the main focus of our (Geller et al., 1997a,b; Kagan, 1997a,b) debate on feasibility of precursors discovery for individual earthquakes.

5. We proposed a model of random defect interaction in a critical stress environment (Kagan and Knopoff, 1981; Kagan, 1982; 1994) which, without additional assumptions, seems to explain most of the available empirical results. In the time domain, the Omori's law of foreshock/aftershock occurrence and, in general, the temporal clustering of earthquake events, has been shown to be a consequence of a Brownian motion-like behavior of random stresses due to defect dynamics. I conjecture that the fractal pattern of earthquake fault geometry is due to the fault self-organization in the conditions of high lithostatic and tectonic shear stresses. Theoretical calculations, simulations and measurements of the rotation of earthquake focal mechanisms (Zolotarev, 1986; Kagan, 1994) suggest that the stress in the earthquake focal zones follows the Cauchy distribution which is one of the stable probability distributions. Therefore, the evolution and self-organized aggregation of defects in the rock medium are responsible for fractal spatial patterns of earthquake faults and rotation of earthquake focal mechanisms.

6. Studying the distribution of earthquake focal mechanisms and modeling of earthquake fault branching opens new avenues for statistical analysis: very little is currently known about the statistics of 3-D rotations, quaternions, and the properties of tensor-valued stochastic processes (Kagan, 1991; 1994; 2000).

How can the results described above influence the development of earthquake science? At present, great technological advances have been made in obtaining earthquake information in increasingly broad temporal, spatial, and frequency domains: space geodetical methods allow us to see the surface strain field over all continents; registration of slow and ultra-slow earthquakes increases our knowledge of deformation in the plastic lithosphere; dense networks of deep borehole seismometers allow identification of new types of seismic signals, etc. However, collection and interpretation of new data should be guided by theoretical understanding, which, in my opinion, is still deficient. The present-day models of earthquake fracture are largely based on experiments with man-made rock surfaces, even though such surfaces have a scale (or scales) imposed by machines. As I discussed above, the empirical and theoretical evidence suggests a scale-invariant nature of the earthquake process. Moreover, these models use mathematical tools going back to the middle of the 19^th century.

As Mandelbrot (1983) points out, the fractal concepts are based on the mathematics of the late 19^th—early 20^th century. Completely new mathematical fields have been developed since that time, and to understand and describe 3-D large-scale deformation of brittle solids we would need to apply new differential geometry, group-theoretical methods, topology of defects in condensed matter, and other modern mathematical tools (see discussion in Kagan, 1999a). I do not expect that even the involvement of pure mathematicians will soon solve the fundamental problems of brittle fracture. In 1992 I argued that the theoretical problems of seismicity description are even more difficult than turbulence theory—the major unsolved problem in modern science, even though such first-rate mathematicians as Kolmogorov, Yaglom, and others investigated turbulence.

Fracture modeling can also help in the theoretical advances. For example, calculations of molecular dynamics demonstrate that the basic properties of solid fracture can effectively be derived from simple physical laws. Similarly, precise laboratory measurements of the fault propagation demonstrate multiple branching of fault surfaces. These simulations reproduce the fractal character of a fracture.

Since the early 19th century the Gaussian (normal) distribution was almost exclusively used for the statistical analysis of data. However, the Gaussian distribution is a special, limiting case of a broad class of stable probability distributions. These distributions which, except for the Gaussian law, have a power-law tail, recently became an object of intense mathematical development (Mandelbrot, 1983; Zolotarev, 1986; Uchaikin and Zolotarev, 1999; Nolan, 2002). The use of these distributions is numerous in physics, finance, and other disciplines. Application of these distributions to the analysis of seismicity and other geophysical phenomena would significantly increase our quantitative understanding of their fractal patterns.

Another important issue is the development of reliable, rigorous, reproducible methods for testing earthquake model hypotheses. Objects of study in geoscience are usually inaccessible to direct observation and experimental validation; the results of inversion are highly non-unique. Thus, in almost all cases the standard techniques of a hypothesis falsification or verification, used in other natural sciences, are not applicable. Therefore geoscientists are not trained to apply such techniques. Hence many geophysical results are not reproducible and often are artifacts due to various defects of the analysis. Rigorous analysis of possible errors including systematic biases is often absent or insufficient. Proposed earthquake occurrence hypotheses are usually not falsifiable (Kagan, 1999a).

It took several decades in medical research to develop techniques (double-blind, placebo-controlled experiments) for reproducible testing. Such methods are not applicable in geophysics. However, in some cases an earthquake occurrence allows for a direct falsification of proposed models: about 10 large (magnitude greater than 7) earthquakes occur annually. Thus, on a global scale or in large seismic regions there is a sufficient number of earthquakes to test most earthquake occurrence hypotheses over a span of a few years.

The problems and challenges facing earthquake science are briefly discussed in a series of my slides, presented at the 2002 workshop in the Minnesota Institute for Mathematics and its Application. See: