Yumeng Hu ,Qi Zhng ,Wenji Zho ,Hito Wng

a Lanzhou Institute of Seismology,China Earthquake Administration,Lanzhou,730030,China

b University of Science and Technology of China,Hefei,230026,China

c Institute of Geology,China Earthquake Administration,Beijing,100029,China

d China Earthquake Networks Center,Beijing,100045,China

Keywords:Transformer Deep learning Seismic P-Wave detection

ABSTRACT Recent years,we have witnessed the increasing research interest in developing machine learning,especially deep learning which provides approaches for enhancing the performance of microearthquake detection.While considerable research efforts have been made in this direction,most of the state-of-the-art solutions are based on Convolutional Neural Network (CNN) structure,due to its remarkable capability of modeling local and static features.Indeed,the globally dynamic characteristics contained within time series data (i.e.,seismic waves),which cannot be fully captured by CNN-based models,have been largely ignored in previous studies.In this paper,we propose a novel deep learning approach,TransQuake,for seismic P-wave detection.The approach is based on the most advanced sequential model,namely Transformer.To be specific,TransQuake can exploit the STA/LTA algorithm for adapting the three-component structure of seismic waves as input,and take advantage of the multihead attention mechanism for conducting explainable model learning.Extensive evaluations of the aftershocks following the 2008 Wenchuan MW 7.9 earthquake clearly demonstrates that TransQuake is able to achieve the best detection performance which excels the results obtained using other baselines.Meanwhile,experimental results also validate the interpretability of the results obtained by TransQuake,such as the attention distribution of seismic waves in different positions,and the analysis of the optimal relationship between coda wave and Pwave for noise identification.

1.Introduction

With the rapid development of seismic monitoring technology,more and more attention has been paid to the efficient detection and differentiation of microearthquakes from massive noise data,such as the intensive aftershock sequences following a destructive earthquake.Indeed,a variety of methods,such as template matching (Frank et al.,2014;Gibbons et al.,2007),similarity searching(Yoon et al.,2015),and higher-order statistics (Saragiotis et al.,2000),have been developed in previous studies for addressing this challenge.Nevertheless,the most widely used technologies for seismic phase detection are developed based on the STA/LTA algorithm and/or its variants,which identify earthquake events when the ratio between the Short-Term Average(STA)and Long-Term Average (LTA) of energy function for seismic waves exceeds the pre-defined threshold (Allen,1982;Guo et al.,2019;Withers et al.,1998).However,the identification performance of all these traditional approaches is usually limited,either due to the trade-off effect between false and missing alarms or due to the that between the computational cost and time sensitivity(Liu and Zhang,2014;Reichstein et al.,2019).

Fig.1.The geographical location of the study region.

Recently,an increasing number of researchers have focusedon developing machine learning,especially deep learning,which provides approaches for enhancing the performance of microearthquake detection.While considerable research efforts have been made in this domain,most of the state-of-the-art solutions are based on Convolutional Neural Network (CNN) structure,due to its remarkable capability of modeling local and static features (Rouet-Leduc et al.,2020;Zhang et al.,2020).For example,Perol et al.(2018) propose an advanced CNN-based approach,CovNetquake,for earthquake detection and location.Based on CovNetquake (Zhu et al.,2019),propose a new approach CPIC for identifying and picking the seismic phase arrivals.However,since CNN is originally designed for pixel-based image processing rather than sequential data modeling,the global and dynamic characteristics of seismic waves cannot be fully captured by CNN-based models,which has been largely ignored in previous studies.

To address the above challenges,in this paper,we propose a novel deep learning approach,TransQuake,for seismic P-wave detection based on the most advanced sequential model in natural language processing(NLP),namely Transformer (Vaswani et al.,2017).To be specific,TransQuake can exploit the STA/LTA algorithm for adapting the three-component structure of seismic waves as input,and take advantage of the multi-head attention mechanism for providing explainable model learning.Extensive evaluations of the aftershocks following the 2008 WenchuanMW7.9 earthquake clearly demonstrates that TransQuake can achieve the best detection performance which excels the results obtained using other baselines.Meanwhile,experimental results also validate the interpretability of the results obtained by TransQuake,such as the attention distribution of seismic waves in different positions,and the analysis of the optimal relationship between coda wave and P-wave for noise identification.

