Multi-parameter algorithm for Earthquake Early Warning

Abstract

Earthquake Early Warning (EEW) system is a real-time earthquake damage mitigation system. It detects, analyzes and transmits information of the impending event at the potential user sites. In the present study, an attempt has been made to develop a multi-parameter-based EEW algorithm for accurate and reliable issuance of EEW. In the present study, not only the individual parameters such as maximum predominant period (), average period (τ_c)_, peak displacement (P_d), cumulative absolute velocity and root sum of squares cumulative velocity have been used to develop an algorithm but also various combinations have been attempted to issue alarm. The criterion for issuing warning is based on the alarm status of combinations of parameters at nearest four stations within selected epicentral distance of the event. It has been concluded that three parameter preference-based approach provides most efficient results at a time window of 4 s. The data-set for estimating various parameters required for multi-parameter algorithm has been taken from K-NET seismic array (Japan). The algorithm has been validated on the data taken from different regions. Based on the validation carried out using Indian and PEER-NGA database, the results revealed better performance of multi-parametric EEW algorithm.

1. Introduction

The recent spurt in seismic activity may not at all be attributed to denser instrumentation with high resolutions increasing the number of events in smaller magnitude seismicity but has been conspicuous by experiencing the damaging earthquakes in last few decades all over the world. This has the support of many of the seismic hazard studies carried out in the last few decades for various regions (Sharma & Lindolhm 2012; Paul & Sharma 2011; Ameer et al. 2005, etc.) Such increase in damaging earthquakes has synergistic effect while combining the development in seismic instrumentation in the recent past resulting in systems like Earthquake Early Warning (EEW) where its importance has been understood by all the seismically active countries to minimize the loss caused by such catastrophic events. Further, the rapid advancement in the field of real-time transmission and processing of seismological data has led to the development of real-time risk mitigation system like EEW system. The available few seconds time before the impending damaging event, allowed by such early warning systems, when wisely used reduced seismic risk by reducing loss of life, property and other resources. The EEW systems utilize the difference in the speed of non-destructive primary waves (P-wave), destructive secondary waves (S-wave) and high-speed transmission of data (electromagnetic waves) for issuing alarm. Depending upon the distance between the earthquake source and the target site, lead time (i.e., the time between the P-wave detection by a seismic network and the S-wave arrival at the target site) is estimated.

The concept first put forth by Dr Cooper around 145 years ago has now shaped into a modern EEW system with state-of-the-art technology being used by many of the countries around the world. A number of earthquake-prone countries are engaged in ongoing efforts to develop EEW systems having fast on accurate real-time operating EEW algorithm. Japan (Nakamura 1988; Odaka et al. 2003; Kamigachi 2009), Mexico (Espinosa-Aranda et al. 1995, 2009), Taiwan (Wu et al. 1999; Wu & Teng 2002, 2007), Romania (Wenzel et al. 1999; Böse et al. 2007) and Turkey (Erdik et al. 2003; Böse et al. 2008; Alcik et al. 2009) have developed EEW system, which proves effective in mitigating seismic risk at short time scales (Allen et al. 2009). Other countries like California (Allen & Kanamori 2003; Allen et al. 2009; Wurman et al. 2007; Cua & Heaton 2007), Switzerland (Cua & Heaton 2007) and Italy (Satriano et al. 2011) are in early stage of EEW system implementation.

High-density seismic networks with real-time data transmission capability are considered to be the primary requisite for an EEW system which is followed by an algorithm for the automatic processing of seismic data, and to retrieve the seismic source parameter information from the processed seismic signal. Thus, the algorithm encompasses the detection of earthquake and determines seismic source parameters to estimate magnitude of the earthquake for issuing warning.

The issuance of warning by EEW system depends upon the accuracy and reliability of predicted parameters used to define the size of the incoming event in real time. Various parameters are used for estimation of magnitude and issuing alarm in EEW system worldwide. In the present study, an attempt has been made to combine various EEW parameters to develop an algorithm, which can issue alarm and estimate magnitude with reliable accuracy in minimal time window. A multi-parameter-based EEW algorithm has been developed and demonstrated for more accurate and reliable issuance of EEW.

