Wavelet transform as a tool for detection of geomagnetic precursors of earthquakes

Phys. Chem. Earth, Vol. 23, No. 9-10, pp. 965-967, 1998 © 1998 Elsevier Science Ltd All rights reserved 0079-1946/98IS--see front matter

Pergamon

PII: S0079-1946(98)00128-1

Wavelet T r a n s f o r m as a Tool for Detection of G e o m a g n e t i c Precursors of E a r t h q u a k e s L. A l p e r o v i c h 1 and V. Z h e l u d e v 2

1Department of Geophysics and Planetary Sciences, Tel-Aviv University, Tel-Aviv, Israel 2School of Mathematical Sciences, Tel-Aviv University, Tel-Aviv, Israel

Received 25 April 1997; accepted 22 December 1997

Abstract. We applied discrete wavelet analysis to US geomagnetic digital network data to reveal ultra low frequency (ULF) waves prior to the occurrence of strong earthquake, such as M=7.0 Loma-Prieta, California earthquake. Via the wavelet approach we s n ~ in discovering the following phenomena: 1. Geomagnetic oscillations of specific scale time within the range 10 nfin- a few hours were generated in the seismoactive region 2 days and 5 hours before the LomaPrieta earthquake. 2. The oscillations were running from the source in the two modes: as a result of earth currents extending from a seismic source (hourly oscillations) and as slow hydromagnetic waves (short packets of 3 + 10rain oscillations).

f Ic

© 1998 Elsevier Science Ltd. All rightsreserved.

1

big. 1. Digital Geomagnetic Network observatories:BRW-Barrow,Alaska;

Introduction

BOU-Boulder,Colorado;BSL-Bay St.Louis;DLR-Del Rio, Texas;FRDFredetikshag, Virginia; FRN-Fresno, California; HON-Honolulu, Hawaii; NEW-Newport, Washington, SIT--Sitka, Alaska; SJG-San Juan, Puerto Rico; TUC-Tueson, Arizona.

The objective of the investigation presented is to search for the low frequency electromagnetic perturbations preceding an earthquake and initiated by it. We used one-minute sampiing rate recordings of the total field from 12 US Digital Geomagnetic Network observatories. The location of these observatories is shown in Fig. 1. In accordance with the 1 min. sampling rate mentioned, the beginning of the earthquake is addressed to the sample number 4323. Simultaneous centered magnetograms of total vector variations from 8 observatories are shown in Fig. 2. The initial point of Fig. 2 is related to 00h00m UT on October 15, 1989. WRhin the 5 day period (15/10/89 -19/10/89), a strong Loma Prieta earthquake had occurred ( San Francisco, October 18, 1989 at 00h 04m UT).

1 FRN, 2 T U C , 3 BOU, 4 DLR, 5 BSL, 8 FRD, 7 SJG 8 NEW, 9000

•

8000

i

7

i

!

!

7000 6000 5000 •

4000

:

3000

! !3! _

2000

2

151

1000

Remarks on mathematical tools employed

0

Up to the recent time the main tools for the search for electric, geomagnetic, ionospheric and other precursors remained

..... ! :

~

; . . . . . . . .". . . ~........ ~........................

:

iii

i

-

Fig. 2. Simultaneous recordings of the geomagnetic field from 8 observatoties

Correspondence to: U Alperovich 965

966

L. Alperovich and V. Zheludev: Wavelet Transform for Detection of Geomagnetic Precursors of Earthquakes