2.Data and method

2.1.Data

The dataset used in this paper is the records of 100 Hz signals from the aftershocks following the May 12th,2008 WenchuanMW7.9 earthquake,which were recorded by 14 permanent earthquake stations (Fang et al.,2017) from July 1stto 31st,2008.The distribution of stations and the main event are shown in Fig.1.Each record of these stations contains 3 dimensions(3D):vertical component,north-south(N-S)component,and east-west (E-W) component.Followed by the settings in (Zhang et al.,2019),we filter the waves using the Bessel filter with bandwidth 2–10 Hz,and generate the negative samples by theFilterPicker(Lomax et al.,2012).Considering the different epicentral distances,we are unable to set a fixed time window that only contains a full P-wave,thus we set the time-window as 50 s to comprehensively analyze all the included waveforms.After data preprocessing,we have obtained 109 719 noise waveforms (negative samples),9 891 aftershock waveforms (positive samples),and the shape of each waveform is structured by 1×3×5 000.Fig.2 shows some examples of two types of waveforms.Note that,instead of directly generating the noise waveforms by background noises,here we.

Select negative samples from those misclassified ones as aftershocks after passing through the filter.Inspired by the Adversarial Machine Learning (Goodfellow et al.,2015),which is based on the Nash Equilibrium of Game Theory,we creatively generate these misclassified waveforms as adversarial samples.In this way,our model is able to capture the critical points of earthquake detection and subsequently extract vitally inherent features of seismic waveforms.Such process overcomes the weakness of traditional detection methods such as overfitting,and provides us a more robust model(Szegedy et al.,2014).

Fig.2.Some examples of labeled waveforms in the dataset.(a) Positive examples;(b) Negative examples.

Fig.3.The diagram of our problem formulation.

Fig.4.The architecture of our model,TransQuake.(a) Data preprocessing;(b) Fully connected and position code layer;(c) The full architecture of our model.

2.2.Problem formulation

Taking 3-D waveforms in a 50 s window as the original input,our model aims to learn and output the possibility of the results of each window:whether it contains an earthquake P-wave arrival or not.Fig.3 shows the diagram of our problem formulation,where U-D orientatesvertically,N–S and E-W orientate horizontally,standing for the northsouth channel and the east-west channel respectively.To formulate the problem,we useDto represent a training dataset which containsiwaveforms,and three components of each waveform correspond to our three channels.Formally,the training dataset is expressed as:

wherewirepresents a waveform and channelimeans a single signal.Eachwicorresponds to a labelli,whereli=｛0，1｝represents a noise whileli=｛1，0｝represents an earthquake event containing P-wave arrival.Finally,we can define the problem as follows:The objective is to learn a binaryclassification modelMfrom the training datasetD,which can classify a new waveform windowwithrough outputri,whereri=1 represents an earthquake event with P-wave andri=0 otherwise.

2.3.Model architecture

Fig.5.The architecture of the encoder in TransQuake.

Our model TransQuake is based on the most advanced sequential model Transformer in the NLP field.The original Transformer contains two parts,namely encoder and decoder (Vaswani et al.,2017).Since the decoder part is designed for generative tasks rather than discriminative tasks,we only use the encoder part to establish our model.The full architecture of our model is shown in Fig.4.Different from the original Transformer,to adapt the continuous waveform data rather than independent words as the input,we first propose a new approach to obtain the input of encoders (as shown in Fig.4(a) and (b)).Then,with the adapted encoders of the Transformer,we further exploit a linear layer and a softmax layer to obtain the classification results,i.e.,whether the input data contain a seismic P-wave or not.The related mathematical notations are summarized in Table 1.

Table 1 Some important mathematical notations.

Table 2 Definition of the confusion matrix for model evaluation.

2.3.1.Input