2. Data-set

The EEW systems are based on the parameters relating the size of the impending event. Such parameters are estimated empirically for their threshold values, thus, increasing the requirement of a good data-set to work upon. In the present study, the data-set has been carefully selected from seismically active regions having good coverage of instrumentation. The basic requirements for the data have been classified into three subsets namely, (1) data required for regression to estimate the parameters for EEW, (2) validation of data in another seismically active region and (3) to check the stability of parameters using worldwide data.

The first data-set has been chosen from K-NET seismic array in Japan. Out of the whole data, 1726 records of 105 events having 5 ≤ M ≤ 7.2 with epicentral distance ≤ 60 km have been selected. To validate the developed algorithm in seismically active region, the Himalayas have been chosen as the study region. About 300 stations of strong motion accelerographs (with inter-station spacing of about 40 km) have been installed in north and north-eastern India by Department of Earthquake Engineering (Kumar et al. 2012) which has recorded a good amount of strong motion data from Himalayan region. All 28 events have been selected from these seismically active regions having 51 digital records from stations having epicentral distance up to 60 km from a range of magnitude varying between 3.3 and 6.8. To check stability of the algorithm, another data-set (not used for regression analysis in the study) has been selected mainly from: Southern California, Taiwan and Turkey. This data-set consists of 219 earthquake records of 14 earthquakes recorded at 174 strong motion stations with magnitude ranging 4.27 ≤ M ≤ 7.62 within 60 km of epicentral distance. The distribution of selected events magnitude with number of records for K-NET, Indian Himalaya and PEER-NGA data-set has been shown in Figure 1. Figures 2, 3 and 4 show the location of selected K-NET, India, southern California, Turkey and Taiwan events and their corresponding recording stations within desired epicentral range of 60 km.

Figure 1. Magnitude and corresponding number of earthquake records used from (a) Japan, (b) Indian and (c) PEER-NGA strong motion data-set.

Figure 2. Map of Japan showing location of earthquakes and stations used in the present study.

Figure 3. Northern Indian region showing location of earthquakes and stations used in the present study with triangles representing the strong motion stations, stars representing the epicenters of the earthquakes and ellipse representing the regions (1, 2, 3 and 4) from where the strong motion data has been selected for validation.

Figure 4. Location of earthquakes and recording stations in (a) Southern California, (b) Turkey and (c) Taiwan. Solid triangles representing the strong motion stations and stars representing the epicenters of the earthquakes.

short-legendFigure 4.

In the present study, the vertical components of the acceleration time series having sampling rate of 100 samples per second have been used. The data having higher sampling rate has been re-sampled using developed MATLAB code. The velocity time series have been obtained by integrating the acceleration strong motion records, and then filtered with a high-pass Butterworth filter having a cut-off frequency of 0.075 Hz to remove the low-frequency drift introduced during the integration process.

3. Methodology

The EEW is generally classified into two approaches called regional warning approach and onsite warning approach. Both types of EEW approaches are designed based on different EEW parameters such as maximum predominant period ( ), average period (τ_c)_, peak displacement (P_d), cumulative absolute velocity (CAV), bracketed cumulative absolute velocity (BCAV), windowed bracketed cumulative average velocity (BCAV-W) and root sum of squares cumulative velocity (RSSCV). is considered to be the first EEW parameter, which utilizes first few seconds of P-wave to calculate the maximum possible predominant period within a selected time window (Nakamura 1984, 1988). Kanamori (2005a) proposed τ_c similar to . While is calculated from the ratio of velocity and acceleration records, τ_c was proposed from ratio of velocity and displacement records. Further, determines predominant period of the P-wave, while τ_c represents initial portion of P-wave over a fixed time window, i.e., the pulse width of P-wave and it varies according to the size of event, thus used for estimating the event magnitude (Wu & Kanamori 2008a). Both the and τ_c parameters deal with the frequency content present in the initial portion of P-wave and when used together to estimate the magnitude, more accurate results were obtained (Shieh et al. 2008). The τ_c parameter was used by Wu and Kanamori (2005a, 2005b) for the development of onsite EEW system in Taiwan.

