A frequency-based parameter for rapid estimation of magnitude

This study introduce a new frequency parameter called sfcwt, which can be used to estimate earthquake magnitude on the basis of the first few seconds of P-waves, using the waveforms of earthquakes occurring in Japan. This new parameter is introduced using continuous wavelet transform as a tool for extracting the frequency contents carried by the first few seconds of P-wave. The empirical relationship between the logarithm of sfcwt within the initial 4 s of a waveform andmagnitude was obtained. To evaluate the precision of sfcwt, we also calculated parameters s p and sc. The average absolute values of observed and estimated magnitude differences (jMest Mobsj) were 0.43, 0.49, and 0.66 units of magnitude, as determined using s p , sc, and sfcwt, respectively. For earthquakes with magnitudes greater than 6, these values were 0.34, 0.56, and 0.44 units of magnitude, as derived using s p , sc, and sfcwt, respectively. The sfcwt parameter exhibited more precision in determining the magnitude of moderateand small-scale earthquakes than did the sc-based approach. For a general range of magnitudes, however, the s p -based method showed more acceptable precision than did the other two parameters.


Introduction
Given the occurrence of destructive earthquakes in seismically intensive regions, the resultant human and financial losses, and the scientific inability to predict earthquakes, investigations into earthquake early warning systems (EEWSs), and their development are critical requirements. Accordingly, researchers have recently proposed many approaches to improving EEWS performance [1,3,4,12,14,15,17,26,27,29]. Earthquake records generally contain valuable information on seismic sources [20]; records of the first waves of an earthquake contain information on magnitude, and those of second waves carry information on energy [2,12]. For this reason, source parameters (magnitude and location) in EEWSs are estimated through an analysis of the initial phases of Pwave arrival. Because the frequency content carried by the initial seconds of a signal does not depend on epicentral distance, earthquake magnitude can be directly estimated using frequency parameters without correcting distance effects [6,24,25,28]. This method enables the derivation of empirical magnitude-scaling relationships in EEWSs on the basis of frequency information from the first few seconds of P-wave arrival [1,2,6,12,16,[21][22][23][24][25]28]. In this regard, two basic parameters, namely, s max P and s c , were introduced for use in EEWSs in California [5,13]. In the s max P -based method proposed by Nakamura [16], s P is measured through acceleration and velocity time series within the initial seconds of receiving P-waves (usually 3-4 s); the calculation proceeds as follows [1]: where i denotes the number of samples, x i is the ground motion velocity within a sampling interval, X i represents the smoothed ground velocity squared, D i is the time derivation of the ground velocity (acceleration), and a

A r h i v e o f S I D
stands for the smoothing regularization constant. Variable a is set to 1 À ð1=F s Þ. Allen and Kanamori [1] considered s max P the maximum value of s P within 2-4 s after P-wave arrival. The s P measured within the initial seconds of Pwave arrival can be disregarded, because this parameter depends largely on pre-signal noise level [18,19] [21]. Frequency parameter s c , proposed by Wu and Kanamori [24] for scaling earthquake magnitude, is calculated by thus: where u and v denote the displacement and velocity records of ground motions, respectively. The main weakness of s max P and s c theories is that they are justifiable for noise-free monochromatic signals [12,18,23,29], yet an earthquake time series encompasses a wide range of frequencies characterized by unavoidable noise. Shieh et al. [21] indicated that s max P is sensitive to the initial noise level of a signal. In an analysis of the underlying theory and characteristics of s max P , Wolfe [23] demonstrated that this parameter is a non-linear function of spectral amplitude and that it assigns considerable weight to high frequencies.
Magnitude estimation grounded in s max P is, therefore, probably associated with significant errors. As Wolfe [23] proposed, wavelet and multi-taper concepts can be used to decrease the weighted dependence of s max P on frequency contents and precisely determine the correlation between s max P and earthquake size. To define a new frequency parameter, namely, the log-average period (s log ), Ziv [29] used the power spectrum of a signal in the initial seconds of P-wave arrival. s log is measured on the basis of the early proportion of the velocity spectrum as follows: (a) the initial few seconds of P-wave arrival are extracted from the total waveform, (b) the effects of a windowed signal on two ending edges of the signal are eliminated using a Hanning window, (c) power spectrum coefficient P i x i ð Þ is calculated on the basis of Fourier transform, and (d) the s log value is measured through the following equation: where * shows that P i x i ð Þ is resampled per 0.1 log unit of frequency within frequency ranges of 0.1-10 Hz. At the end of the calculation, the logarithm of the average s max P , s c , and s log measurements of each event are plotted against earthquake magnitudes, and the best line is fitted using the least-squares method as follows: Earthquake early warning systems use this logarithmiclinear relationship in estimating earthquake magnitude. In this relationship, a and b denote linear regression coefficients, and s and M denote the frequency parameter and magnitude, respectively. We introduced a new magnitudescaling parameter called s fcwt , which is based on timefrequency continuous wavelet transform (CWT) as a combination of spectral methods and wavelet transform.