the direct comparison of signals recorded on a net of observatories or of the signals subjected to a Fourier processing. However, the Fourier analysis is poorly suited for problems concerned with detection of low-magnitude short-living events, especially in the absence of an a priori information on the duration and frequency content of perturbation sought. To the contrary, the Wavelet Analysis (WA) and its extensions, such as wavelet packet analysis, for example, provide relevant means for such cases. WA appeared primarily in geophysics (Goupilland et al., 1984) and now it is a powerful tool in numerous areas of fundamental and applied research. In some sense, it is a drastic extension and generalization of the Fourier analysis. In the WA, signals studied are expanded into a unit of waveforms well localized both in time and frequency domains, unlike sine (cosine) waves used in the Fourier analysis which may be looked upon as a limit case of WA. An immense collection of "testing" waveforms supplied with efficient computational procedures is available in the modern WA. Moreover, current techniques of WA provide means for flexible adaptation to a signal analyzed (Wickerhauser, 1994). This yields relevant tools to search for localized perturbations shadowed by the noise background. Our basic mathematical tool was the wavelet packet analysis. First, by primary investigation of the class of signals available, we choose a family of the analyzing wavelet packets. This family is ranged by scales (durations), blocks, which correspond to frequency content, and, within a single block, by time location. Then we computed the correlation coefficients of a signal studied with these testing waveforms. The components of the signal most matched to some waveforms from the family contributed the large coefficients to the corresponding waveforms whereas random incoherent noise background produced small coefficients with any waveform. The next step was the operation of thresholding that is elimination from the table of coefficients those ones whose magnitudes did not exceed some prescribed threshold level. In doing so, we obtained comparatively few significant coeffidents associated with various scales, blocks and locations, which present a peculiar signature of the signal. To observe changes of components of a signal from one time scale to another, to trace those in the space, to compare with each other ones, it was helpful to implement reconstruction from several blocks. Such reconstructed blocks comprised contribution to the signal from corresponding testing waveforms. It allowed to extract isolated events and to follow their propagation from one observatory to another.

1 FRN, 2 TUC, 3 BOU, 4 DLR, 5 BSL, 6 FRD, 7 SJG 8 NEW,

te

......

12

Ti

,

...d, ,~

~ ,

. . . . !~I....... tl,~ . .llj,

. . . .

.....

......

,

,1

i

. . . . . .

.

' "' '

4 "":4 "'i . . . . '

,I .....

: ....

.....

'

'i'''

" II

:

6,519)

4,851) '

i

:

8.92e)

' 3,412) 2,31s)

Ill

! o

1ooo

,~

3ooo ~ ,ooo ~ DecompotJtlon scale - J-6; Block -2; S'p'nmletS, Nothre~hold

,ooo

Pig. 3. Time dependenceof waveletcoefficientsof the signals shown in Fig. 2. Analyzingwaveletis Symmlet8with time-scaleabout 1 hour. Thresholdingwas performedat the level 130

3

Description tained

o f t h e i n v e s t i g a t i o n a n d the results ob-

After testing numerous classes of wavelet types, we chose for our analysis the Symmlets of 8-th order and Spline wavelets of 4-th order (Wickerhanser, 1994) as most relevant to the problem being investigated. Fig. 3 shows time-dependence of wavelet-p.acket coefficients at the 2rid block from 6th scale of decomposition oF signals from 8 observatories for the Symmld; 8 wavelet packets. Fig. 4 displays the reconstruction of the block. We point out that in the picture, as well as in subsequent ones, the first number inside brackets corresponds to the index of an observatory, whereas the second one indicates the magnitude scaling factor. The observatory nearest to the site of the earthquake is Fresno observatory (FRN). Fig. 5 displays results of processing the signals by the stationary wavelet transform. We emphasize on the fact that with the threshold level 1 FRN, 2 TUC, 3 BOU, 4 DLR, 5 BSL. B FRD, 7 SJG 8 NEW, (8,170) (L77) S.88) 5,1t~) 4.10e)

To determine exactly locations of events of interest we used, beside the wavelet packets, the stationary wavelet transform (Nason and Silverman, 1995). Here the set of testing waveforms was arranged by scales and locations. Waveforms associated with a certain scale were constructed as an appropriate dilation and shifts along the time axis of a single basic wavelet. It is important that, unlike the conventional WA, the adjacent waveforms at each scale were shifted relative to each other at one sampling point.

3,tl9) (2,,~) (1,se)

7000

R~on=truction from scale - J-6, block -2; SymmkltS, No threshold

Pig. 4. Geomagneticperturbationsreconstructedfrom the waveletcoefficientsdisplayedin Fig.3

967

L. Alperovich and V. Zheludev: Wavelet Transform for Detection of Geomagnetic Precursors of Earthquakes I FRN, 2 TUC, 3 B O U , 4 DLFL 5 n s L , 6 FRD, 7 SJG 8 NEW,

FI::IN.2 TUC, 3 BOU, 4 DLR, 5 B~L, ElFRD, 7 S,,IGI8 NEW,

oi ,

16

.

i

i

S,4459)

•

Y I

6,1720)

I

12

Iu I

ill

1o

!

5,2437)

12 •

- -

•

ll:i i

lO

'

I

:

.

| :

:w

I:

r

s,1657)

I I

I 2,1640) !