Wu and Kanamori (2005a, 2005b) and Wu et al. (2006, 2007) projected another EEW parameter which is based on the maximum amplitude of the displacement time history calculated from the vertical strong ground motion component in a time window of 3 s after P-onset. P_d parameter is used to estimate the peak ground velocity, i.e., the damageability of the impending earthquake at a site. Kennedy and Reed (EPRI 1988) introduced CAV parameter. The relation between the CAV and the local intensity was noticed and it is used as a parameter to distinguish a potentially damaging event from a non-damaging one at nuclear power plant sites regardless of the response spectral values (EPRI 1988). EPRI (1991) modified the CAV calculation by removing the non-damaging accelerations values from the records. The time histories having the maximum acceleration values greater than the predetermined acceleration value (0.025 g) within a selected time domain were considered for calculating modified CAV values and are termed as standardized cumulative absolute velocity (BCAV). Alcik et al. (2009) modified BCAV parameter as BCAV-W. The calculation for BCAV-W was performed by windowing BCAV calculation on a broader window length (W) basis. BCAV-W parameter along with PGA parameter is used in Istanbul EEW system with more accurate alarms and less false alarms. Later, a new EEW parameter called RSSCV (Bhardwaj et al. 2013a) has been considered, in which the energy content in the incoming time series is assumed to be present in terms of velocity components in relation to other component of ground motion. Therefore, not only the predominant velocity in the time series is considered but also the cumulative effect is included in terms of its energy component.

The above-said parameters have been used individually or in combinations by various groups for EEW systems. However, it has been found that each parameter has different correct alarm capability than other. Such capability depends on the parameter reflecting and reacting to specific property of the ground motion defining the size of the event in the form of the energy present at that moment. Therefore, in the present study, an attempt has been made to develop a multi-parameter-based EEW algorithm which provides more accurate and reliable EEW. The mathematical formulation of the EEW parameter considered in the designing of multi-parameter EEW algorithm is listed in Table 1.

Table 1. EEW parameters considered for the development of multi-parameters-based EEW algorithm.

EEW parameters	Numerical formula
(Nakamura 1984)
τc (Wu & Kanamori 2005a, 2005b)
Pd (Wu & Zhao 2006)	Maximum value of high-pass filtered vertical displacement amplitude
CAV (EPRI, 1988)
RSSCV (Bhardwaj et al. 2013a)

In EEW system, the first step is automatic detection of P-onset in a continuously streaming data. Allen's (1978) automatic P-onset detection algorithm is widely used in EEW systems, which are based on short-term average and long-term average ratio. However, RSSCV parameter (Bhardwaj et al. 2013a) is preferred for automatic P-onset picking in real-time RSSCV parameter as it has its own merits of providing estimate of the energy in the signal which can be used for magnitude estimation of the event too. The present algorithm has included the module for location of the events and presently uses the locations of the earthquakes known apriori. The parameters estimated using the initial portion of incoming time series have been analyzed for retrieving the information about earthquake source parameters which are empirically regressed with the magnitude of the event. Five parameters namely , τ_c, P_d, CAV and RSSCV have been used in the present study for which the threshold values are estimated using the recorded data.

To obtain empirical relationships between EEW parameter and magnitude, , τ_c and P_d have been calculated at different time windows (1–5 s). The threshold magnitude (considered as 6.0 in the present study) is then used to define the thresholds for these parameters to issue warning in future. However, in case of CAV and RSSCV parameters, the threshold values for issuing warning have been determined iteratively and manually using K-NET data-set that comprises of 1726 records of 105 earthquakes. All the five EEW parameters , τ_c, P_d, CAV and RSSCV are analyzed at different time windows (1–5 s) and specified time domain amplitude levels have been estimated for issuing warning for M ≥ 6. Estimation of impending earthquake magnitude helps in determining the amount of damage it can cause to target site which is utmost important for issuing warning. Simultaneously, with each addition of second in processing time, location and intensity of the event are also estimated for issuing more accurate and reliable warning by an EEW system. Since the initially computed magnitudes are not always updated when more and more data is recorded, therefore the magnitude information is given qualitatively only. For issuance of warning, various magnitude estimation approaches have been used for EEW around the world. Japan uses combination of energy and growth rate measure for magnitude estimation, while in Mexico maximum displacement amplitude has been used for magnitude estimation (Kamigaichi et al. 2009). In the present study, Brune's model has been used to estimate the magnitude in real time (Bhardwaj et al. 2013b).