Learning from the earthquake real-time monitoring system of the National Digital Seismic Network,we extract the 150 time-steps STA and the 2000 time-steps LTA of the waveforms (Guo et al.,2019).Using STA/LTA,the data are able to generate better waveform features and model performance.It indicates that machine learning methods combining traditional physical models can complement each other and achieve better performance (Reichstein et al.,2019).Next,as the Transformer is proposed for NLP,its original input is a batch of words,while our input is 3-D continuous-time data.To adapt to this difference,we propose a module to make corresponding changes in input processing.In NLP tasks,each word is embedded into a vector whose size equals the hidden units.Then the position code (PC) is added.As shown in Fig.4,we reshape the waveform and follow a fully connected (FC)layer to explore the relationship between each position of vectors,and then add the PC to get the input of encoders.The processing details are as follows.Considering the combination of 3 components in each position and the suitable size of the attention model,we reshape each waveform from 3 × 5000 intok× 3T,wherek=（3 × 5000）/ 3T.Then we parallelly feed it to FC and PC layers.In the FC layer,each neuron cell connects all new cells as:

whereymis the value of a new cell withm∈｛0，1，…，d｝ andddenotes the dimension of hidden units;Wnis the wight that should be trained;xnis the value of the original cell.After the FC layer,the last dimension of data equals the hidden units.In the PC layer,we embed each position number of reshaped data with the Xavier initializer to a vector whose length also equals hidden units(Glorot et al.,2010).As FC and PC are the same shapes,we are able to add two matrices corresponding elements together and take the new matrix as our input.

2.3.2.Encoders

Before feeding the new input into the encoders,we use a dropout layer.The dropout layer randomly drops out a certain percentage(P)of model neurons in training,which effectively avoids the model from relying on certain neurons and enhances the generalization ability.After all preparations,data are fed into the encoders of our model.The detailed process of the first encoder is shown in Fig.5.Using the layernormalization after each residual connection layer and working together with feed-forward layers,it prevents the extreme situation,balances the data,and avoid model deterioration.However,compared with CNN,the Transformer model fully relies on the attention mechanism,which is the real workhorse of our model and can model the sequential information well.Specifically,the mission of the self-attention layers is to find out the dependency between different positions in the input that can better encode the information on the current position.For example,in the NLP field,the word“he”in the sentence“The boy didn’t go to school,because he felt sick.”may pay the highest attention to the word“boy”.In this task,the P-wave may pay the highest attention to the S-wave.Furthermore,the multi-head attention mechanism operates multiple attention in parallel and enables the model to explore more relationships than single attention which greatly strengthens the ability to receive information on different positions and decreases the training time.Take the attention mechanism as mapping a query and corresponding key-value pairs to an output,we have the core formula as:

whereQ，K，Vare matrices packed respectively by a batch of queries,keys,values (Alammar,2018).Note that we explore the internal relationship among waveforms;thus,Q=K=Vis set before they go through their own linear layers and begin the attention calculation.

2.3.3.Output

To turn the vectors into a classification result,we use a linear layer,which uses the Xavier initializer andbias=0.1.We then then add a softmax layer to get the possibility of each label.Finally,our model will output the result label which has a higher probability.As we mentioned before,ri=0 represents a noise wave andri=1 represents an earthquake event with P-wave arrival.

3.Results and discussion

In our experiments,we split the dataset into 2 contiguous parts:the first 5/6 samples are collected in July for training,and the last 1/6 are prepared for validation(1/12)and test(1/12).Since the dataset is severely imbalanced(i.e.,much more negative samples than positive samples),the data-driven model will prefer the majority class,which greatly decays the performance of machine learning models.To balance the training data,we produce additional earthquake waveforms by adding zero-mean Gaussian noise(i.e.,set signal-to-noise ratio randomly between 20 and 80 dB) to the original positive waveform(Perol et al.,2018).After that,the number of positive and negative samples in the training dataset is equal(Zhou and Liu,2005).