equal to 1500 at the FRN observatory, two areas with the highest coefficients are not observed at any other observatory (see Fig.). The first area is located near the 1000th minute, and the second one is just before the earthquake (4323th minute). The perturbations are most apparent within the chosen range of the wavelet durations at the FRN (with tlareslaold level at 1500). It is interesting, that the perturbation is observed practically at all observational points for thresholds in the range 1000,1200, excluding DLR, but the coefficients attain its maximal intensity at FRN. In the ease when geomagnetic perturbation caused by processes of preparation of the earthquake propagated from the seismoactive area at a slow velocity, the WA should produce significant wavelet coefficients for different observatories at different times, corresponding to the times when the perturbation reached observational points. We applied the procedure of total comparison to the same set of stationary wavelet coefficients as before. Namely, the stationary wavelet coefficients are ranged accordingly to time moments. At any moment corresponding coefficients from all observatories are compared. Only one whose magnitude exeee~ all other's magnitudes is captured and addressed to the related observatory. All except this coefficient are eliminated. At the next moment the procedure is repeated and an observatory with maximal coefficient gains a new term. Such a procedure is helpful for tracing the wave propagation from one site to another. Fig. 6 displays the result of the total comparison. From Fig. 6 we recognize distinct anomalies within the time interval 1000 + 2000 rain., immediately prior to the earthquake at t=4000 and just after that, at t=4500. We carried out the procedures of WA for the whole 5 day interval containing the the moment of earthquake using various families of wavelet-packet testing functions. The conclusion about an elevated geomagnetic activity captured by comparatively long wavelets ( 3hour) nearby the epieentral zone (FRN, TUC) still stands and does not depend on the waveletpacket form.

P

i I,,

i

t7.17~)

I

i ! ,

!

6,172O1 II(5,2437) 4.1787)

i I Ill

3,1657)

I ,

I: • . I r ,4OO0

le.44sg)

,

{.ill,

i

•

V

I

Iii i i I Y

5OOO

dIComperisorl; Decompc~#ion scale - J-8; Spline4

Spline4, Threshold ~1500

Fig. 5. Thresholdedstationarywavelettransformof the signalsshownin Fig.2. Thresholdingwas performedat the level1500

'

i

•

I 7OOO

c)Decomposition scale - J - B ;

L

-|

I

t:

I

I

i

•

~V &

i I i |,

t

i :

1,1~2)

I

'

il

I I

mv'a

4,1787)

:

I

'

7,1712)

I

14

le

6OOO

2,1~45)

1.1842)

I 7000

Fig. 6. Resultof the operationof totalcomparisonappliedto the stationary wavelet~ansformof the signals.

4

Conclusions

The procedure of total comparison revealed an increase of the geomagnetic activity prior two days to the earthquake and immediately at 5 hours before that at the observatories nearest to the earthquake. The block-wise reconstruction allowed to study the thin structure of the anomalies discovered. It was found that time-scales of the perturbations are within the ran.ge from 10 rain to a few hours. Spatial features of the oscillations with maximal intensity at the low-latitude region near the epicentral area (obs. FRN and TUC) and the direction of its propagation from that area differentiate them from the magnetospheric perturbations. Maximal intensity of the latter ones tends to high latitudes and the oscillations propagate from the high to middle and low latitudes. The investigation revealed at least 3 kinds of wave packets: 1. Wave packets comprised the hourly oscillations running with extremely low speed at V ~ 0.1 + 0.5kra/see. Most likely, propagation of 5 hours perturbations could be thought of as a result of spreading the electrical current in the ground with finite conductivity from a local seismic source. For the simple homogeneous ground model with apparent resistivity p = 103Ohm • ra and perturbation with the time-scales ~" = 2. 104see we had V ,,, 1 V / ~ / r ,,, 0.5kin~see. 2.The second type of the revealed perturbations were short packets of 3+10min oscillations whose horizontal velocities were about l+2km/s. 3.The 3-rd type were fast oscillations with the horizontal velocities at V >> 2 k m / s . References

R Goupilland,-A.Groc;smannand J. Morlet, Cycle octave and related t~ansformationin seismicsignal analysis,Geoerploration, 23, 85-102, 1985/84. W.V.Wickerhauser,AdaptedWaveletAnalysisfrora Theoryt o Software, AK Peters,Wellesley,MassaehuseUs,1994. G. R. Nason and B. W. Silverman,The StationaryWaveletTransform and someStatisticalApplications,/nWaveletsand statistics, (A. Antoniadis and G. Oppenheirneds), SPRINGER-VERLAG,N.Y., pp. 281300, 1995.