Figure 5 shows the steps involved in designing of EEW algorithm starting from pre-processing of vertical component of time series such as baseline correction, and filtering followed by P-onset picking. After P-onset detection, EEW parameters are calculated at four nearest stations for a selected time window. The calculated EEW parameters are then compared with their preset threshold values. If the calculated parameter has value higher than the corresponding threshold value, an alarm status is set. The alarm is issued for the event only when a defined number of stations give alarm status. Meanwhile, Brune's model-based magnitude estimation approach (Bhardwaj et al. 2013b) is simultaneously estimating the magnitude of the event at all the triggered stations at a selected time window after P-onset. The average of event magnitude estimated at the triggered station gives the final magnitude of the event. It may be noted that the magnitude estimated has been used as information only, and it is not participating in issuance of alarm.

Figure 5. Steps involved in designing of EEW algorithm. P1, P2, P3 and P4 represent considered EEW parameters, THR represents the preset threshold values and M represents minimum number of alarm required to issue warning along with estimated magnitude for the event.

The algorithm is designed in a way such that it may issue alarms based on individual parameter and/or a combination of parameters. The parameters used in EEW system reflect the basic properties of the incoming time series in terms of their amplitude with respect to real time, ground motion in terms of acceleration, velocity and displacement. While increase in amplitude with respect to time is used to mark the P-onset, the other properties are used for estimating the size, incoming energy and the damageability. Based on the type of source (thrust/strike-slip faults), azimuthal direction (radiation pattern), focal depth (shallow/deep), distance (local/regional) and the magnitude (micro-earthquakes/moderate earthquakes/strong earthquakes), these parameters react differently. Therefore, based on the prevailing seismotectonics and the physical sources present in the region, a region specific combination of EEW parameters is required to be worked out for an efficient EEW system. The next section contains an exercise to find the appropriate combination of these parameters for realistic issuance of EEW.

4. Analysis and results

In the present study, K-NET strong motion data-set has been used for the analysis due to its high spatial density of the network in the source region. The analysis carried out at five windows (1, 2, 3, 4 and 5 s) to search for the accurate and minimal time window is shown in Figure 6 where the observed , τ_c and P_d values at each station have been plotted against the catalogue magnitude at each event. The least square fit to the average values of , τ_c and P_d with magnitude have been shown by different lines. Figure 6 reveals an increase in slope and improvement in magnitude estimate with each additional second of data stabilizing at 4 s. The scattering effect of at lower magnitudes gets reduced with increase in window length which verifies Brown et al. (2009) observation. In case of τ_c at 1 s, data results in smaller slope to the scaling relation leading to systematic underestimation of magnitude for large events due to saturation at higher magnitude (Kanamori 2005), whereas an increase in slope and improvement in magnitude estimate for larger events has been observed at higher time windows. In case of P_d with the increase in time window length, same improvement in magnitude estimation has been observed. In the present study, separate scaling relations for each time window has been calculated. The least square fit to the data at different time windows produces scaling relations as shown in Table 2. The CAV and RSSCV values have been calculated from 1726 records of 51 earthquakes using equations (1 and 2).

Figure 6. Effect of number of seconds of P-wave data on scaling relations between (a) , (b) τ_c and (c) P_d with magnitude. Crosses represent observation at individual stations and different types of lines are linear best fit scaling relations to crosses at different time windows. Combined effect of all the five time windows (1, 2, 3, 4 and 5 s) shows the rise in slope with each additional second of data.

Table 2. Scaling relation obtained between , τ_c and P_d with magnitude at different window lengths.