In the training process,considering the possibility of label error and the robustness of our model,we use the label smoothing algorithm for our model(Szegedy et al.,2016).Specifically,given a smoothing rate ε,a uniform distribution appears over labels,the label is reformulated as

whereKis the number of label channels.Here,since the original labelliis equal to ｛0，1｝ or ｛1，0｝,K=2,and ε=0.1;thus,the reformulated positive label=｛0.95，0.05｝,while the negative one is equal to｛0.05，0.95｝.Besides,we utilize thel2-regularized cross-entropy loss function of λ=10-3to quantify the gap between label and output (Ng,2004).To narrow the gap,we optimize the model parameters using Adam Optimizer with β1=0.9,β2=0.98,and ε=10-8(Kingma et al.,2014).Meanwhile,the dropout rate is also set to 0.1.

3.1.Evaluation metrics

In our experiments,we evaluate our approach with four widely-used evaluation metrics,namely Accuracy,Precision,Recall,and F1 Score.In particular,the definition of each evaluation metric is listed as follows:

where the definitions of TP,TN,FN,and FP are shown in Table 2.Indeed,the Accuracy quantifies the percent of true samples;the Precision focuses on how many positive samples our model judges are correct;the Recall focuses on the percentage of positive samples that can be detected correctly;and the F1 score aims to strike a balance between Precision and Recall,which encourages us to find the earthquake examples accurately and completely.

Table 4 The overall performance.

Table 5 Performance with different time window.

3.2.Baselines

To better validate the effectiveness of our model,we have compared it with numerous state-of-the-art machine learning models that are widely used in the seismic field.Specifically,we first select some traditional machine learning models includingLogistic Regression,SupportVector Machine(Cortes et al.,1995),andRandom.

Forests(Breiman,2001).Then,since TransQuake is a deep learning-based sequential model,we also select some advanced deep learning methods as baselines.

• Recurrent Neural Network(RNN)is a classic neural network structure for sequential data modeling (Rodriguez et al.,1999).In ourexperiments,we choose the most advanced RNN-based model Long short-term memory(LSTM)as the baseline.

• CNN structure is widely applied for seismic detection,due to its remarkable capability of modeling local and static features.Among existing CNN models,we choose two advanced models,ConvNetquake and CPIC (Zhu et al.,2019),as baselines.Specifically,ConvNetquake is the first CNN model proposed for earthquake detection.In our experiments,we add the number of 1d-convolution layers from 8 to 11,and other parts are set following the criterion proposed by Perol et al.(2018).

• MSDNN and MSDNN+Multi-task Learning are two advanced deep learning models for earthquake detection.Particularly,MSDNN combines the advantage of both CNN and LSTM,while MSDNN+Multi-task Learning can take advantage of the homologous earthquakes to enhance MSDNN with additional clues (Zhang et al.,2019).

3.3.Overall performance

After using the validation dataset for fitting our model,we can obtain the best settings of model parameters as shown in Table 3.

The overall performance comparisons of different models on the test dataset are shown in Table 4.The results suggest the following conclusions.First,our approach TransQuake consistently outperforms all baselines.Particularly,in terms of the F1 score,the performance of TransQuake is 259% higher than RNN,131% higher than Random Forest,12.7% higher than CPIC,12.4% higher than ConvNetquake,4.7%higher than MSDNN,and 2.1% higher than MSDNN+Multi-task learning,respectively.Notably,as mentioned before,since the negative samples in our dataset are noises misclassified as aftershocks by low thresholdFilterPickerrather than the simple background noises,the difficulty of our task is much higher than previous studies(Lara-Cueva et al.,2017;Perol et al.,2018;Rouet-Leduc et al.,2019),reasonably resulting in slightly lower performance with weaker baselines compared with those obtained by the original paper (using background noises as negatives).In turn,it validates the difficulty of our task and the discrimination of our model.

Table 3 The setting of model parameters in our experiments.

Fig.6.The process of generating an average attention-simplified view.(a) Average attention-simplified view;(b) 8-head attention-simplified view.

Fig.7.Average attention-simplified views of different waveforms.

