- Email: [email protected]
Premonitory seismicity patterns near Vancouver Island, Canada
Tectonophysics, Elsevier
167 (1989) 299-312
Science Publishers
299
B.V., Amsterdam
Premonitory
in The Netherlands
seismicity
D. BROWN
patterns near Vancouver Canada
‘, Q. LI ‘, E. NYLAND
Island,
’ and D.H. WEICHERT
of Physics, lJniuersit_v of Alberta Edmonton,
’ Department ’ Geologml
- Printed
2
Alta. T6G 2JI (Canada)
Survey of Canada, Pacific Geoscience Centre, Box 6000, Sidney, B.C. V8L 4B2 (Canada) (Received
November
4, 1987; revised version
accepted
August
18. 1988)
Abstract Brown,
D., Li, Q., Nyland,
Canada.
E. and Weichert,
In: M.J. Berry (Editor),
The search
1989. Premonitory
Hazard
Assessment
seismicity
patterns
for large seismic events requires
and Juan
de Fuca
plates
M z 6.5. The
strongest
earthquakes
three
complete
generalized
for M 2 2.5 are preceded
sense depicted
by the algorithm
been preceded
by activation,
occupy
30% of the total
about
significant
similarity
by activation “M8”
suggested
have
of the earthquake
time-space
domain
the earthquake
of 5 years and 100 km. Such anticipation
could justify
attempts
to make short-term
using patterns
There is evidence
is a natural
consequence
predictions.
Analysis
for the coupling
of seismicity
of the model of a lithosphere
Introduction We seek to diagnose the approach of the strongest earthquakes in that part (Fig. 1) of the Pacific North American plate boundary where the Farallon plate fragments into the Explorer and the Juan de Fuca plates which then subduct under the
them. In the
catalog
of the world
future
strong
some precautionary of quiescence
where
earthquakes measures
restrict
has
been
range, In the 6.0 may have
(TIP) of a strong alarms”.
to be
appears
One event of magnitude
regions
as an assemblage
to define
flow in the lower magnitude probability
I!land,
167: 299-312.
the seismicity
are no “false
in the subduction
Vancouver
this threshold
there
to be tested, of anticipating
accuracy earthquake.
and
flow here and in other
which remains
since
by Soviet scientists.
considered
threshold
Canada)
occurred
near
Tectonophysics,
an a-priori Island,
but two others were not. The times of increased
between
tested. There is a possibility,
(the area of Vancouver
(M z 6.0) which
patterns
and Prediction.
area of the Explorer acceptably
for premonitory
D.H.,
Earthquake
This
algorithm
earthquake implies
a
MX was
in this area with an
and could be usef 11 for
the location
of the anticipated
zone with that above it. Such inter2 ction
of interacting
blocks.
the present state of seismicity analysis as a means of describing the phenomenon of earthquake flow and of suggesting hypothetical regularities. Empirical data are necessary to formulate such hypotheses,
and they are as insufficient
in ..he area of
in-
Vancouver Island as they probably ,ire everywhere. As a consequence we create hypotheses which might apply in diverse regions and use them
northwest of the of the southwest
to study the similarities in the situations before strong earthquakes. The approach se.lrches for
coast of Canada extending to the north end of Vancouver Island. For convenience in this paper we refer to this zone as the Vancouver Island area. Since there is no comprehensive theory which would define a priori an algorithm of earthquake prediction, a-posteriori data fitting is necessary in
those features of seismicity which are I-easonably independent of the tectonic environment. There is a quite surprising number that can be identified world wide. The approach used here considers tectonic information only in so far as seismicity data can be crudely associated with major tectonic
North
American
cludes United
a part States,
0040-1951/89/$03.50
plate.
The region
of the Pacific and that part
considered
0 1989 Elsevier Science
Publishers
B.V
iO(,
Once such a preliminary
features.
the approach
searches
for anomalous
to large earthquakes
without
strictions
related
intended
in this paper
merely comment
to tectonic
expressed
and CN, developed Institute
Moscow
based
specific
on
hypotheses.
re-
It is not
the patterns;
we
and reproducible
a relatively
wide
perceived
subregions patterns
free parameters applications
MS
that
For
world
lower thresholds equally
various
calculated.
This set of components
demands
a
set of criteria
set of premonitory
phenomena and makes effective use of the data available for the Vancouver Island area. Both algorithms recognize that the process of evolution to a large earthquake may be very different in different cases. Both algorithms are open in that the inco~oration of new data constitutes part of the aigorithm. Of the two we have chosen M8 as
a vector
and
describing
are a conse-
in time
window
of the
that
Marly of the
are fixed for all
P with
characteristics
of the catalog is divided
among
four groups.
For every
window
between
a threshold
in-
the components
M8 deals with seven characteristics
for each component
desig-
spaced
groups.
values
than
8 for which it was originally
every
OII
wide. in spite of the fact that
stants
Earth.
no hearing
in these regions. in the algorithm
M8 is now is use for much ned.
have
at the 0. Yu.
et al.. in prep.) reasoning
geographic
the magnitude
in the algorithms
for the Physics
of pragmatic
relatively
prior
artificial
by investigators
(Gabrielov
quence
patterns
making
to explain
is made
that they exist.
The hypotheses Schmidt
selection
is
into
divided
can be chosen
of the vector so that pS of at1
the start of a catalog
and the first
strong earthquake in the window, or the end of the catalog if there are no strong earthquakes in the window, are above this limit. All such values
the simplest.
are designated as extreme. We now examine the vector calculated at time t from t - 3 yrs to t and identify the components that were extreme at least
The M8 algorithm searches for process activation in a seismicity catalog and considers objects
once in this period. A group is extreme if at least one of its characteristics is extreme. We calculate
defined
by time moments
t and windows
centered
the number
of such extreme
characteristics,
on a small number of regularly spaced points in the territory to be analysed. These windows can be
the number of extreme TIP (times of increased
of various shapes and are a type of automatic division of the seismicity into possibly overlapping
(f. I + T) if for two consecutive of such extreme characteristics
Fig. 1. The study area.
h, and
groups, g, and declare a probability) for a period times the number is 2 6 and the
PREMONITORY
number
SEISMICITY.
VANCOUVER
of extreme
these parameters
groups
a pattern
applied
geodynamic
observations
time
record
can
series,
in the Vancouver
that indicate
phenomena, may
whose collective of increased
Although
(Press
et
1979;
Dobrovolsky
al.,
1980a, 1988: Gabrielov
MS is desig-
ers who
the approach
consult
wish
those references,
here. The problem
100 km and
a time accuracy
of a few
identify
for this algorithm
from a study of very large earthquakes
world wide
of a large
toring
We introduce functions
Keilis-Borok this
et
may
the essentials
earthquake
of time
in a region.
s and
in the catalog
consider by calcu-
lating
strong earthquakes, sufficiently complete
time t and all events in the interval [t-s. t]. First we separate the aftershocks from
M, and the start time of a catalog, tO. The algo~thm is
normalized to a fixed rate of production of mainshocks and can be applied in a very wide variety of situations if we accept the hypothesis of the self-similarity, world wide, of the earthquakegenerating process. The results of this study enhance our understanding of the underlying similarity of seismicity world
wide.
previously nate the
If such self-similarity
in fact exists,
unsuspected relationships may domiphysics of earthquake generation. In
addition to the study of self-similarity the analysis clarifies the way seismic volumes such as subducting oceanic slabs and overriding lithosphere may interact. Finally, the results aid in the prediction of destructive seismic events. A suggestion of a region in which there is an increased probability
of the catalog
that depend
than considered in the California and Knopoff. As a consequence
the
maximum
distance
at which
considered an aftershock the mainshock-aftershock
area by Gardner we increased the an event
could
be
by 50% (Table 1). Once separation has been
made by the application of magnitude-~tependent space-time windows (Table 1) we can consider the
TABLE
1
if a mainshock occurring
has magnitude
within
tude less than
T days
after
M or less, any event within the event and having
M is an aftershock
2.5
29
6
result.
3.0
33
11.5
3.5
39
22
Definitions
4.0
45
42
The basic assumption that underlies the pattern recognition approach is the assumption that the earthquake flow includes a large chaotic component. As a result, the specific symptoms of an approaching strong earthquake may be different each time. Accordingly the algorithm uses a set of
on
mainshocks using a modified form of the generally accepted thresholds (Gardner and Knopoff, 1974). Inspection of the events in the data that we used clearly indicated aftershocks at larger distances
of a large seismic event could be used to initiate more detailed studies using either other methods or more time-dependent data to refine the MX
of the functions of M8
to
by moni-
reduce the number of free parameters for the study to two. They are the threshold for defining
characteristics
al.,
and read-
recognition
flow of earthquakes a period
1975;
approach
we discuss
of time defined
is fully
Briggs,
is to use pattern
the TIP of a strong the total
and
et al., in prep.)
to implement
accuracy
derived
(TIP)
inbut
of the
of the algorithm
elsewhere
with a positional
The parameters
can be diagnostic
probability
the basis
described of
behavior
informative
earthquake.
but
Island
each of which taken
be insufficiently
earthquake
of a strong of about
patterns
be
of quantitative
since 1952. The algorithm
ned to identify
years.
longest
premonitory dividually time
approach
that can be used with M8 is the record
seismicity
of
of most of
earthquake.
recognition
to any multicomponent
the chronologically area
activation
In any case no TIP is raised for
3 years after a strong Such
is 2 4. The choice
guarantees
the characte~stics.
301
ISLAND
M
R
T
4.5
52
5.0
60
150 290
X3
5.5
70
6.0
81
510
6.5
91
790
7.0
105
915
7.5
121
960
8.0
141
9X5
R
a magni-
flow of mainshocks { m,
) where
34,.
B,(e))
mainshock.
m,
t is the origin and
that occurred
17,.
time such that I~. , > I,.
parameter
k
those
is the depth.
is the number
of
mainshocks
exercise
that occurred
normalize
in an interval
of
.s before a
greater
than a
set designated
mainshock.
of mainshock
of simply
summing
the events
which we designate
Z, and is as close to 20/yr
aftershocks
measure
as N,.
as possible
I.,.
which we consider that
occur
This function.
depends
on the choice
quakes
of
magnitude
as :V,. I_, . and
and
for another
Zz. The seventh
in M8 is the count
within referred
two
days
of
of
a
to as B below.
of -. M only in that earthless than M ._ are not
mainshocks. We can now consider the earthquake stream to be characterized at a time I by a vector of func-
If ,0 = b/3, then S is a measure of the rate at which the radius of circular fractures is generated, and if /I = 2h/3, then S is the rate at which fracture area is generated. Here. h is the usual coefficient in the energy magnitude relationship log(E) = A + bhf.(The parameter a is inserted to normalize the functions). We carry out this process for all events whose magnitude IV, is between M and an upper limit R (R= M, - 0.1) and the result as S(r; M, a,
s, (Y. p). This pair
of functions can be combined to produce a measure of the way in which seismicity concentrates in space. _M, M, s, (Y. p)
tions
deduced
from the preceding
R, .T)y
defines the ratio of the average radius of a fracture (S/N) to the average event separation for all events between M and R, if they are uniformly distributed in theregion (N”‘). The rate at which mainshocks are produced is not uniform. If we assume that on average changes in a linear fashion. the function -M, s) = N(r: M. t-f,,)
[ - 10 I -
2,) -
s
it
time period
.r.
The problem is to deduce. by pattern recognition, whether such a time is part of a TIP. The duration T of such a TIP. the time within which a strong earthquake is expected. is taken to be the 5 years used in all recent applications of M8 world wide. In order to declare a TIP at time f we require that over a period +, prior to and including I. three of the groups ( N,. N, }, ( L,. L, }. { Z,, Z, ). and { B ) contain functions that have extreme values
and that at least four of the functions N,. Z,. Zz. and B have extreme values. N2, L,. I.?. In this sense extreme values are those in the upper 10% of the range achieved by the function from the start of the catalog
S(K M, R, .s, a, P)
- N(r; M, I-ff,,-s)
by adjust-
set of functions
the number
activ-
I&
on the
to count
ity calculates a weighted form of this sum where the weights are proportional to the size of the
[/v@; ly, s) - A+;
Ilnear trend.
explicitly
the result and this is achieved
A more sophisticated
=
depend
o/’
M. If we expect to compare results with -. in other areas we must find some wav to
function
Z(r:
from a long-term
three functions
lower limit M. We refer to this count as N( I; _M. s).
define
the gencratmn
ing -M so that the average yearly rate of occmrence of mainshocks is in one case IO;‘vr for one
time I and which had a magnitude
earthquake. Instead we calculate:
deviates
to which
in the first e days after
the mainshock. It is a simple
the degree
mainshocks These
H(e)
number
measures
of the
+ is the longitude.
M is the magnitude
of \:ectors
( I,. A,. 9,.
has six components
where i is the sequence
X is the latitude, aftershocks
to bc a sequence
to the time of occurrence
of the first large (M > MO. the prediction threshold) event. Once this condition occurs the interval (I. t + T] is declared to be a TIP. There is no straightforward means of deriving this algorithm from more fundamental statements. It follows from the plausible idea that patterns should exist in the earthquake stream and that these patterns should, as much as possible, be physically realistic and flexible. In the final analysis its justification must rest on the pragmatic statement that it works in many areas of the world. An area such as the southwest coast of Canada can be considered as a collection of perhaps over-
PREMONITORY
lapping define
SEISMICITY,
regions
defined
latitude
is then
such proposed
region
a side of
of latitude.
Each
to check
by the catalog
regions
occur
per
that do not
per year are discarded
and
the M8 algorithm
is applied
determine
a TIP exists in any of them.
whether
to the remainder
to
once
of regionalization
a small
a large
is possiearthquake
If such an area is divided it is possible
or quiescent
event
a quiescent
to occur within similar
can
activation
quiescence
that may precede
of the latter. A
be seen on a 2-dimensional
of event number with
area
area the large event is likely
or on the boundary
process
to classify
in the 6 years prior
(Fig. 2). If an activated
against
time (Fig. 11). A pronounced activation
a large (M
in which
of the event
to the seismic contains
coupled
that contains
area
them as activated
more than one region for each region
an
of location
histogram
and a TIP does not always
to produce
into small cells of size of the order of the accuracy
The regions can overlap, and usually do. As a consequence the same large event can appear in appear
refinement
might occur is identified.
seismic
ten mainshocks
A further
of ble
to have
ten mainshocks
year. Those proposed contain
it may be possible
defined
is then examined
of at least
the regions
area of high probability.
The square region centered
over the time span covered
an average
We
the area
by equal intervals
R = e Mclp5-6 in degrees
whether
seismicity.
by first diving
cells defined
and longitude.
on this event
303
ISLAND
by their
such regionalization
into rectangular
length
VANCOUVER
of smaller scenario
magnitude
and
lack of large events events
is the
(Lamoreaux,
1982)
large shocks.
> M) earthquake. In the restrospective analyses of a catalog such results suggest geographic areas
The seismicity of the Vancouver Island area
that are sensitive to the onset of large earthquakes and may suggest areas that should be monitored with more care. By considering the intersection of
The seismicity of the west coast of Canada has been described by many authors (e.g., Milne et al.,
Elf. ISI
cl
Unknown seismic state Activated seismic state Quiescent seismic state AQ boundary
Fig. 2. A schematic
diagram
to illustrate
the idea of activation
quiescence.
52' N 51' N 50°N 49'N 48' N
46'N 45' N 44'N 32' W30° W 128°W 126' W 124' W 122'W
2o" w 118 W
Fig. 3. The seismicity of the study area. 2 00 !____!__--!----! ‘841
t
1
t86i
1
I
1881
1
1
(931
1
300
I
1
((341
1
1
1951
% 1
1955
1
1959
1 .t
l
1963
I
1
1967
1
1
197t
1
1
1975
1
1
1979
f
t
i9e3
1
1
._ _!_-__a
--_-,_.
-_
1 i
; i?
1
192q
6.00
5.00
4.00 __!____I____!
4 42 45
to
4
1
12 17 22 37 22
10 5 9 9
;
2 2 id 2
i 15
5 6 4 9 3 6 24 14 15 6 13 17 9 11 152 12 3
6 1
37 49 30 14 ; 61 28 24 7 19 33 3 23 17 15 18 R 19 13 3 26 26 20 14 4 39 52 24 16 6 22 13 5 15 4 9 9 3 4 7 9 8 5 25 11 42 15 40 4 18 7 78 43 71 42 23 55 67 7 25 129 72 10 24 5 77 48 2 9 3 1 2 23 ,____!____!____!____!_-__!_---!----!----!----!----!----!----
84 83
Fig. 4. A count table indicating the number of mainshocks
4 ; 5 4 4 i 1 2 1 3 5 2 2 1 4 1 7 1 13 1
3
4 4 8
2 2 2 1 2
;
1 2 I
i
; 1 1 1 4 ;
2
i
2
1
i
i
in a particular magnitude interval during a particular time interval for the entire catalog.
PREMONITORY
SEISIMICITY.
VANCOUVER
4.2 4.5 3.3
2.1
4.3 3.9 3.0
51°N. 50’N. 49-N .
3.8 6.3 5.1 4.6
48’N.
5.0 5.2 4.3 3.3 3.0
47*N -
4.6 4.4
46’N I 49N
2.6 5.1 5.0 3.4 4.8 R
3.6
4.8 3.9 4.6 4.3 4.8 4.3 4.3
4.4
3.7
5kT.6
4.1 4.5
3.7 I I
132; W 130&W 128”W
126&W 124”W
Fig. 5. The largest earthquakes
305
ELAND
122&W 120&W It%” W
in 1 o x 1 o rectangles.
1978; Basham et al., 1982; Basham et al., 1985). Figure 3 shows the geographic distribution of epicenters for events of magnitude 4 or greater in the area which we studied. In this area, from 118” to 130 o W and from 45 o to 51” N, we merged the seismicity catalog of the Geological Survey of Canada with the catalog from the U.S. National Earthquake Information Service. Particularly for the earlier years the Canadian catalog contained many significant earthquakes as far south as 45 o N. The duplicate events were removed manually. Usually, NEIS data were accepted in the United States and GSC data were accepted in Canada. If only one agency reported a magnitude the report of that agency was accepted. Events for which no magnitudes were reported were assigned a zero magnitudes while if a preferred mag~tude was given it was used. For these data this is usually Mr. We also removed those events since May 1980 that were clearly related to the Mt. St. Helen’s Volcanic Eruption. The data are summarized in Fig. 4 as a time-magnitude histogram and in Fig. 5 as a map showing the largest magnitude earthquakes in regions measuring 1” X 1 O. Patterns in space and time can be distorted by changes in the observation networks and associated changes in the degree of completeness of the seismicity catalogs at lower magnitudes. Evolution of the station configuration in the Vancouver Island area, and the associated increase in complete-
ness of the catalog at small magnitudes, has been discussed by Basham et al. (1982). The pattern recognition approach used in this work does not require such a stringent completeness condition as discussed there. Thus Fig. 4 indicates that since 1952 a sufficient number but certainly not all of the magnitude 2.5 earthquakes are present in the catalog for some parts of the test area. The effect of missing small magnitude events seems to predominate in the function that counts af-:ershocks and probably reduced its effectiveness as a diagnostic characteristic in Vancouver Island. SimiIar considerations apply to the time variable location accuracy of the catalog. The spatial resolution of the analysis in this paper is about 0.5 ‘, reflection epicentral accuracy. Focal depths are the least certain, and in some subareas nonexistent. However, the only application of focal depths is as an exercise to determine whether patterns in shallow seismicity precede an intermediate earthquake. Data analysis In the analysis of the subcatalog discussed above we primarily used the standard parameter set developed for M8 by an analysis of a large number of regions world wide. We adjusted the threshold for prediction A+, and the start time of the catalog to to satisfy completeness considerations (Fig. 4) as required to normalize the Vancouver Island earthquake flow to others used in MS and applied local judgment to adjust the windows used to eliminate aftershocks. Table 2 shows the nine strongest earthqu~es of the area and indicates that the threshold for prediction should be M, = 6.5. This choice is by no means unique but is a stable choice; small variations do not seriously affect our conclusions. In order to select earthquakes for prediction we chose, before beginning the prediction process, all mainshocks of magnitude greater than 6 that occurred since 1957. Since a mainshock of magnitude 6.4 that occurred in 1971 near 50 ON and 128O W was followed within 1 year by one of magnitude 6.0 near 49” N and 129O W we noted these events as possibly being predictive. The 6.0 magnitude is substantially lower than our threshold
306
I .4BLF. 2 ‘I he mne strongest earthquakes
m the region h
Year
Month
Da>
I872
12
15
5
1917
12
22
1918
12
1926
~-
min
-_
!.
Lat.
Long.
-._-__-_ I>epth Magmrude _____.__.-...__ 44
.7?
4X.6
121.4
15
4x
50.0
128.0
6.5
6
X
41
S.8
49.6
125.9
7.0
11
1
1
39
1x.0
4X.X
128.5
6.0
1939
2
x
6
39
25.X
49.1
12x.0
6.5
1946
6
23
17
13
25.X
49.8
125.3
7.3
1946
7
18
6
6
5x.5
49.5
129.7
6.5
1965
4
29
15
28
44.0
47.4
122.3
59
6.5
1976
12
20
20
33
12.0 -
49.0
128.7
18
6.7
for prediction so we place no great emphasis on the results for this event. We excluded two large aftershocks of magnitude 6.4 and 6.1 which occurred in 1972 and 1978, and two large events of magnitude 6.3 that occurred in 1954 and 1956 are not likely to be predictable since the preceding seismicity catalog is inadequate. Two events of magnitude 6.0 that occurred in 1958 and 1961 were ignored on the grounds that their magnitudes were both small and their time of occurrence was uncomfortably close to 1957, the beginning of a useful prediction dataset. We therefore set as our goal the prediction of one event near Seattle and
three west of the northern end of Vancouver Island. Figure 6 shows the TIPS that resulted from the analysis of all the seismicity data available to us. The heavy lines under the time scale indicate TIPS, the light lines indicate the center of the region in which the TIP occurred. The three events west of the northern end of Vancouver Island, of magnitude 6.0, 6.4 and 6.7, are associated with three regions in which TIPS occurred. Although the prediction for the magnitudes 6.0 and 6.4 should strictly be considered false alarms, we accept them as being within the uncertainty of
Fig. 6. The centers of the regions analyzed with M8. A heavy line on the pointer to a location indicates the occurrence of a TIP in the time interval read from the scale in the
upper right of
the diagram. In this analysis we used a magnitude
time in the catalog of 1952.
threshold of 6.5 and a start
PREMONITORY
SEISMICITY.
magnitude preceded
VANCOUVER
determinations.
The magnitude
by a TIP in two regions
the event, in another
region.
data
due
from
a region
for both regions The magnitude at a depth ceded
two years before
east of the event.
6.4 event
that contain
resulted
The
was seen
not seen in the region
east. The magnitude
by
It was not seen in the
6.0 event just to the northwest
in three TIPS and
6.7 was
and preceded
a 5 year long TIP, which expired
magnitude
307
ISLAND
to the
in a TIP
it.
mation The
regions
6.5 event near Seattle shallow.
their overlap
activated
at the bottom and
considered.
I
I
I
I
I
I
.6:0’,‘,‘,‘,‘,‘,‘,‘,‘,‘, \\\\\\\\\\ ,,,,,,,,,,,,,
by of
The space-
time volume is area (km2) multiplied
,,,,,,,,,, \\\\\\\\\\
common
by time (yrs).
I
are
fairly
speculation
the
intersections
common This
to all TIPS is shown
but
to explore tectonic
both the 6.5
does not prove
of the separate
the predictions,
In order
in that part of the
domains
it is certainly
more closely
activity
that
the
~~11 in fact sugges-
the potential
as a tool to predict
deeper events we attempted so separate the catalog into essentially crustal earthquakes and events related to subduction. Depth determirations in this area are not adequate for careful separation SO we
chose to discriminate by considering only events of depth known to be less than 25 km. This
I
Fig. 7. The regions covered by overlapping
large
whether
to each TIP. The part
of Fig. 7 and contained
the 6.7 event.
of shallow
/\1
regions
occurred
area. The quality of prediction is estimated space-time volume of all TIPS in percentage
here
ated with the TIPS all occurred
TIP did not appear argues strongly for coupling between the shallow and deep seismicity in this
volume
considered
of the sum of their area defines
that pro-
Nevertheless
esti-
region of high risk. We note that the events associ-
refine tive.
space-time
in Fig. 6 gives a quality
(Fig. 7) and it is an interesting
the existence of four successful TIPS, one almost successful TIP, and only one region in which a
total
shown
of 30%.
of the time domains
of 59 km but the activation
the event was largely
The analysis
TIPS and one time intervals
common
to all TIPS for each event.
308
1960
Fig. 8. The results
of an analysis
using
only
events
whose
depths
explanations
is clearly a set of Iithospheric events but does not include all of them. The results of an analysis of this subset are shown in Fig. 8; there are insufficient events in the remainder of the set to allow application of M8. As expected, since most TABLE
have
1970
been determined
1990
1960
IO be 25 km or less. See Fig. 6 for
of layout.
seismicity near Seattle is shallow, the TIPS around the Seattle event are unchanged. One region to the south of the 6.7 event is lost due to insuf~~ient events in its s&catalog but the others remain viable TIPS.
3
A summary
of the TIPS (times of increased
probability)
for the Vancouver
TIPS
Center of window rp
x
49*N
129ew
Island area
Strong earthquakes
Time
Date
48”N
128” W
128* W
~970,6.20-1977.6.2~
1970,6.20-1978,6.20
i 1971.3.13
6.7 6.4 I
1971.125
6.0
1976,12.20
6.7
i 19?1,3.13
6.4
1971 J2.5
6.0
( 1970,6.20-1975q6.20
72
32
1971,125
6.0
6.0
1983,6.20-1989.6.20
?
?
72
19654.29
6.5
72
48’N
124* W
1959,4.29-1%7,4.29
49ON
123O W
1957,10.29-196410.29
48*N
122OW
1960,4.29-1965‘4.29
49ON
121”W
1957.10.29-1964.10.29
47“N
121°W
( 1962,4.29-1967.4.29 1982,4.29-198924.29
* Time period
M
of
expectations (months)
1976.12.20
SOoN
Tie
from start time of a TIP to a strong
event of M > Me.
84 1%5,4.29
6.5
60 84
19654.29
6.5
36
?
?
84
*
PREMONITORY
SEISMICITY,
VANCOUVER
309
ISLAND
Center of window: 49N 129W .*****A * * * * ~*,***~,~*,“,*,* *I************* I***h*F******* ***********~* .************** .*****~******~************* *****A********* **********~***************** .***********I** I*************** .************** .**********~@t.&************ **************I .*********I**** A************* r*# * A * * * * * *r”*a***** ******** 1
Nl N2 Ll L2 Zl 22 B 1965
1970
1975
1970
1975
Center of window: 50N 128W
1965 Center of window: 48N 124W Nl N2 Ll L2
.* ** * ** **** * **-??**************** .********* * * * * * * * * I) *%.********** .**************** *******n******** .**************** a ‘** * * * * * * * I) * * * .*******n******** ********J******* .***************I **************** .**************** * * * * * * * * *2” A * * * * * ;c**1****** * * * *y&#.&* * * * *************** * * * * A* * **************** ************I*******~* * * **’ ‘* .q#* n * * * * .******** *** 1
Zl 22 B 1955
*
1960
1965
Center of window: 48N 122W Nl N2 Ll L2 Zl 22 B
‘*********3h**L**** .********* *F * .* *I) ** ** ** ******* * * * ***** .*********’ ********* ?* * * * * * * * * ********* ********I ********** .********* ********** ‘******y&+******* .**I** * * * * I * *?,**** .*I******* ************#r ‘L. * * * * * * * *** x
1955
1960
1965
Center of window: 47N 121W Nl N2 Ll L2 Zl 22 B 1955
u.********, * ****** * * * I) I ******. * ********I * * * * * *I * . ********* .*********. ****~*****~*******. ********* ***~*********** * * . I)********. **A****** *********. ******xArte******. ********* *I)*******. ********* ****~*~* ***** * * * * . 1960
1965
Fig. 9. The part played by each function in establishing the TIP. The shaded area are the time intervals in which the activation was sufficient to establish a TIP (referred to here as a “ vote”).
The analysis that exploits the entire dataset is associated with two TIPS for which large events have not yet occurred.
In one case the seismicity
the anticipated earthquake. We and others found it helpful to examine activatu~n on a more detailed
scale and
of
we show some of the regions
subdivided
than 6.5 began in the middle
of
cells. We have examined
1982 and will expire in 1991. If past experience
is
significant
activity
years prior
to the large event.
If no such activity
occurred
the cell was identified
as quiescent.
collection
of quiescent
data
indicate
magnitude a guide
that
greater this TIP
a TIP
is associated
end of Vancouver
Island.
ated with the northwestern States.
It began
expires
in 1991.
Some
traits
of an
earthquake
quiescence
h;t\e and
with the northern
The other TIP is associcorner
in the middle
of the United
of 1983 and
also
of the catalog
contribute
to the the
time intervals that have created TIPS for five windows. The asterisks indicate the occurrence of extreme values for the functions on the left at times
indicated
at the bottom
of the plots.
The
regions are the time intervals in which the activation was sufficient to establish a
TIP. Table
may
strong
TIPS and others do not. In Fig. 9 we illustrate
shaded catalog
cells
3 is a another
means
of summarizing
these results. As indicated in Fig. 6 the M8 algorithm is not particularly effective in specifying the location of
Cantor ot wlndow UN 12oW
Canter of wlndow: SON 128W
the
event
events
on
in it durrng
cells bounded the
will occur. plot
into 256 iI
each cell to determine
occurred
indicate
I-‘& It!
in
the
area
the o A
by activated
in which
the
next
In each of the four cases boundary
between
quiescent and activated areas, but the regions the northwest show a better defined area
the in of
quiescence than those in the southeast. The same process can be used to refine the time of occurrence of the anticipated earthquake. In time, the quiescence when
an earthquake
histogram time for
will help us to refine the time will occur.
Figure
11 is a
of event number against magnitude and the magnitude 6.7 and 6.5 events. A
quiescent period occurs before the magnitude 6.7 event and a courageous analyst would have predicted the event when the activity began to increase. The histogram of magnitude not show such an effect as clearly.
6.5 event does
Illllflllllllll!l
Conclusions
It is scarcely surprising that the seismicity of Vancouver Island bears a strong similarity to that found elsewhere. It would have been surprising if
Centor
of wlndow:
4ON 122W
Cmtw ot wtmlow 48N lZ4W
llllllllI_ccI
q ahsceflt
Vancouver Island prior to at least four large seismic events. This suggests the possibility of
cdl
El Ouioscont boundwy .
this area were radically different from the twenty areas previously tested by Keilis-Borok and his colleagues (Gabrielov et al.. in prep.). Apparently identifiable patterns exist in the seismicity of
Loutton ot strong wrthqwk*
Fig. 10. Plots of activation
and quiescence
events in the catalog.
for the four large
defining time and space intervals that indicate an increased probability of large earthquakes. Three successful retrospective predictions alone are not convincing, but the agreement of this method with many cases world wide suggests that the current TIPS in the northwest and the southeast parts of the study area should be taken seriously. One earthquake of magnitude 6.0 is included in the TIP that precedes a larger event, but two others are not. This indicates that seismicity of
PREMONITORY
SEISMICITY.
VANCOUVER
311
ISLAND
Center of window: 5ON 12’7W 4.00
2.00 3.00 !____!____!____!___ 69
1
70
1
I
0
f
3
0
1
71110 72
1
73
1 IO
74
1
75
1
i I
1 0
1 1 1
10 10
76
1
10
77
1
10
4 5 4 8 8
78110
---
!----!
3 9 166
6.00
5.00
.- !----!----!____!--__!____!
2 7 1 3 1 1
1 3 1 3 2 2
2
; 1 2 1 2 1 6 6 1
2 ; 1 2 2 2 3 2 1
i
; 1
; 2
ii
2 3
2
1
3 1
;
i 2 1
2 2 I
i : i 1 3 1
:
:
:
1 1
;
: I
: 1
i 1 1 2 1 1
:
7.00 _!____!___-
1
i
:
2 1
: 1
:
:
8 13 7 12 6 8 4 2 5 5 6 9 6 15 18 12 14 16 0
-!__--!____!--__!___-I___-!___-!
4 1
3 3
1 2
1 6
2
2
CVENTS
Center of window: 47N 123W 2.00 3.00 !____!____!____!-_ 58
1
1 0
59
1
10
60
1
I
61
1
10
62
1
1 0
63
1
1 0
0
13 10 9 8 8 9 6 6 5 4 2
64110
7 3 4 6 I 2 7 3 1 3 2 1 1
6.00 -!----!
-_-
-_-
5 3 2 6 2 1 1 2 1 1
:
1
0
26 16 15 24 II 13 14 13 7 5 6 5 2 4 2 3 0
1 1
1
1
65
1
1 0
i 1
2
:
8 2
4 4
2 6
66110
166
3
4
3
2
1
0
0
EVENTS
Fig. 11. Magnitude-time
histograms
of events in regions
centered
on two large earthquakes
(see text).
magnitude 6.0 is evidence of activation for events of larger magnitude. but that the activation phenomena required to predict events of magnitude 6.0 by MS do not exist in this catalog. The two events of magnitude 6.0 that are outside the TIP were not followed by large events but the one that was within the TIP was present there because the TIP was initially generated for an event of M = 6.4. Since this M = 6.0 event was foflowed by an event of M = 6.7 and the other two were not followed by large events, these smaller mainshocks were correctly interpreted by the algorithm. The failure to predict them is not a failure of the algorithm; in this area it appears to work only for magnitudes in excess of 6.4. In addition we observe quite clear patterns in shallow seismicity prior to the deep large event at 47O30’N, 122” W depth 59 km). This suggests that the subducting lithosphere may be coupled to shallower structures in some as yet unexplained way. This in turn implies that other near-surface geophysical measurements, such as strain and geochemical studies, may have potential as detectors of precursory patterns to strong deep earthquakes. At this time we make no statement on the potential for a megathrust earthquake in Cascadia. The interval of adequate seisrnicity catalogs may not be sufficiently large to identify this event. The only evidence that may have a bearing is the strangely low number of mainshocks in the area and the anomalous difficulty we encountered in assessing aftershocks. AcknowIe%lgements The analysis reported here was directed by V.I. Keilis-E3orok while he visited the University of Alberta in July and August of 1987. It benefitted greatly from discussions with Garry Rogers of the Pacific Geoscience Centre, Ge&gica.l Survey of Canada. We also acknowledge our debt to a legion of workers who developed the coding for the M8 algorithm.
Keilis-Borok’s visit was supported t’rom the Fund for the Future at the University of Alberta. The Natural Sciences and Engineering Research Council of Canada. the Department of Physics. and the Departments of Computer Services and Computing Sciences at the University of Alberta also aided this research in various ways.
References Basham, P.W., Weichert, D.H., Anghn, F.M. and Berry, M.J., 1982. New probabilistic strong seismic ground motion maps of Canada: a CompiIation of earthquake source zones, metbods and results. Minist. Energy Mines Reaourc., Earth Phys. Branch Gpen File 82-33. Basham, P.W., Weichert, D.H., &$n, FM. and Berry, M.J., 1985. New probabiBstic strong seismic ground motion maps of Canada. BuU. Seismol. Sot. Am., 75: 563-595. Dobrovobky, I.R., Zubkov, S.I. and Miachkin, V.I., 1979. Estimation of the size of earthquakes preparation zones. Pure Appl. Geophys.. 117: 1025-1044. Gabrielov, A.M., Dimitreava, O.E., Keihs-Borok, V.I., Kosobokov, V.G., Kuznetsov, I.V., Levshina, T.A., Mitroev, K.M., Molcban, G.M., NegmatuBaev, SKh., Pisarenko, V.F., Prozeroff, A.G., Rinehart, W., Rotvain, I.M., Shebalin, P.N., Shnirman, M.G. and S&eider, SYu., in prep. Cent. Reg. SismoI. Am. Sud., Lima. Gardner, J.K. and Knopoff, L., 1974. Is the sequence of earthquakes in Southern California, with aftershocks removed, Poissonian? BuII. Seismol. Sot. Am., 64: 1363-1367. Keihs-Borok, V.I., Knopoff, L. and Rotvain, I.M., 1980a. Bursts of aftershocks, long-term precursors of strong earthquakes, Nature, 283: 259-263. KeiIis-Borok, V.I., Knopoff, L., Rotvain, I.M. and Sidorenko, T.M., 1980b. Bursts of seismicity as long-term precursors of strong earthquakes. J. Geophys. Rea., 85: 803-811. Keilis-Borok, V.I., Knopoff, L., Rotvain, LM. and Allen CR., 1988. Intermediate-tetm prediction of occurrence times of strong earthquakes. Nature, 335: 690-694. Lamoreaux, R., 1982. Cluster patterns in seismicity. Ph.D. Thesis, Univ. AIber@ Edmonton. M&e, W.G., Rogers, G.C., Riddihougb, RP., Hyndman, R.D. and McMecban, G.A., 1978. Seismic&y of western Canada. Can. J. Earth Sci., 15: 1170-1193. Press, F. and Briggs, P., 1975. CbandIer wobble, ear&quakes, rotation and geomagnetic changes. Nature, 256: 270-273.
167 (1989) 299-312
Science Publishers
299
B.V., Amsterdam
Premonitory
in The Netherlands
seismicity
D. BROWN
patterns near Vancouver Canada
‘, Q. LI ‘, E. NYLAND
Island,
’ and D.H. WEICHERT
of Physics, lJniuersit_v of Alberta Edmonton,
’ Department ’ Geologml
- Printed
2
Alta. T6G 2JI (Canada)
Survey of Canada, Pacific Geoscience Centre, Box 6000, Sidney, B.C. V8L 4B2 (Canada) (Received
November
4, 1987; revised version
accepted
August
18. 1988)
Abstract Brown,
D., Li, Q., Nyland,
Canada.
E. and Weichert,
In: M.J. Berry (Editor),
The search
1989. Premonitory
Hazard
Assessment
seismicity
patterns
for large seismic events requires
and Juan
de Fuca
plates
M z 6.5. The
strongest
earthquakes
three
complete
generalized
for M 2 2.5 are preceded
sense depicted
by the algorithm
been preceded
by activation,
occupy
30% of the total
about
significant
similarity
by activation “M8”
suggested
have
of the earthquake
time-space
domain
the earthquake
of 5 years and 100 km. Such anticipation
could justify
attempts
to make short-term
using patterns
There is evidence
is a natural
consequence
predictions.
Analysis
for the coupling
of seismicity
of the model of a lithosphere
Introduction We seek to diagnose the approach of the strongest earthquakes in that part (Fig. 1) of the Pacific North American plate boundary where the Farallon plate fragments into the Explorer and the Juan de Fuca plates which then subduct under the
them. In the
catalog
of the world
future
strong
some precautionary of quiescence
where
earthquakes measures
restrict
has
been
range, In the 6.0 may have
(TIP) of a strong alarms”.
to be
appears
One event of magnitude
regions
as an assemblage
to define
flow in the lower magnitude probability
I!land,
167: 299-312.
the seismicity
are no “false
in the subduction
Vancouver
this threshold
there
to be tested, of anticipating
accuracy earthquake.
and
flow here and in other
which remains
since
by Soviet scientists.
considered
threshold
Canada)
occurred
near
Tectonophysics,
an a-priori Island,
but two others were not. The times of increased
between
tested. There is a possibility,
(the area of Vancouver
(M z 6.0) which
patterns
and Prediction.
area of the Explorer acceptably
for premonitory
D.H.,
Earthquake
This
algorithm
earthquake implies
a
MX was
in this area with an
and could be usef 11 for
the location
of the anticipated
zone with that above it. Such inter2 ction
of interacting
blocks.
the present state of seismicity analysis as a means of describing the phenomenon of earthquake flow and of suggesting hypothetical regularities. Empirical data are necessary to formulate such hypotheses,
and they are as insufficient
in ..he area of
in-
Vancouver Island as they probably ,ire everywhere. As a consequence we create hypotheses which might apply in diverse regions and use them
northwest of the of the southwest
to study the similarities in the situations before strong earthquakes. The approach se.lrches for
coast of Canada extending to the north end of Vancouver Island. For convenience in this paper we refer to this zone as the Vancouver Island area. Since there is no comprehensive theory which would define a priori an algorithm of earthquake prediction, a-posteriori data fitting is necessary in
those features of seismicity which are I-easonably independent of the tectonic environment. There is a quite surprising number that can be identified world wide. The approach used here considers tectonic information only in so far as seismicity data can be crudely associated with major tectonic
North
American
cludes United
a part States,
0040-1951/89/$03.50
plate.
The region
of the Pacific and that part
considered
0 1989 Elsevier Science
Publishers
B.V
iO(,
Once such a preliminary
features.
the approach
searches
for anomalous
to large earthquakes
without
strictions
related
intended
in this paper
merely comment
to tectonic
expressed
and CN, developed Institute
Moscow
based
specific
on
hypotheses.
re-
It is not
the patterns;
we
and reproducible
a relatively
wide
perceived
subregions patterns
free parameters applications
MS
that
For
world
lower thresholds equally
various
calculated.
This set of components
demands
a
set of criteria
set of premonitory
phenomena and makes effective use of the data available for the Vancouver Island area. Both algorithms recognize that the process of evolution to a large earthquake may be very different in different cases. Both algorithms are open in that the inco~oration of new data constitutes part of the aigorithm. Of the two we have chosen M8 as
a vector
and
describing
are a conse-
in time
window
of the
that
Marly of the
are fixed for all
P with
characteristics
of the catalog is divided
among
four groups.
For every
window
between
a threshold
in-
the components
M8 deals with seven characteristics
for each component
desig-
spaced
groups.
values
than
8 for which it was originally
every
OII
wide. in spite of the fact that
stants
Earth.
no hearing
in these regions. in the algorithm
M8 is now is use for much ned.
have
at the 0. Yu.
et al.. in prep.) reasoning
geographic
the magnitude
in the algorithms
for the Physics
of pragmatic
relatively
prior
artificial
by investigators
(Gabrielov
quence
patterns
making
to explain
is made
that they exist.
The hypotheses Schmidt
selection
is
into
divided
can be chosen
of the vector so that pS of at1
the start of a catalog
and the first
strong earthquake in the window, or the end of the catalog if there are no strong earthquakes in the window, are above this limit. All such values
the simplest.
are designated as extreme. We now examine the vector calculated at time t from t - 3 yrs to t and identify the components that were extreme at least
The M8 algorithm searches for process activation in a seismicity catalog and considers objects
once in this period. A group is extreme if at least one of its characteristics is extreme. We calculate
defined
by time moments
t and windows
centered
the number
of such extreme
characteristics,
on a small number of regularly spaced points in the territory to be analysed. These windows can be
the number of extreme TIP (times of increased
of various shapes and are a type of automatic division of the seismicity into possibly overlapping
(f. I + T) if for two consecutive of such extreme characteristics
Fig. 1. The study area.
h, and
groups, g, and declare a probability) for a period times the number is 2 6 and the
PREMONITORY
number
SEISMICITY.
VANCOUVER
of extreme
these parameters
groups
a pattern
applied
geodynamic
observations
time
record
can
series,
in the Vancouver
that indicate
phenomena, may
whose collective of increased
Although
(Press
et
1979;
Dobrovolsky
al.,
1980a, 1988: Gabrielov
MS is desig-
ers who
the approach
consult
wish
those references,
here. The problem
100 km and
a time accuracy
of a few
identify
for this algorithm
from a study of very large earthquakes
world wide
of a large
toring
We introduce functions
Keilis-Borok this
et
may
the essentials
earthquake
of time
in a region.
s and
in the catalog
consider by calcu-
lating
strong earthquakes, sufficiently complete
time t and all events in the interval [t-s. t]. First we separate the aftershocks from
M, and the start time of a catalog, tO. The algo~thm is
normalized to a fixed rate of production of mainshocks and can be applied in a very wide variety of situations if we accept the hypothesis of the self-similarity, world wide, of the earthquakegenerating process. The results of this study enhance our understanding of the underlying similarity of seismicity world
wide.
previously nate the
If such self-similarity
in fact exists,
unsuspected relationships may domiphysics of earthquake generation. In
addition to the study of self-similarity the analysis clarifies the way seismic volumes such as subducting oceanic slabs and overriding lithosphere may interact. Finally, the results aid in the prediction of destructive seismic events. A suggestion of a region in which there is an increased probability
of the catalog
that depend
than considered in the California and Knopoff. As a consequence
the
maximum
distance
at which
considered an aftershock the mainshock-aftershock
area by Gardner we increased the an event
could
be
by 50% (Table 1). Once separation has been
made by the application of magnitude-~tependent space-time windows (Table 1) we can consider the
TABLE
1
if a mainshock occurring
has magnitude
within
tude less than
T days
after
M or less, any event within the event and having
M is an aftershock
2.5
29
6
result.
3.0
33
11.5
3.5
39
22
Definitions
4.0
45
42
The basic assumption that underlies the pattern recognition approach is the assumption that the earthquake flow includes a large chaotic component. As a result, the specific symptoms of an approaching strong earthquake may be different each time. Accordingly the algorithm uses a set of
on
mainshocks using a modified form of the generally accepted thresholds (Gardner and Knopoff, 1974). Inspection of the events in the data that we used clearly indicated aftershocks at larger distances
of a large seismic event could be used to initiate more detailed studies using either other methods or more time-dependent data to refine the MX
of the functions of M8
to
by moni-
reduce the number of free parameters for the study to two. They are the threshold for defining
characteristics
al.,
and read-
recognition
flow of earthquakes a period
1975;
approach
we discuss
of time defined
is fully
Briggs,
is to use pattern
the TIP of a strong the total
and
et al., in prep.)
to implement
accuracy
derived
(TIP)
inbut
of the
of the algorithm
elsewhere
with a positional
The parameters
can be diagnostic
probability
the basis
described of
behavior
informative
earthquake.
but
Island
each of which taken
be insufficiently
earthquake
of a strong of about
patterns
be
of quantitative
since 1952. The algorithm
ned to identify
years.
longest
premonitory dividually time
approach
that can be used with M8 is the record
seismicity
of
of most of
earthquake.
recognition
to any multicomponent
the chronologically area
activation
In any case no TIP is raised for
3 years after a strong Such
is 2 4. The choice
guarantees
the characte~stics.
301
ISLAND
M
R
T
4.5
52
5.0
60
150 290
X3
5.5
70
6.0
81
510
6.5
91
790
7.0
105
915
7.5
121
960
8.0
141
9X5
R
a magni-
flow of mainshocks { m,
) where
34,.
B,(e))
mainshock.
m,
t is the origin and
that occurred
17,.
time such that I~. , > I,.
parameter
k
those
is the depth.
is the number
of
mainshocks
exercise
that occurred
normalize
in an interval
of
.s before a
greater
than a
set designated
mainshock.
of mainshock
of simply
summing
the events
which we designate
Z, and is as close to 20/yr
aftershocks
measure
as N,.
as possible
I.,.
which we consider that
occur
This function.
depends
on the choice
quakes
of
magnitude
as :V,. I_, . and
and
for another
Zz. The seventh
in M8 is the count
within referred
two
days
of
of
a
to as B below.
of -. M only in that earthless than M ._ are not
mainshocks. We can now consider the earthquake stream to be characterized at a time I by a vector of func-
If ,0 = b/3, then S is a measure of the rate at which the radius of circular fractures is generated, and if /I = 2h/3, then S is the rate at which fracture area is generated. Here. h is the usual coefficient in the energy magnitude relationship log(E) = A + bhf.(The parameter a is inserted to normalize the functions). We carry out this process for all events whose magnitude IV, is between M and an upper limit R (R= M, - 0.1) and the result as S(r; M, a,
s, (Y. p). This pair
of functions can be combined to produce a measure of the way in which seismicity concentrates in space. _M, M, s, (Y. p)
tions
deduced
from the preceding
R, .T)y
defines the ratio of the average radius of a fracture (S/N) to the average event separation for all events between M and R, if they are uniformly distributed in theregion (N”‘). The rate at which mainshocks are produced is not uniform. If we assume that on average changes in a linear fashion. the function -M, s) = N(r: M. t-f,,)
[ - 10 I -
2,) -
s
it
time period
.r.
The problem is to deduce. by pattern recognition, whether such a time is part of a TIP. The duration T of such a TIP. the time within which a strong earthquake is expected. is taken to be the 5 years used in all recent applications of M8 world wide. In order to declare a TIP at time f we require that over a period +, prior to and including I. three of the groups ( N,. N, }, ( L,. L, }. { Z,, Z, ). and { B ) contain functions that have extreme values
and that at least four of the functions N,. Z,. Zz. and B have extreme values. N2, L,. I.?. In this sense extreme values are those in the upper 10% of the range achieved by the function from the start of the catalog
S(K M, R, .s, a, P)
- N(r; M, I-ff,,-s)
by adjust-
set of functions
the number
activ-
I&
on the
to count
ity calculates a weighted form of this sum where the weights are proportional to the size of the
[/v@; ly, s) - A+;
Ilnear trend.
explicitly
the result and this is achieved
A more sophisticated
=
depend
o/’
M. If we expect to compare results with -. in other areas we must find some wav to
function
Z(r:
from a long-term
three functions
lower limit M. We refer to this count as N( I; _M. s).
define
the gencratmn
ing -M so that the average yearly rate of occmrence of mainshocks is in one case IO;‘vr for one
time I and which had a magnitude
earthquake. Instead we calculate:
deviates
to which
in the first e days after
the mainshock. It is a simple
the degree
mainshocks These
H(e)
number
measures
of the
+ is the longitude.
M is the magnitude
of \:ectors
( I,. A,. 9,.
has six components
where i is the sequence
X is the latitude, aftershocks
to bc a sequence
to the time of occurrence
of the first large (M > MO. the prediction threshold) event. Once this condition occurs the interval (I. t + T] is declared to be a TIP. There is no straightforward means of deriving this algorithm from more fundamental statements. It follows from the plausible idea that patterns should exist in the earthquake stream and that these patterns should, as much as possible, be physically realistic and flexible. In the final analysis its justification must rest on the pragmatic statement that it works in many areas of the world. An area such as the southwest coast of Canada can be considered as a collection of perhaps over-
PREMONITORY
lapping define
SEISMICITY,
regions
defined
latitude
is then
such proposed
region
a side of
of latitude.
Each
to check
by the catalog
regions
occur
per
that do not
per year are discarded
and
the M8 algorithm
is applied
determine
a TIP exists in any of them.
whether
to the remainder
to
once
of regionalization
a small
a large
is possiearthquake
If such an area is divided it is possible
or quiescent
event
a quiescent
to occur within similar
can
activation
quiescence
that may precede
of the latter. A
be seen on a 2-dimensional
of event number with
area
area the large event is likely
or on the boundary
process
to classify
in the 6 years prior
(Fig. 2). If an activated
against
time (Fig. 11). A pronounced activation
a large (M
in which
of the event
to the seismic contains
coupled
that contains
area
them as activated
more than one region for each region
an
of location
histogram
and a TIP does not always
to produce
into small cells of size of the order of the accuracy
The regions can overlap, and usually do. As a consequence the same large event can appear in appear
refinement
might occur is identified.
seismic
ten mainshocks
A further
of ble
to have
ten mainshocks
year. Those proposed contain
it may be possible
defined
is then examined
of at least
the regions
area of high probability.
The square region centered
over the time span covered
an average
We
the area
by equal intervals
R = e Mclp5-6 in degrees
whether
seismicity.
by first diving
cells defined
and longitude.
on this event
303
ISLAND
by their
such regionalization
into rectangular
length
VANCOUVER
of smaller scenario
magnitude
and
lack of large events events
is the
(Lamoreaux,
1982)
large shocks.
> M) earthquake. In the restrospective analyses of a catalog such results suggest geographic areas
The seismicity of the Vancouver Island area
that are sensitive to the onset of large earthquakes and may suggest areas that should be monitored with more care. By considering the intersection of
The seismicity of the west coast of Canada has been described by many authors (e.g., Milne et al.,
Elf. ISI
cl
Unknown seismic state Activated seismic state Quiescent seismic state AQ boundary
Fig. 2. A schematic
diagram
to illustrate
the idea of activation
quiescence.
52' N 51' N 50°N 49'N 48' N
46'N 45' N 44'N 32' W30° W 128°W 126' W 124' W 122'W
2o" w 118 W
Fig. 3. The seismicity of the study area. 2 00 !____!__--!----! ‘841
t
1
t86i
1
I
1881
1
1
(931
1
300
I
1
((341
1
1
1951
% 1
1955
1
1959
1 .t
l
1963
I
1
1967
1
1
197t
1
1
1975
1
1
1979
f
t
i9e3
1
1
._ _!_-__a
--_-,_.
-_
1 i
; i?
1
192q
6.00
5.00
4.00 __!____I____!
4 42 45
to
4
1
12 17 22 37 22
10 5 9 9
;
2 2 id 2
i 15
5 6 4 9 3 6 24 14 15 6 13 17 9 11 152 12 3
6 1
37 49 30 14 ; 61 28 24 7 19 33 3 23 17 15 18 R 19 13 3 26 26 20 14 4 39 52 24 16 6 22 13 5 15 4 9 9 3 4 7 9 8 5 25 11 42 15 40 4 18 7 78 43 71 42 23 55 67 7 25 129 72 10 24 5 77 48 2 9 3 1 2 23 ,____!____!____!____!_-__!_---!----!----!----!----!----!----
84 83
Fig. 4. A count table indicating the number of mainshocks
4 ; 5 4 4 i 1 2 1 3 5 2 2 1 4 1 7 1 13 1
3
4 4 8
2 2 2 1 2
;
1 2 I
i
; 1 1 1 4 ;
2
i
2
1
i
i
in a particular magnitude interval during a particular time interval for the entire catalog.
PREMONITORY
SEISIMICITY.
VANCOUVER
4.2 4.5 3.3
2.1
4.3 3.9 3.0
51°N. 50’N. 49-N .
3.8 6.3 5.1 4.6
48’N.
5.0 5.2 4.3 3.3 3.0
47*N -
4.6 4.4
46’N I 49N
2.6 5.1 5.0 3.4 4.8 R
3.6
4.8 3.9 4.6 4.3 4.8 4.3 4.3
4.4
3.7
5kT.6
4.1 4.5
3.7 I I
132; W 130&W 128”W
126&W 124”W
Fig. 5. The largest earthquakes
305
ELAND
122&W 120&W It%” W
in 1 o x 1 o rectangles.
1978; Basham et al., 1982; Basham et al., 1985). Figure 3 shows the geographic distribution of epicenters for events of magnitude 4 or greater in the area which we studied. In this area, from 118” to 130 o W and from 45 o to 51” N, we merged the seismicity catalog of the Geological Survey of Canada with the catalog from the U.S. National Earthquake Information Service. Particularly for the earlier years the Canadian catalog contained many significant earthquakes as far south as 45 o N. The duplicate events were removed manually. Usually, NEIS data were accepted in the United States and GSC data were accepted in Canada. If only one agency reported a magnitude the report of that agency was accepted. Events for which no magnitudes were reported were assigned a zero magnitudes while if a preferred mag~tude was given it was used. For these data this is usually Mr. We also removed those events since May 1980 that were clearly related to the Mt. St. Helen’s Volcanic Eruption. The data are summarized in Fig. 4 as a time-magnitude histogram and in Fig. 5 as a map showing the largest magnitude earthquakes in regions measuring 1” X 1 O. Patterns in space and time can be distorted by changes in the observation networks and associated changes in the degree of completeness of the seismicity catalogs at lower magnitudes. Evolution of the station configuration in the Vancouver Island area, and the associated increase in complete-
ness of the catalog at small magnitudes, has been discussed by Basham et al. (1982). The pattern recognition approach used in this work does not require such a stringent completeness condition as discussed there. Thus Fig. 4 indicates that since 1952 a sufficient number but certainly not all of the magnitude 2.5 earthquakes are present in the catalog for some parts of the test area. The effect of missing small magnitude events seems to predominate in the function that counts af-:ershocks and probably reduced its effectiveness as a diagnostic characteristic in Vancouver Island. SimiIar considerations apply to the time variable location accuracy of the catalog. The spatial resolution of the analysis in this paper is about 0.5 ‘, reflection epicentral accuracy. Focal depths are the least certain, and in some subareas nonexistent. However, the only application of focal depths is as an exercise to determine whether patterns in shallow seismicity precede an intermediate earthquake. Data analysis In the analysis of the subcatalog discussed above we primarily used the standard parameter set developed for M8 by an analysis of a large number of regions world wide. We adjusted the threshold for prediction A+, and the start time of the catalog to to satisfy completeness considerations (Fig. 4) as required to normalize the Vancouver Island earthquake flow to others used in MS and applied local judgment to adjust the windows used to eliminate aftershocks. Table 2 shows the nine strongest earthqu~es of the area and indicates that the threshold for prediction should be M, = 6.5. This choice is by no means unique but is a stable choice; small variations do not seriously affect our conclusions. In order to select earthquakes for prediction we chose, before beginning the prediction process, all mainshocks of magnitude greater than 6 that occurred since 1957. Since a mainshock of magnitude 6.4 that occurred in 1971 near 50 ON and 128O W was followed within 1 year by one of magnitude 6.0 near 49” N and 129O W we noted these events as possibly being predictive. The 6.0 magnitude is substantially lower than our threshold
306
I .4BLF. 2 ‘I he mne strongest earthquakes
m the region h
Year
Month
Da>
I872
12
15
5
1917
12
22
1918
12
1926
~-
min
-_
!.
Lat.
Long.
-._-__-_ I>epth Magmrude _____.__.-...__ 44
.7?
4X.6
121.4
15
4x
50.0
128.0
6.5
6
X
41
S.8
49.6
125.9
7.0
11
1
1
39
1x.0
4X.X
128.5
6.0
1939
2
x
6
39
25.X
49.1
12x.0
6.5
1946
6
23
17
13
25.X
49.8
125.3
7.3
1946
7
18
6
6
5x.5
49.5
129.7
6.5
1965
4
29
15
28
44.0
47.4
122.3
59
6.5
1976
12
20
20
33
12.0 -
49.0
128.7
18
6.7
for prediction so we place no great emphasis on the results for this event. We excluded two large aftershocks of magnitude 6.4 and 6.1 which occurred in 1972 and 1978, and two large events of magnitude 6.3 that occurred in 1954 and 1956 are not likely to be predictable since the preceding seismicity catalog is inadequate. Two events of magnitude 6.0 that occurred in 1958 and 1961 were ignored on the grounds that their magnitudes were both small and their time of occurrence was uncomfortably close to 1957, the beginning of a useful prediction dataset. We therefore set as our goal the prediction of one event near Seattle and
three west of the northern end of Vancouver Island. Figure 6 shows the TIPS that resulted from the analysis of all the seismicity data available to us. The heavy lines under the time scale indicate TIPS, the light lines indicate the center of the region in which the TIP occurred. The three events west of the northern end of Vancouver Island, of magnitude 6.0, 6.4 and 6.7, are associated with three regions in which TIPS occurred. Although the prediction for the magnitudes 6.0 and 6.4 should strictly be considered false alarms, we accept them as being within the uncertainty of
Fig. 6. The centers of the regions analyzed with M8. A heavy line on the pointer to a location indicates the occurrence of a TIP in the time interval read from the scale in the
upper right of
the diagram. In this analysis we used a magnitude
time in the catalog of 1952.
threshold of 6.5 and a start
PREMONITORY
SEISMICITY.
magnitude preceded
VANCOUVER
determinations.
The magnitude
by a TIP in two regions
the event, in another
region.
data
due
from
a region
for both regions The magnitude at a depth ceded
two years before
east of the event.
6.4 event
that contain
resulted
The
was seen
not seen in the region
east. The magnitude
by
It was not seen in the
6.0 event just to the northwest
in three TIPS and
6.7 was
and preceded
a 5 year long TIP, which expired
magnitude
307
ISLAND
to the
in a TIP
it.
mation The
regions
6.5 event near Seattle shallow.
their overlap
activated
at the bottom and
considered.
I
I
I
I
I
I
.6:0’,‘,‘,‘,‘,‘,‘,‘,‘,‘, \\\\\\\\\\ ,,,,,,,,,,,,,
by of
The space-
time volume is area (km2) multiplied
,,,,,,,,,, \\\\\\\\\\
common
by time (yrs).
I
are
fairly
speculation
the
intersections
common This
to all TIPS is shown
but
to explore tectonic
both the 6.5
does not prove
of the separate
the predictions,
In order
in that part of the
domains
it is certainly
more closely
activity
that
the
~~11 in fact sugges-
the potential
as a tool to predict
deeper events we attempted so separate the catalog into essentially crustal earthquakes and events related to subduction. Depth determirations in this area are not adequate for careful separation SO we
chose to discriminate by considering only events of depth known to be less than 25 km. This
I
Fig. 7. The regions covered by overlapping
large
whether
to each TIP. The part
of Fig. 7 and contained
the 6.7 event.
of shallow
/\1
regions
occurred
area. The quality of prediction is estimated space-time volume of all TIPS in percentage
here
ated with the TIPS all occurred
TIP did not appear argues strongly for coupling between the shallow and deep seismicity in this
volume
considered
of the sum of their area defines
that pro-
Nevertheless
esti-
region of high risk. We note that the events associ-
refine tive.
space-time
in Fig. 6 gives a quality
(Fig. 7) and it is an interesting
the existence of four successful TIPS, one almost successful TIP, and only one region in which a
total
shown
of 30%.
of the time domains
of 59 km but the activation
the event was largely
The analysis
TIPS and one time intervals
common
to all TIPS for each event.
308
1960
Fig. 8. The results
of an analysis
using
only
events
whose
depths
explanations
is clearly a set of Iithospheric events but does not include all of them. The results of an analysis of this subset are shown in Fig. 8; there are insufficient events in the remainder of the set to allow application of M8. As expected, since most TABLE
have
1970
been determined
1990
1960
IO be 25 km or less. See Fig. 6 for
of layout.
seismicity near Seattle is shallow, the TIPS around the Seattle event are unchanged. One region to the south of the 6.7 event is lost due to insuf~~ient events in its s&catalog but the others remain viable TIPS.
3
A summary
of the TIPS (times of increased
probability)
for the Vancouver
TIPS
Center of window rp
x
49*N
129ew
Island area
Strong earthquakes
Time
Date
48”N
128” W
128* W
~970,6.20-1977.6.2~
1970,6.20-1978,6.20
i 1971.3.13
6.7 6.4 I
1971.125
6.0
1976,12.20
6.7
i 19?1,3.13
6.4
1971 J2.5
6.0
( 1970,6.20-1975q6.20
72
32
1971,125
6.0
6.0
1983,6.20-1989.6.20
?
?
72
19654.29
6.5
72
48’N
124* W
1959,4.29-1%7,4.29
49ON
123O W
1957,10.29-196410.29
48*N
122OW
1960,4.29-1965‘4.29
49ON
121”W
1957.10.29-1964.10.29
47“N
121°W
( 1962,4.29-1967.4.29 1982,4.29-198924.29
* Time period
M
of
expectations (months)
1976.12.20
SOoN
Tie
from start time of a TIP to a strong
event of M > Me.
84 1%5,4.29
6.5
60 84
19654.29
6.5
36
?
?
84
*
PREMONITORY
SEISMICITY,
VANCOUVER
309
ISLAND
Center of window: 49N 129W .*****A * * * * ~*,***~,~*,“,*,* *I************* I***h*F******* ***********~* .************** .*****~******~************* *****A********* **********~***************** .***********I** I*************** .************** .**********~@t.&************ **************I .*********I**** A************* r*# * A * * * * * *r”*a***** ******** 1
Nl N2 Ll L2 Zl 22 B 1965
1970
1975
1970
1975
Center of window: 50N 128W
1965 Center of window: 48N 124W Nl N2 Ll L2
.* ** * ** **** * **-??**************** .********* * * * * * * * * I) *%.********** .**************** *******n******** .**************** a ‘** * * * * * * * I) * * * .*******n******** ********J******* .***************I **************** .**************** * * * * * * * * *2” A * * * * * ;c**1****** * * * *y&#.&* * * * *************** * * * * A* * **************** ************I*******~* * * **’ ‘* .q#* n * * * * .******** *** 1
Zl 22 B 1955
*
1960
1965
Center of window: 48N 122W Nl N2 Ll L2 Zl 22 B
‘*********3h**L**** .********* *F * .* *I) ** ** ** ******* * * * ***** .*********’ ********* ?* * * * * * * * * ********* ********I ********** .********* ********** ‘******y&+******* .**I** * * * * I * *?,**** .*I******* ************#r ‘L. * * * * * * * *** x
1955
1960
1965
Center of window: 47N 121W Nl N2 Ll L2 Zl 22 B 1955
u.********, * ****** * * * I) I ******. * ********I * * * * * *I * . ********* .*********. ****~*****~*******. ********* ***~*********** * * . I)********. **A****** *********. ******xArte******. ********* *I)*******. ********* ****~*~* ***** * * * * . 1960
1965
Fig. 9. The part played by each function in establishing the TIP. The shaded area are the time intervals in which the activation was sufficient to establish a TIP (referred to here as a “ vote”).
The analysis that exploits the entire dataset is associated with two TIPS for which large events have not yet occurred.
In one case the seismicity
the anticipated earthquake. We and others found it helpful to examine activatu~n on a more detailed
scale and
of
we show some of the regions
subdivided
than 6.5 began in the middle
of
cells. We have examined
1982 and will expire in 1991. If past experience
is
significant
activity
years prior
to the large event.
If no such activity
occurred
the cell was identified
as quiescent.
collection
of quiescent
data
indicate
magnitude a guide
that
greater this TIP
a TIP
is associated
end of Vancouver
Island.
ated with the northwestern States.
It began
expires
in 1991.
Some
traits
of an
earthquake
quiescence
h;t\e and
with the northern
The other TIP is associcorner
in the middle
of the United
of 1983 and
also
of the catalog
contribute
to the the
time intervals that have created TIPS for five windows. The asterisks indicate the occurrence of extreme values for the functions on the left at times
indicated
at the bottom
of the plots.
The
regions are the time intervals in which the activation was sufficient to establish a
TIP. Table
may
strong
TIPS and others do not. In Fig. 9 we illustrate
shaded catalog
cells
3 is a another
means
of summarizing
these results. As indicated in Fig. 6 the M8 algorithm is not particularly effective in specifying the location of
Cantor ot wlndow UN 12oW
Canter of wlndow: SON 128W
the
event
events
on
in it durrng
cells bounded the
will occur. plot
into 256 iI
each cell to determine
occurred
indicate
I-‘& It!
in
the
area
the o A
by activated
in which
the
next
In each of the four cases boundary
between
quiescent and activated areas, but the regions the northwest show a better defined area
the in of
quiescence than those in the southeast. The same process can be used to refine the time of occurrence of the anticipated earthquake. In time, the quiescence when
an earthquake
histogram time for
will help us to refine the time will occur.
Figure
11 is a
of event number against magnitude and the magnitude 6.7 and 6.5 events. A
quiescent period occurs before the magnitude 6.7 event and a courageous analyst would have predicted the event when the activity began to increase. The histogram of magnitude not show such an effect as clearly.
6.5 event does
Illllflllllllll!l
Conclusions
It is scarcely surprising that the seismicity of Vancouver Island bears a strong similarity to that found elsewhere. It would have been surprising if
Centor
of wlndow:
4ON 122W
Cmtw ot wtmlow 48N lZ4W
llllllllI_ccI
q ahsceflt
Vancouver Island prior to at least four large seismic events. This suggests the possibility of
cdl
El Ouioscont boundwy .
this area were radically different from the twenty areas previously tested by Keilis-Borok and his colleagues (Gabrielov et al.. in prep.). Apparently identifiable patterns exist in the seismicity of
Loutton ot strong wrthqwk*
Fig. 10. Plots of activation
and quiescence
events in the catalog.
for the four large
defining time and space intervals that indicate an increased probability of large earthquakes. Three successful retrospective predictions alone are not convincing, but the agreement of this method with many cases world wide suggests that the current TIPS in the northwest and the southeast parts of the study area should be taken seriously. One earthquake of magnitude 6.0 is included in the TIP that precedes a larger event, but two others are not. This indicates that seismicity of
PREMONITORY
SEISMICITY.
VANCOUVER
311
ISLAND
Center of window: 5ON 12’7W 4.00
2.00 3.00 !____!____!____!___ 69
1
70
1
I
0
f
3
0
1
71110 72
1
73
1 IO
74
1
75
1
i I
1 0
1 1 1
10 10
76
1
10
77
1
10
4 5 4 8 8
78110
---
!----!
3 9 166
6.00
5.00
.- !----!----!____!--__!____!
2 7 1 3 1 1
1 3 1 3 2 2
2
; 1 2 1 2 1 6 6 1
2 ; 1 2 2 2 3 2 1
i
; 1
; 2
ii
2 3
2
1
3 1
;
i 2 1
2 2 I
i : i 1 3 1
:
:
:
1 1
;
: I
: 1
i 1 1 2 1 1
:
7.00 _!____!___-
1
i
:
2 1
: 1
:
:
8 13 7 12 6 8 4 2 5 5 6 9 6 15 18 12 14 16 0
-!__--!____!--__!___-I___-!___-!
4 1
3 3
1 2
1 6
2
2
CVENTS
Center of window: 47N 123W 2.00 3.00 !____!____!____!-_ 58
1
1 0
59
1
10
60
1
I
61
1
10
62
1
1 0
63
1
1 0
0
13 10 9 8 8 9 6 6 5 4 2
64110
7 3 4 6 I 2 7 3 1 3 2 1 1
6.00 -!----!
-_-
-_-
5 3 2 6 2 1 1 2 1 1
:
1
0
26 16 15 24 II 13 14 13 7 5 6 5 2 4 2 3 0
1 1
1
1
65
1
1 0
i 1
2
:
8 2
4 4
2 6
66110
166
3
4
3
2
1
0
0
EVENTS
Fig. 11. Magnitude-time
histograms
of events in regions
centered
on two large earthquakes
(see text).
magnitude 6.0 is evidence of activation for events of larger magnitude. but that the activation phenomena required to predict events of magnitude 6.0 by MS do not exist in this catalog. The two events of magnitude 6.0 that are outside the TIP were not followed by large events but the one that was within the TIP was present there because the TIP was initially generated for an event of M = 6.4. Since this M = 6.0 event was foflowed by an event of M = 6.7 and the other two were not followed by large events, these smaller mainshocks were correctly interpreted by the algorithm. The failure to predict them is not a failure of the algorithm; in this area it appears to work only for magnitudes in excess of 6.4. In addition we observe quite clear patterns in shallow seismicity prior to the deep large event at 47O30’N, 122” W depth 59 km). This suggests that the subducting lithosphere may be coupled to shallower structures in some as yet unexplained way. This in turn implies that other near-surface geophysical measurements, such as strain and geochemical studies, may have potential as detectors of precursory patterns to strong deep earthquakes. At this time we make no statement on the potential for a megathrust earthquake in Cascadia. The interval of adequate seisrnicity catalogs may not be sufficiently large to identify this event. The only evidence that may have a bearing is the strangely low number of mainshocks in the area and the anomalous difficulty we encountered in assessing aftershocks. AcknowIe%lgements The analysis reported here was directed by V.I. Keilis-E3orok while he visited the University of Alberta in July and August of 1987. It benefitted greatly from discussions with Garry Rogers of the Pacific Geoscience Centre, Ge&gica.l Survey of Canada. We also acknowledge our debt to a legion of workers who developed the coding for the M8 algorithm.
Keilis-Borok’s visit was supported t’rom the Fund for the Future at the University of Alberta. The Natural Sciences and Engineering Research Council of Canada. the Department of Physics. and the Departments of Computer Services and Computing Sciences at the University of Alberta also aided this research in various ways.
References Basham, P.W., Weichert, D.H., Anghn, F.M. and Berry, M.J., 1982. New probabilistic strong seismic ground motion maps of Canada: a CompiIation of earthquake source zones, metbods and results. Minist. Energy Mines Reaourc., Earth Phys. Branch Gpen File 82-33. Basham, P.W., Weichert, D.H., &$n, FM. and Berry, M.J., 1985. New probabiBstic strong seismic ground motion maps of Canada. BuU. Seismol. Sot. Am., 75: 563-595. Dobrovobky, I.R., Zubkov, S.I. and Miachkin, V.I., 1979. Estimation of the size of earthquakes preparation zones. Pure Appl. Geophys.. 117: 1025-1044. Gabrielov, A.M., Dimitreava, O.E., Keihs-Borok, V.I., Kosobokov, V.G., Kuznetsov, I.V., Levshina, T.A., Mitroev, K.M., Molcban, G.M., NegmatuBaev, SKh., Pisarenko, V.F., Prozeroff, A.G., Rinehart, W., Rotvain, I.M., Shebalin, P.N., Shnirman, M.G. and S&eider, SYu., in prep. Cent. Reg. SismoI. Am. Sud., Lima. Gardner, J.K. and Knopoff, L., 1974. Is the sequence of earthquakes in Southern California, with aftershocks removed, Poissonian? BuII. Seismol. Sot. Am., 64: 1363-1367. Keihs-Borok, V.I., Knopoff, L. and Rotvain, I.M., 1980a. Bursts of aftershocks, long-term precursors of strong earthquakes, Nature, 283: 259-263. KeiIis-Borok, V.I., Knopoff, L., Rotvain, I.M. and Sidorenko, T.M., 1980b. Bursts of seismicity as long-term precursors of strong earthquakes. J. Geophys. Rea., 85: 803-811. Keilis-Borok, V.I., Knopoff, L., Rotvain, LM. and Allen CR., 1988. Intermediate-tetm prediction of occurrence times of strong earthquakes. Nature, 335: 690-694. Lamoreaux, R., 1982. Cluster patterns in seismicity. Ph.D. Thesis, Univ. AIber@ Edmonton. M&e, W.G., Rogers, G.C., Riddihougb, RP., Hyndman, R.D. and McMecban, G.A., 1978. Seismic&y of western Canada. Can. J. Earth Sci., 15: 1170-1193. Press, F. and Briggs, P., 1975. CbandIer wobble, ear&quakes, rotation and geomagnetic changes. Nature, 256: 270-273.