Time window (s)	Obtained scaling between EEW parameter and magnitude	Conversely, in terms of magnitude
1
2
3
4
5

CAV is computed by integrating acceleration record a(t) to time duration (t_max) of record (1)

The RSSCV is an EEW parameter which includes the cumulative effect (amplitude and time) of ground motion duration. RSSCV is also an integral EEW parameter same as CAV (Bhardwaj et al. 2013a), and is defined as (2) where v_i is the velocity vector calculated by taking integral of strong motion records for n number of samples in the selected window.

The threshold values for five time windows have been calculated by manually observing the best fit on the records as shown in Figures 7 and 8. The individual station records are marked by cross and average of event is marked with plus. The solid horizontal line represents the specified threshold value for the particular time window best obtained after analysis. In case of CAV and RSSCV values, scatter gets reduced with increasing time window with a rise in correct alarm rate. It has been found that with each additional second, scattering in the data gets reduced and better estimation of magnitude for large events can be achieved. This verifies the finding of Allen and Kanamori (2003) that the initial magnitude estimate using 1 s data is a minimum estimate and it gets improved with each additional second of data, but with a reduction in warning time.

Figure 7. Specified threshold values of CAV for all the five time windows (1, 2, 3, 4 and 5 s). 1726 records are marked with cross and plus represents the average CAV value of 105 earthquake events. Solid horizontal line represents best fit threshold value and dashed vertical line is marked at M = 6. The CAV values scatter gets reduce with increasing time window with a rise in correct alarm rate.

Figure 8. Specified threshold values of RSSCV for all the five time windows (1, 2, 3, 4 and 5 s). 1726 records are marked with cross and plus represents the average RSSCV values of 105 earthquake events. Solid horizontal line represents best fit threshold value and dashed vertical line is marked at M = 6. The RSSCV values scatter gets reduce with increasing time window with a rise in correct alarm rate.

The threshold values for CAV and RSSCV have been calculated for five time windows by observing the best fit on the records. CAV describes the cumulative effect of generated ground signal in the selected time, thus used to differentiate between a damaging and non-damaging event. RSSCV is another EEW parameter which includes the cumulative effect (amplitude and time) of ground motion duration (Bhardwaj et al. 2013a). RSSCV helps in enhancing the SNR and reduces the standard error in observational seismology. In the present study, RSSCV proved to be an important EEW parameter as it has been used for automatic P-pick as well. A comparison of all the estimated threshold values for the five EEW parameters at different time windows have been listed in Table 3 which reveals that the amplitude of all the threshold values increases with the increase in time window.

Table 3. Calculated threshold values of different EEW parameter for issuing warning for an event having M ≥ 6 at different time windows.

	Time window
Parameters	1 s	2 s	3 s	4 s	5 s
(s)	0.95	1.00	1.06	1.10	1.14
τc (s)	1.02	1.17	1.20	1.42	1.55
Pd (cm)	0.13	0.27	0.51	0.95	1.38
CAV (cm/s)	3.00	8.00	10.00	23.00	41.00
RSSCV (cm/s)	0.30	1.00	1.70	5.20	10.00

The algorithm considers four nearest stations which get triggered on the onset of the event within 60 km epicentral range. If the calculated EEW parameter value at the considered station exceeds the preset threshold value, it is considered as a vote for the event. If out of four nearest stations, three stations vote within the selected time window, then the EEW algorithm issues alarm for that event. For the check in the efficacy of the EEW parameters and the developed algorithm, the alarms are classified in four categories as: (a) correct alarm (CA) means issuing alarms where alarm needs to be issued and have been issued; (b) false alarm (FA) means no alarm needs to be issued but issued; (c) correct all clear (CAC) means no alarm needs to be issued and also not issued; and (d) missed alarm (MA) means issuing alarm where alarm need to be issued but not issued. The concept is described in Figure 9 where cross represents the calculated EEW parameter values; for example, RSSCV value calculated at 5 s time window and solid horizontal line represents the specified threshold value. The CA, FA, CAC and MA zones are represented by Quarter-I, Quarter-II, Quarter-III and Quarter-IV, respectively, in Figure 4. Further, CA and CAC are together termed as the correct detection (CD) while the FA and MA are the undesired alarms and together called as incorrect alarm (ICA). The efficiency of individual EEW parameter has been calculated in terms of number CD and ICA. To provide reliable and accurate warning results in minimal time length, the search for the best possible combination of EEW parameters has been carried out.