Fig.8.Two views of average attention map of all TP examples and all TN examples.(a)attention-full view of TP;(b)attention-simplified view of TP;(c)attention-full view of TN;(d) attention-simplified view of TN.

3.4.Discussion of the length of the time window

As shown in Table 5,we set 4 different time windows from 20 s to 50 s.Generally,metrics,especially F1,become better with the increasing of the time window length,indicating that the information besides the P wave also contributes to detection.Given the efficiency,we set the time window as 50 s in the paper instead of increasing the length.In addition,short time windows of our model still perform better than most state-of-the-art models.Thus,it is possible to use our model in situations with different time requirements.

3.5.Case study:the visualization of multi-head attention

As a crucial part of our model,the attention mechanism can help the model pay more attention to the important information and avoid noise interference,providing us a way to understand the judgment logic behind the model.Here,we useT=50,k=100 positions,which means each position corresponds to a 0.5 s time window of the waveform.To appropriately show their attention for different purposes,we make two views of attention maps:one consists of 100 rectangles and briefly presents the focuses (attention-simplified view),the other is a 2-D matrix with size 100×100 which helps us further explore the inner relationship(attention-full view).Using attention-simplified views,Fig.6 shows an example of the distribution of attention scores.The darker the color,the higher the attention score (i.e.,more important for model judgment).It can be seen that our model spontaneously focuses on the P wave and S wave which are important for earthquake waveforms.Specifically,Fig.6(b) shows how our multi-head attention mechanism works.The multi-head attention uses different attention distributions to focus on different important positions(e.g.,Head 1 focuses on S-wave and Head 3 focuses on P-wave).This kind of structure can make the attention mechanism more precise and comprehensive.Then,Fig.7 shows some attention maps of different waveforms,where we can find that our model always accurately focuses on the key positions(i.e.,P-wave and S-wave),no matter how the waveforms change.

Furthermore,in order to analyze the effect of the attention mechanism on waveforms of different types in detail,we plot the attention of all TPs (600 aftershocks) and that of all TNs (8 967 noises) respectively in Fig.8.In attention-simplified views (i.e.,Fig.8(b) and (d)),they both highly focus on the P-wave and exhibit a slight difference in the attention distribution of the S-wave.In attention-full views(i.e.,Fig.8(a)and(c)),they not only present the high attention position but also tell the positions are paying attention to it (i.e.,the latent correlation between different positions of waveforms).In the attention-full view,each square denotes the attention weight ofYtoX.As shown in Fig.8(c),TN(noise)has more diagonal focuses,which means it usually changes its focus for different waves rather than builds stable relationships with other positions.Another important difference is shown in the red boxes in Fig.8(a)and(c):earthquakes own their most stable relationship between the coda wave and P-wave while noises do not (it meaninglessly focus on the diagonal of their first motion).Besides,the high weights of S-wave and coda wave paying to P-wave in Fig.8(a) further validate that the existence of the S-wave and coda wave contributes to identifying the P-wave,especially the coda wave.

4.Conclusion and prospection

In this paper,we propose a novel deep learning approach,Trans-Quake,for enhancing the performance of seismic P-wave detection.Specifically,TransQuake is designed on the basis of the most advanced deep sequential model,namely Transformer,integrating both STA/LTA algorithms and multi-head attention mechanisms.Compared with other deep learning models in relevant studies,TransQuake can effectively model the globally dynamic information contained within seismic waves.In addition,it shows high interpretability for result investigation.We have conducted extensive experiments on the aftershocks following the Mw7.9 Wenchuan earthquake using a number of state-of-the-art baselines.The experimental results clearly validate the effectiveness of TransQuake in terms of both seismic P-wave detection and result interpretation,such as the attention distribution of seismic waves in different positions,and the analysis of the optimal relationship between coda wave and P-wave for noise identification.Indeed,the architecture of TransQuake can be easily adapted to multidimensional data from various fields.Therefore,the future works will focus on the integration of different data resources(e.g.,geomagnetic and geothermal data)into our model for seismic detection.