Data
The data catalogue treated in this work consisted of vertical component accelerograms associated with 72 earthquakes with magnitudes (M JMA ) between 3 and 8 and epicentral distances smaller than 100 km. The data are also related to the borehole accelerometers recorded in the KiK-net strong-motion seismograph networks of Japan's National Research Institute for Earth Science and Disaster (NIED) Prevention from 2000 to 2016 (Fig. 1). Vertical borehole records were used to eliminate near-surface effects. To reduce undesirable effects from background noise, records with high signal-to-noise ratio ([ 5) were chosen among the vertical components. This ratio indicates the power of the main signal-to-background noise power. It is clear that a greater amount of this ratio represents the high quality of the seismic signal. Moreover, the earthquakes chosen were those occurring from inland and oceanic rim events with focal depths of less than 50 km. The data set likewise

Method based on s fcwt
Similar to Ziv [29], we used the power spectrum of a signal within the first few seconds of P-wave arrival to determine the frequency features of the signal. However, Ziv [29] used Fourier transform to determine the power spectrum, similar to Atefi et al. [3], we adopted wavelet transform coefficients C(k; sÞ for such calculation [3]: in which w k;s ðtÞ 1 ffiffi The wavelet transform of a time-frequency f t ð Þ (here, the initial few seconds after P-wave arrival) is defined as the inner product of f t ð Þ and w (wavelet family), resulting in wavelet coefficient C. This time-frequency transform of f t ð Þ was formulated as Eq. (8) [7][8][9][10][11]. In this study, the average of all the samples of the wavelet transform coefficients was used to estimate the power spectrum at each scale of k 0 . In general, s fcwt values are measured as follows. (1) After the removal of a trend and mean, a Butterworth bandpass filter with cut-off frequencies of 0.1 and 10 Hz is applied to acceleration signals, and then, 4-s time windows are eliminated after obtaining P-waves from records. (2) To eliminate sudden cutting effects on two ends of a windowed signal, a cosine taper of 0.05 is applied on the signal. (3) Windowed signals are decomposed using wavelet transform (Eq. 8), and C(k; sÞs are calculated in the 1024 scale of k. Daubechies-12 is used for the frequency analysis of the signals. Subsequently, (4) a pseudofrequency is calculated using the f = F c /(k 9 D) relation, in which F c is the central frequency of the wavelet function, and D represents the sampling periods. The average power spectrum for each pseudo-frequency is calculated using the following equation: where N denotes the number of Cðf i ; sÞ coefficients. (5) Finally, s fcwt is the inverse of the pseudo-frequency, which has the highest value of Pðf i Þ. To obtain the empirical magnitude-scaling relationship, the logarithm of the averaged s fcwt measurements of each earthquake record is plotted against magnitude. Then, the best line is fitted to the measurements using the least-squares method. Atefi et al. [3] introduced a wavelet-based proxy, log-averaged scale klog, using weighted average power spectrum of different scales to obtain empirical magnitude-scaling relationship. In this study, we use frequency values which have maximum power spectrum. To better cooperation, we use the same database used in Atefi et al. [3]. Figure 2 presents the results obtained on the basis of s max P , s c , and s f cwt in the 4-s time window. The figure shows the logarithm of the s max P , s c , and s f cwt measurements against magnitude (gray circles). Obtaining a suitable empirical magnitude-scaling relationship necessitates the use of the average frequency parameters related to records of each earthquake. Correspondingly, the average values of the frequency parameters are plotted using black circles in Fig. 2. The best regression line was fitted to the logarithm of the event-average values of s max P , s c , and s f cwt and earthquake magnitude through the least-squares method. The constant coefficients of a and b associated with the best fitted line (Eq. 7) are listed in Table 1. The mean and standard deviation of the difference between observed and estimated magnitude jM est À M obs j residuals were used to verify the accuracy of the proposed s f cwt relationship. Figure 3 and Table 1 present the residuals against the observed magnitudes of earthquakes. The results showed that the mean residuals were 0.45 ± 0.43, 0.49 ± 0.36, and 0.66 ± 0.59 magnitude units in the 4-s time window, as determined by s max P , s c , and s f cwt , respectively. To further evaluate the accuracy of the results for earthquakes with M\ 5.5, the mean and standard deviation of the residuals were separately derived ( Table 1). The s max pbased method was more precise for earthquakes with M\ 5. The s f cwt -based approach exhibited more acceptable precision in determining the magnitudes of mediumand small-scale earthquakes than did the conventional s cbased method.

Conclusion and discussion
Earthquake early warning systems use the frequency information of ground motion within the first few seconds of P-wave arrival in estimating earthquake magnitude. To this end, s max P -and s c -based methods were proposed. However, these parameters were formulated on the basis of assumptions that are acceptable for monochromatic and noise-free signals. In other words, the parameters tend to be

www.SID.ir A r h i v e o f S I D
biased toward high frequencies. Therefore, magnitude estimations, particularly for areas with high noise levels, are associated with unavoidable errors. To improve the accuracy of magnitude estimation on EEWSs, we developed a new frequency-based parameter called s f cwt using CWT as a tool for extracting the frequency contents carried by the first few seconds of P-wave arrival. For comparison, s max P and s c parameters were also measured using the same database for s f cwt . In general, three empirical magnitudescaling relationships were extracted using s max P , s c , and s f cwt at a 4-s time window. These relationships were as follows: To compare the precision of relations (10)- (12), the average residuals for earthquakes with M JMA \ 5.5 were calculated ( Table 1). The s max P parameter appeared to be more sensitive to magnitude ranges of earthquakes than were s fcwt and s c (Fig. 3, Table 1). The values calculated using s f cwt within 4 s were plotted against the calculated values at time windows of 2 and 3 s for comparison (Fig. 4). As expected, the values were close to one another for small values of s f cwt (frequency content of small

A r h i v e o f S I D earthquakes).
For larger values of s f cwt (frequency content of large earthquakes), the calculated values clearly differed (particularly those calculated within 2 and 4 s). We can conclude that s f cwt was positively correlated with magnitude. Therefore, given that the magnitude ranges of the earthquakes studied reached 8.0, relations for estimating the magnitude obtained at 4 s were proposed (relations 10-12). The s max P -and s c -based methods calculate the frequency content carried by a signal using two components of records, whereas the s f cwt -based approach uses only a single component of acceleration. Nevertheless, processing is more complex and time-consuming with the use of s f cwt (because of its inherent CWT). This challenge can be overcome by upgrading hardware and using powerful processing systems. In general, the precision of this parameter can be evaluated through early warning systems across the world.

Data and resources
The used records are obtained from KiK-net on-line databases. Figures are prepared using Generic Mapping Tools ( [22]; last accessed August 2016). Data processing is performed using Seismic Analysis Code (https://ds.iris.edu/ds/ nodes/dmc/software/downloads/ sac, last accessed May 2014).