Figure 9. Four different alarm zones as CA, FA, CAC and MA with cross represents the EEW parameter (RSSCV) value calculated at individual stations. The horizontal line represents the specified threshold level and the vertical line represents the line at M = 6.

The EEW parameters give different efficiencies based on the number of CA, MA, CAC and FA. The more will be the number of CA and CAC, the better will be the efficiency of the algorithm. MA and FA are the undesired alarms and reduce the efficiency. Thus, in a reliable and accurate EEW algorithm, it is required that CD number should be more than ICA. In the selected data-set, there are 105 earthquakes out of which there are 24 earthquakes having M ≥ 6 for which warning need to be issued. There are 81 earthquakes having M ≤ 5.9 for which warning need not to be issued because for the present study warning threshold has been set at M ≥ 6. The number of CA and MA are calculated out of 24 earthquakes and CAC and FA are calculated out of 81 earthquakes while, CD and ICA values are calculated from total 105 earthquakes. Time taken by the algorithm for the accurate and correct warning is also considered as an important point of concern in an EEW algorithm. Thus, the efficiency of the algorithm has been checked at five different time windows (1– 5 s) to find out the appropriate time window for issuing reliable warning.

The analysis of all the five EEW parameters (, τ_c, P_d, CAV and RSSCV) considered in the present study have their own rate of CA and ICA as discussed above at different time windows as shown in Table 4. In starting 2 s, the results obtained by the EEW parameters are not consistent and reliable, but at 3 s and above, the results gets stabilized. It has been observed that at 4 s and above, the performance level of all the five EEW parameters attains its reliable and fairly accurate state. Further, in the present study, the next step is to search the best possible combination of EEW parameter that gives most accurate and reliable results in all the five time windows.

Table 4. Number of CD and ICA achieved out of 105 events using considered EEW parameters.

			τc		Pd		CAV		RSSCV
Time window (s)	CD	ICA	CD	ICA	CD	ICA	CD	ICA	CD	ICA
1	98	7	100	5	90	15	75	30	52	53
2	79	26	98	7	93	12	78	27	53	52
3	97	8	98	7	94	11	76	29	67	38
4	98	6	99	6	102	3	89	16	95	10
5	100	5	103	2	103	2	98	7	97	8

5. Search for the best EEW parameters combination

Different possible combinations of considered EEW parameters have been tested for the designing of accurate and fast warning EEW algorithm. The Four possible combinations which have been analyzed in the present study for the search of best EEW algorithm are:

• (C1) Three EEW parameters preference-based approach – if for an event three out of five EEW parameters have similar alarm status, the event has been marked with the same status.

• (C2) Four EEW parameters preference-based approach – if for an event four out of five EEW parameters have similar alarm status, the event has been marked with the same status.

• (C3) Five EEW parameters preference-based approach – if for an event all the five EEW parameters have similar alarm status, the event has been marked with the same status.

• (L1) Logic combination of EEW parameters preference-based approach – comprised of the combination of similar types of EEW parameters and is logically represented as: (3)

The sets in L1 are selected based on the physical properties being represented by the EEW parameter in its own sub-set. In Table 5, the number of CD and ICA obtained using the selected EEW parameter preference-based approaches have been compared. The number of CD in case of L1-based approach is better than C3-based approach, comparable with C2-based approach and poorer than C1-based approach. Thus, it has been concluded that on comparing the number of CD obtained by all the possible four approaches, it has been found that at a time window of 4 s, C1 approach provides the efficient, accurate and reliable results.

Table 5. Comparison in number of CD and ICA achieved out of 105 events using various EEW parameter preferences-based approaches.

	EEW parameter preference based-approaches
	C1		C2		C3		L1
