V. G. Gitis
Institute for Information Transmission Problems, Russian Academy of Sciences, Moscow, Russia
Zhang Zchaocheng, Wang Guixuan, and Qin Xinxi
Center for Analysis and Prediction, State Seismological Bureau, Beijing, China
Retrospective predictions of Tangshan (July 28, 1976, M = 7.8) and Datong (October 19, 1989, M = 6.1) earthquakes are considered as an example. Daily time series obtained from 10 stations since 1972 for Tangshan and from 10 stations since 1981 for Datong were processed to model the process of earthquake preparation and to detect the earthquake precursors. The suggested space-time approach is shown to be effective for research on earthquake prediction.
The obstacles to the solution of the latter problem have several main causes:
Nonetheless, we may presume on the basis of a vast experience accumulated in many countries, that precursors were actually observed in a number of cases. On the observation station, their form can be approximated by several types of signals; i.e., the appearance of distortion of the trend, the baylike anomaly, the unipolar or double-polar impulse. It is also assumed that the precursors can occur and disappear repeatedly in the course of preparation of an earthquake.
A great amount of research was dedicated to the determination of the methods of recognition of earthquake precursors from seismological, geophysical, hydrogeological, geochemical, etc., data. We shall mention here but a few of the summarising works [Ma Zonglin, et al., 1989; Continental Earhquakes, 1993; Sobolev, 1993].
The representativity and completeness of data on the processes of earthquake preparation are determined primarily by the density of the measurement network and by the time longevity of synchronous observations. Apparently, the prognostic observations on intracontinental earthquakes accumulated in China are the most impressive by their volume and the time periods they cover. In contrast to the observations collected in the Benioff zones, in China it is possible to analyze the measurements of the stations surrounding the epicenter of an earthquake.
In the present paper we suggest an approach to the analysis of geodynamic processes from multidisciplinary data of geophysical monitoring and discuss the results of the tests with the methods of short-term earthquake prediction; we also make preliminary estimates (within the scope of the available experimental data) of the significance of geophysical anomalies preceding the Tangshan (north-eastern China, July 28, 1976, M=7.8) and Datong (October 19, 1989, M=6.1) earthquakes.
Table 1 gives a total list of the series of diurnal values of geophysical, hydrogeological and hydrogeochemical parameters, which were used in the process of research.
The stations 1 to 10 were selected for the analysis of the Tangshan earthquake; the stations 7 to 16 were chosen for the analysis of the Datong earthquake.
In the course of the analysis, two major principles were strictly observed: I) we included data limited by one date of 24 hours before the earthquake to prevent using later data; 2) we applied observation series of equal length. In this way, the analysis of variations before the Tangshan earthquake was based on the observation period ranging from the beginning of 1972 till July 27, 1976 and from the beginning of 1981 till October 18, 1989 for the Datong earthquake.
When processing the time series, we used the simplified additive model of the signal. It was presumed that the measured signal X(t) is composed of the response to the geodynamic process S(t), periodic components of the seasonal rhythm U(t), and the high-frequency noise with the zero average V(t):
X(t) = S(t) + U(t) + V(t). | (1) |
Figure 1 (a)
shows an example of the time series of diurnal values of tilts
at station FS. Intensive seasonal variations of the annual period are apparent;
they are complicated by greater high-frequency peaks; the beginning of the latter
in July coincides with the periods of intensive rainfall, as indicated by meteorological
data. The spectral analysis of the initial data shows that all seasonal variations have
maximums of the annual period. Figure 1 (c)
gives an example of the changes in the
water level in the bore-hole at station TS. At this point the periodicity of the
seasonal variations, though less distinct, still shows numerous one-day peaks, most
of which were caused, apparently, by technogenic interference or random errors of
the operators. In order to exclude the peaks for all signals, we conducted the
median smoothing out with a 5-days interval.
It is assumed that after elimination of seasonal rhythms in the absence of earthquake preparation the dynamic fields are inhomogeneous in space, but quasistationary in time. The appearance of a precursor, which occupies a certain related subset of elements of the raster, violates the stationarity of the process.
In order to determine the non-stationarity, the current observation interval is divided into two consecutive sub-intervals T1 and T2, of which the former is several times longer than the latter. The duration of both sub-intervals are determined by the researcher depending on the statistical properties of the background signal within the interval without precursor and on the expected duration of anomalous changes. The problem is reduced to testing the hypothesis of statistical homogeneity of two random sample sets. The hypothesis about the coincidence of parameters of the sample set distributions is tested by a certain statistics, which depends on the selected statistical model. In the case when several dynamic fields of different physical nature are jointly analyzed (the multidisciplinary approach), their values at raster points are multidimensional vectors, and multidimensional statistical models are used to verify the hypothesis.
In our case, the vector dynamic fields, computed from the time series of different types of the measured parameters, are not representative owing to the restricted number of stations. We, therefore, resort to a compromise by using the GEOTIME system to obtain a unified scalar dynamic field for all temporal series recorded at all available stations regardless of the physical nature of the measured values.
For this purpose, the signals from all stations should be first of all standardized (normalized by variances).
After standardization, the data processing is carried out by three stages, i.e., I - cleaning the time series of seasonal rhythms; II - calculation and transformation of dynamic field of precursors; III - calculation and analysis of dynamic fields of precursors' significance.
The relatively short series of diurnal observations before the Tangshan and Datong earthquakes of about 5 and 9 years duration, correspondingly, did not allow us to analyze the long-term precursors, and we analyzed the possible appearance of only middle- and short-term precursors. The rather large and complicated by form variations of seasonal rhythms also provided little opportunity for identification of precursors with duration of more than a year. Consequently, at the present state of research, the analysis was aimed at finding anomalous changes of from one to several months duration.
The pattern of seasonal rhythm was estimated by the time series interval 2 years long. It was assumed that the seasonal rhythm remains unchanged during the next year (the year of prognosis). For example, the seasonal rhythm for the Tangshan earthquake was estimated by data of 1972 and 1973 and predicted for 1974. This procedure was repeated for every following year. In order to retain for analysis the whole observation series before the Tangshan earthquake since 1972, the image of the seasonal rhythms, based on the material of 1973-74 and 1974-75, was predicted backwards, i.e., for 1972 and 1973, respectively. The seasonal patterns computed for all these years were then subtracted from the initial time series.
Figure 1 (b, d) shows examples of signals obtained after subtraction of seasonal rhythms. It is obvious that these realizations contain not only anomalies associated with incomplete removal of the seasonal trend, but also other deviations of unknown nature and probably the earthquake precursors among them.
For different physical parameters the precursors can be both positive and negative anomalies, as illustrated, in particular, by the given examples. In 1975-76, at one of the stations, a lowering of the signal's level was recorded (Figure 1 (b)), whereas at another station (Figure 1 (d)) the level has risen, and both anomalies can be regarded as the possible precursors of the Tangshan earthquake.
![]() | (2) |
where Z(l, j, t) is interpolated value at point (l, j, t), Yn(t) is value of signal of n-station at moment t, rn is the distance from n-station to the point (l, j, t), with rn = 0, Z(l, j, t) = Yn(t).
![]() | (3) |
where
![]() |
![]() |
T1 = 73, T2 = 6, which at 5-days step of the time co-ordinate corresponds to 365 and 30 days.
It is obvious that with H0 hypothesis the mathematical expectation of statistics G(a, j, t) is close to 0. But where as the counts down of the Z(a, j, t) signal are strictly correlated, the estimation of mean square deviations of the difference of selected averaged values, if H0 hypothesis is applied, meets with great difficulties. Let us assume that with H0 hypothesis the variance s2(a, j, t) of the signal is a slowly varying function of t and so it is the same in both windows T1 and T2. Then the upper limit of the root mean square of random variable G(a, j, t) is equal to 2s(a, j, t). It is convenient to normalize the G statistics by the statistical estimate of root mean square:
![]() | (4) |
where
![]() |
The u(a, j, t) statistics evaluates the degree of deviation from the stationary. With H0 hypothesis the statistics has the expectation close to 0 and the root mean square less than 1.
The analysis of the maps reveal anomalies in July-August of some years owing to incomplete removal of the seasonal rhythm. The periods from September through June are free of this drawback, and their anomalies, apparently, had other causes including those of geodynamic character.
Unlike the two previous years, in September of 1975, a powerful anomaly with a
maximum occurred over almost the entire studied area of northeastern China. It gradually
attenuated and practically disappeared in April 1976. Only two anomalous stations remained,
TS and SQ, in the region of the future Tangshan earthquake. In May 1976, the anomaly in
the area of these stations grew intensively, as shown in
Figure 2 (1)
representing the
map of that period. At the same time, a less intensive anomaly appears in the Beijing
area. For easier orientation, Figure 2 (1)
also shows tectonic faults and epicenters of
the future Tangshan earthquake and of the two strongest aftershocks with magnitudes 6.9
and 7.0. There is a striking irregularity in the distribution over the territory of
recording stations that supplied the initial observation series.
The next map in Figure 2 (2)
represents the period from May 14 till June 13.
At that time the anomaly in the Beijing area reached its maximum, and that in Tangshan
was steadily growing. The map in Figure 2 (3)
(May 29-June 28, 1976) shows the maximum
of the anomaly in Tangshan region; after that time it gradually attenuates till the
moment of the earthquake (Figure 2 (4) ).
The dynamics of the development of anomalies
along the A-B profile of latitudinal strike (Figure 2 (4))
can be followed in Figure 3
.
The curves 1, 2, 3, 4 correspond to the periods of maps in
Figure 1 (1, 2, 3, 4).
For comparison, Figure 4
shows 6 maps of the previous three years covering the
period of the beginning (May 4-June 3) and the maximum (May 39-June 28) of the Tangshan
anomaly. It is apparent that during all these years, in the Tangshan area, no significant
anomalous changes were recorded in the field of the complex of studied parameters.
An analysis of these maps shows that the dynamics of the field of the studied
parameters largely differ from the dynamics before the Tangshan earthquake. This is
evident from the series of maps in Figure 5.
During the period from July 21 till
August 20, 1989 (Figure 5a),
in the region of Tangshan, an anomaly appears, which
then grows in intensity and spreads west and east. The anomaly reaches its maximum
in the period from August 20 till September 19 (Figure 5c),
and then reduces (Figure 5d).
In this case, the irregularity of distribution of measurement stations becomes apparent.
In the southwest, there is only one TY station, and there are no data for the region of
the Datong earthquake.
For comparison, a series of maps for the previous six years is demonstrated in
Figure 6
(1983-1988) for August - September; during that period in 1989, the anomaly
in the studied region reached its maximum. Apparently there no considerable anomalous
changes in intensity and area at that time.
The dynamics of the development of the anomaly in 1989 before the Datong earthquake
are shown in Figure 7.
The plots 1, 2, 3, 4 along the latitudinal profile A-B correspond
to the periods of corresponding maps in
Figure 5a, b, c, d.
The migration of the process
is clearly traced in the western direction; it is less distinct in the eastern direction
where the region is confined to the location of CL station. Unfortunately, the data of JX
station hear Haichen were not available to the anthars for the period of the Datong
earthquake.
The epicentral character of manifestation of precursors was discussed in [Chen Yong, et al., 1992; Zhang Zhaocheng, et al., 1992], where retrospective data are given on the dependence between the relative number of stations that recorded the precursor and the distance to the earthquakes epicenters. For the earthquakes with magnitude 7, these data are given in Table 2.
Let us assume that at distances of the order of 1000 km, the precursors cease to bee observed and the amplitude of a precursor diminishes approximately by the same law as the number of stations that recorded the earthquake precursors. In this case, Table 2 can be approximated by the attenuation function
F(R) = A(exp(-0.02R) - 0.27) | (5) |
where A is the amplitude of the signal, R is the distance to the source of the signal in km.
In this part of the paper we shall make an attempt to model a generalized precursor of the Tangshan earthquake by using the attenuation model (5) to localize the future source.
Let us indicate the time series obtained on observation stations n = 1, 2, ..., N by u(n, t). We shall assume that, for any fixed moment of time, the center of the precursor's signal coincides with one of the points of the raster k = 1, 2, ..., K and that the values of the time series at the observation stations are expressed by relationship
un = AF(Rnk) + xn | (6) |
where A is the amplitude of the signal, k is the number of the raster point, in which the center of the signal is located; xn are the independent Gauss values with zero mathematical expectancy and the mean square deviation s.
![]() | (7) |
With known k =
and if the signal propagates from point
,
the estimation of the amplitude A(
)
is determined by expression
![]() | (8) |
By sorting out all the possible values of k = 1, 2, ..., K, we obtain different
values of A() and residual sums
of squares
![]() | (9) |
In order to evaluate the center of the signal - the K point-, let us analyse the criterion of relation of square estimation of the measured signal to the square of discrepancy (signal/noise ratio)
![]() | (10) |
where
![]() | (11) |
We should note that the vectors of the signal F = (F(R1k), F(R2k), . . . , F(Rnk)) and of the discrepancy, obtained by the least squares method, are orthogonal and in sum compose the vector of observations u = (u1, u2, . . . , uN).
Therefore,
![]() | (12) |
and
![]() | (13) |
If we suppose that at the analysed moment of time the signal is absent, that is un = xn, then a(k) is the relation of the Gauss square value with zero mathematical expectancy and dispersion s2 to the sum of squares N-1 of Gauss values with the same parameters.
It is convenient, instead of a(k), to consider the monotonously associated with it value
![]() | (14) |
Value b(k) acquires the magnitude ranging from 0 to 1.
The value b(k) = 1 corresponds
to the situation when un = AkF(Rnk),
i.e., the signal/noise ratio is equal to .
At b(k) = 0, the
signal is absent. At b(k)
0.5,
the discrepancy in parameter A evaluation is greater
than the amplitude evaluation. Therefore, independently of the value of amplitude
Ak estimation, we can admit that, in the studied region, the signal corresponding to the
assumed model was not observed.
On the basis of these considerations, at b(k) 0.5,
the estimation of the position of the signal
center we shall choose for the condition of the maximum of function b(k),
k = 1, 2,..., K, i.e.,
![]() | (15) |
and (
)
value shall be assumed as the estimation of the amplitude of the signal.
The first stage of processing, as in part 3, included standartization of series,
suppression of the seasonal rhythms, and squaring. At the next stage, the time series
of deviations for stationarity were calculated by formulas (3) and (4) with substitution
in them of coordinates l, j
by the number of the stations is at the value of parameter
T1 = 365 days and T2 = 30 days.
The last stage was the transition to the space-time model
by estimation of parameters of the model
and
(
)
with 5 days interval.
Figure 7
shows plots of changes of criterion of the presence of signal
b(),
estimation of its amplitude
(
)
and the distance r(
, C)
between the center of signal
and the real earthquake epicenter. It is apparent that during almost the whole
interval up to August 1975 b(
) < 0.5.
This means that, according to our supposition,
the epicentral precursor signal was either absent, or did not comply to the model (5).
A year before the earthquake, since August-September 1975, the b(
)
values sharply increased and exceeded 0.5. In other words, at that time interval, the suggested
model conforms best with the set of real data. Concurrently, an essential increase
of the signal's amplitude was recorded. Finally, since November 1975, the center of
the modelled anomaly was shifted to the epicentral area of the Tangshan earthquake
with deviations from the future epicenter reaching about 50 km.
Two approaches were realized:
The following preliminary conclusions can be drawn from the analysis of the initial data and maps at our disposal showing the complex of the studied geophysical, hydrogeological and geochemical parameters.
From May till July 1976, an anomaly was recorded in the region of the future Tangshan earthquake; its amplitude, probably, exceeded the random deviations of the background noise. The coincidence of the maximum of the anomaly with the Tangshan earthquake and the absence of such anomaly in the analyzed previous years does not contradict the hypothesis about its geotectonic origin. If that supposition is true, then the anomaly is a medium-short -term precursor that appeared in the vicinity of the earthquake's focus. Let us recall that precursors with duration longer than a year or less than a few days are not analyzed in this paper on account of restricted observation series.
On the other hand, the Datong earthquake of October 19, 1989, occurred against the background of regional anomaly covering a large part of north-eastern China and localized near the epicenter. We believe, that the earthquake proper is only indirectly connected with the observed anomaly in the field of geophysical parameters. Both events could have been caused by regional changes in the stress state of the environment the nature of which is as yet unknown.
We wish to emphasize that our conclusion is open to further revision, because it is based on obviously incomplete experimental data and extremely irregular observation network.
Other conclusions of the present research are as follows.
It is necessary to eliminate all possible errors of external origin in the initial data.
It is desirable to add to the analysis the available observation series from a larger set of stations preferably covering the area with greater regularity.
It seems expedient to carry out a formalized comparison of the available series of prognostic observations with similar series of certain meteorological parameters and, primarily, the air and soil temperature, the amounts of precipitation, and atmospheric pressure.
It is important to include into the processing the data on seismicity for better insight into the geodynamic process.
The prognostic features could be more readily identified by creation of models of space-time precursors formation based on the physics of the earthquake focus and on the experience of experimental observations.
Meanwhile, the elaborated variant of the GEOTIME system allows us to study the geodynamic process using a complex of heterogeneous data of prognostic observations. The experience of the present research has indicated the direction in which further development of the system should be continued.
The research is supported by grants No. 97-05-65906 and No. 97-07-90326 from Russian Foundation for Basic Researches.
2. Continental Earhquakes. Seismol. Press, Beijing, 1993, p. 576.
3. Sobolev G. A., Fundamentals of Earthquake Prediction. Nauka, Moscow, 1993, p. 311 (in Russian).
4. Gitis V. G., Osher B. V., Pirogov S. A., Ponomarev A. V. , Sobolev G. A., Jurkov E. F., A System for Analysis of Geological Catastrophe Precursors. Journal of Earthquake Prediction Research, Vol. 3, no. 4, 1994, pp. 540-555.
5. Gitis V. G., Osher B. V., Pirogov S. A., Ponomarev A. V. , Sobolev G. A., Jurkov E. F., Dynamic Fields Analysis System. Cahiers du Centre Europeen de Geodynamique et de Seismologie, Vol. 9, 1995, pp. 129-140.
6. Chen Yong,Wang Wei, Zhu Yueqing and Ji Ying., Multidisciplinary approach used in expert system for earthquake prediction in China. Journal of Earthquake Prediction Research, Vol. 1, no. 1, 1992, pp. 107-113.
7. Zhang Zhaocheng, Zheng Dalin, Luo Yongsheng and Jia Qing., Studies on earthquake precursors and the multidisciplinary earthquake prediction in China mainland. Journal of Earthquake Prediction Research, Vol. 1, no. 2, 1992, pp. 191-205.
Load: | [HTML document with graphics for local use] |