Time window (s)	CD	ICA	CD	ICA	CD	ICA	CD	ICA
1	99	6	69	36	38	67	69	36
2	96	9	71	34	42	63	70	35
3	101	4	72	33	54	51	71	34
4	104	1	93	12	75	30	91	14
5	104	1	98	7	88	17	96	9

Estimation of magnitude to issue warning for an event within a few seconds after it originates is integral task of an EEW algorithm. After issuing the warning for M ≥ 6, the algorithm further explores the data for more accurate magnitude determination. Brune's model has been used to fit on initial P-wave data after P-onset, at different time windows starting from 1 to 10 s for estimating spectral parameters (Bhardwaj et al. 2013b). Bhardwaj et al. (2013b) uses variable length of time windows which results in a diversity of measurement of earthquake size. Best results with minimum window length that can be used for magnitude estimation is found to be 5 s with an uncertainty of ±0.3 units. After P-onset, a time window of 5 s may be selected for calculating the spectral parameters (Ω_o, f_c and f_max) by fitting Brune's model in the selected source spectra. These parameters are then used for estimating the magnitude of the earthquake. The estimated magnitude is found to be in good agreement with the catalogue magnitude with a minimum uncertainty which is due to the inherent error in M_w estimation (Das et al. 2011; Das 2013).

The analyses of K-NET data-set revealed that most reliable performance of the EEW algorithm is obtained by using three parameter preference-based approach at 4 s. Further, in case of magnitude estimation at 5 s, good magnitude estimation with small deviation from catalogue magnitude has been obtained (Bhardwaj et al. 2013b). Thus, obtained final EEW algorithm is shown in Figure 10. The proposed EEW algorithm starts initially with vertical component of time series pre-processing which includes baseline correction by fitting polynomial of first order and filtering of data by using fifth-order high-pass Butterworth filter having cut-off frequency of 0.075 Hz. The P-onset is marked and EEW parameters: , τ_p, P_d, CAV and RSSCV are calculated at four nearest stations for a time window of 4 s after P-onset. The calculated EEW parameters are then compared with their preset threshold values such as 1.1 s, 1.42 s, 0.95 cm, 23 cm/s and 5.2 cm/s, respectively. If the calculated parameter has value higher than the corresponding threshold value, an alarm status is set. The alarm for the event is issued only if three parameters out of five have alarm status.

Figure 10. Proposed EEW algorithm, using three parameter preference-based approach at 4 s time window for issuance of warning. , τ_c, P_d, CAV and RSSCV parameters are used in the warning analysis. Brune's model-based magnitude estimation approach is used for simultaneously estimating the magnitude of the event within a time window of 5 s.

Brune's model-based magnitude estimation approach is used with selected time window of 5 s for simultaneously estimating the magnitude of the event at all the triggered stations. After pre-processing, cosine tapered window has been applied on the baseline corrected time series. Further, the windowed time series has been Fourier transformed to obtain velocity and displacement time series from acceleration time series. The source displacement and acceleration spectra are analyzed for calculating three spectral parameters namely Ω_o, f_c and f_max within 5 s window to be used for calculating source parameter viz., M_o and in turn, M_w. The average of magnitudes estimated at various triggered stations give the final magnitude of the event. It may be noted that the magnitude estimated has been used as information only and it is not participating in issuance of alarm.

6. Validation using Indian and PEER-NGA data-set

The developed EEW algorithm has been tested using Indian and PEER-NGA database for validating the proposed EEW algorithm. Indian and PEER-NGA database have been tested at specified threshold values and combination of parameters for validation of the developed algorithm. Indian strong motion data-set comprised of 51 records of 28 earthquakes having magnitude range varying from 3.3 to 6. Out of 28 events, there were 27 events having magnitude less than 6 for which no warning need to be issued. However, there has been a single event of magnitude 6.8 for which warning need to be issued. At an epicentral distance of 60 km, Indian strong motion data-set have small number of records available.

On applying the proposed EEW algorithm for five different time windows, the number of alarms obtained is shown in Table 6. It has been observed that in case of Indian data-set, three EEW parameter preference-based approaches have found to be quite effective from 2 s onwards. The algorithm has also been explained in terms of numbers of CD and ICA in Table 6. At 1 s, 27 events have been correctly detected with only one ICA. The results showed that the Indian data-set is in strong agreement with proposed three parameter preference-based approach at 4 s.

Table 6. Number of CD and ICA achieved out of 28 events using C1-based approach.

Time window (s)	CD	ICA
1	27	1
2	28	0
3	28	0
4	28	0
5	28	0

The PEER-NGA data-set which consist of strong ground motion records of various countries for e.g., Southern California, Turkey and Taiwan. There were total of 219 records of 14 earthquakes having magnitude range of 4.27–7.62 up to 60 km epicentral distance. In 14 events, there are equal numbers of events, i.e., 7 each, for which warning need to be issued (M ≥ 6) and not to be issued (M < 6).

The algorithm has been explained in terms of number of CD and ICA in Table 7. At 1 s, there are 9 events which are correctly detected, while 5 events are wrongly interpreted. At 2 s, the rate of CD increases from 9 to 10 and one point drop in wrong estimation is observed. At 3 and 4 s, similar numbers of CA and ICA alarms are obtained. However, at 5 s, the three parameter preference based algorithm achieves its maximum efficiency with 14 CA and none ICA. At calculated specification of three parameter preference based approach and at time window of 4 s, the algorithm got noticeable efficiency and thus, validation using PEER-NGA database has been achieved.

Table 7. Number of CD and ICA achieved out of 14 events using C1-based approach.

Time window (s)	CD	ICA
1	9	5
2	10	4
3	13	1
4	13	1
5	14	0

Due to advancement in computers, the time taken is relatively very small in computations than in transmission and issuance of warning. The algorithm is designed to issue alarm at a place where the data from the participating station is stream i.e., warning and distribution centre and alarm is issued along with the information of lead time and magnitude. The developed EEW algorithm takes 4 sec to compute EEW parameters and to make decision about the alarm with 1 s of transmission delay and 1 s of processing delay.

7. Conclusions

In seismically active regions encompassing areas either developed or developing, EEW system has been considered as one of the disaster-mitigating strategies. In this context, it is necessary to exploit the incoming time series for EEW parameters and use them appropriately recognizing their physical significance in relation to the size of the impending earthquake event. An algorithm based on multi-parameter approach has been developed in the present study. After the P-onset detection, thresholds for issuing warning have been estimated for individual parameters namely , τ_c, P_d, CAV and RSSCV at different time windows by regressing these parameters with apriori known catalogue magnitude using a data-set from K-NET, Himalaya (India), southern California, Taiwan and Turkey. The warning has been pre-meditated based on crossing of threshold levels for an event at least at three stations out of the four considered nearest stations. Based on the accurate, reliable and fast warning estimation results, the best possible combination of EEW parameter has been selected i.e., three parameter preference-based approach at a time window of 4 s. For the purpose, the warning has been classified into four different alarm types viz., CA, missed alarm (MA), CAC and FA which is further classified in two types: CD and ICA.

The parameters considered in the present algorithm are designed to capture different attributes related to time domain, frequency domain, amplitude domain and power domain. Depending on the properties and the physical characteristics of the source mechanics and the medium involved in the wave propagation, the damageability, i.e., the size of the event, may be reflected in any of the four domains represented as five parameters in EEW system. The endeavour in such cases is always to capture the early birds and issue alarm as soon as possible. Therefore, in the present study, any three out of five parameters were tested and found to give relatively better results within stimulated time of 4 sec. It is possible that in the present case P_d has behaved relatively better than other but this may not be the case when the algorithm from different settings of the geological and tectonic environment of the region. Thus, to be more stable and accurate in issuance of alarm in the present study, we have selected more than one and two parameters.

Based on the validation carried out using Indian and PEER-NGA database, the results revealed better performance of multi parametric EEW algorithm in comparison to classical use of individual parameters in EEW systems.

Disclosure statement

No potential conflict of interest was reported by the authors.