基于多屬性的P波到達(dá)時間自動拾取方法*

辛 平， 李志華

(江南大學(xué) 物聯(lián)網(wǎng)工程學(xué)院 計算機(jī)科學(xué)系，江蘇 無錫 214122)

0 引 言

在地震震相檢測中，首先要做的是精確拾取震相到時，即精確拾取P波到時。地震波的準(zhǔn)確拾取在許多與地震相關(guān)的研究中起著關(guān)鍵作用。目前，自動拾取方法已占據(jù)主導(dǎo)地位[1]。文獻(xiàn)[2]用簡單的振幅和能量閾值來確定P波到時;由于該算法只用了振幅和能量是否大于閾值判定P波到時，所以其拾取的精確度不高；文獻(xiàn)[3]采用長短時平均(short term average/long term average,STA/LTA)方法確定地震震相到時,可以快速拾取震相到時，但其拾取震相到時的精確度受噪聲影響比較大，高噪音時出現(xiàn)誤拾甚至拾取不到P波。此外經(jīng)驗閾值的選取直接影響STA/LTA拾取的準(zhǔn)確度；文獻(xiàn)[4]采用最大似然(maximum likelihood,ML)估計方法來確定地震P波到時，雖然能夠比較精確確定P波到時，但ML方法認(rèn)為噪聲干擾是一成不變的，忽略了突發(fā)人為或者其他非地震因素所引起的噪聲，這顯然是不合理的；文獻(xiàn)[5]通過計算P波偏振度來確定P波到時,該方法提高了識別可靠度，但計算量大；文獻(xiàn)[6]將STA/LTA、偏振分析、AR分析以及瞬時頻率分析這4種方法結(jié)合起來確定P波到時，雖然精確度上有了極大的提高，但算法的計算量很大；文獻(xiàn)[7]采用自回歸方法來拾取P波到時，能夠比較精確確定P波到時；文獻(xiàn)[8]提出了一種基于局部極大值分布的P波到時自動識別方法，算法將能量和頻率特征結(jié)合起來，能夠在高噪音時比較精確拾取P波。

本文提取地震信號的15維特征作為數(shù)據(jù)樣本，并通過特征選擇來對數(shù)據(jù)樣本進(jìn)行降維處理，最后利用支持向量機(jī)(support vector machine,SVM)進(jìn)行分類預(yù)測，確定P波到時。

1 基于改進(jìn)的局部線性嵌入和核分類方法的P波拾取

1.1 多維特征的構(gòu)造

為了獲取豐富的地震信號特征，本文通過對實測地震信號進(jìn)行變換和處理，從時域和頻域兩方面綜合選取15維特征利于數(shù)據(jù)深度挖掘。15維特征包括震相、振幅、偏振度、分形維數(shù)、STA/LTA，ML、振幅比、能量比、曲線長度比、峭度、偏斜度、均方根振幅、平均能量、振幅峰態(tài)，頻域參數(shù)中選取頻率。具體如下:

1)震相到時由STA/LTA確定。

2)振幅從信號中直接讀取。

3)偏振度F采用文獻(xiàn)[5]中的方法計算，F(xiàn)=m1/m2,m1和m2分別為由地震信號數(shù)據(jù)形成的3×3階矩陣的最大特征值和次最大特征值。其矩陣構(gòu)造為

(1)

式中x，y，z為地震信號3分量數(shù)據(jù)，cov(x,y)為x和y的協(xié)方差

(2)

式中μx和μy為x和y的平均值。

4)分形維數(shù)D的計算為D=1-S,S為不同步長和其對應(yīng)的曲線長度擬合的直線斜率。

5)i時刻振幅比為

(3)

式中x為地震數(shù)據(jù)，m為地震數(shù)據(jù)時窗的起點，n為地震數(shù)據(jù)時窗的終點。

6)i時刻能量比為

(4)

7)i時刻曲線長度比為

(5)

式中 Δt為地震數(shù)據(jù)的采樣間隔。

8)峭度為

(6)

式中s為時窗大小。

1.1.1 改進(jìn)的STA/LTA方法

傳統(tǒng)STA/LTA計算為

(7)

式中m為STA窗口大小,n為LTA窗口大小。

為提高震相拾取的準(zhǔn)確度，對STA/LTA方法進(jìn)行改進(jìn),即

(8)

1.1.2 改進(jìn)的ML方法

由于震相初至前后地震信號的振幅有很大的變化，因此可以用ML估計來確定震相的初至?xí)r間。由于地震信號呈隨機(jī)分布，所以本文假設(shè)地震信號變化前后的分布都服從正態(tài)分布，則ML計算為

(9)

式中σ1為震相初至前噪音方差，σ2為震相初至后方差。

(10)

改進(jìn)的ML計算方法概括如下：

3)采用式(10)計算窗口內(nèi)的似然估計值L。

1.2 基于改進(jìn)的LLE特征抽取方法

LLE計算序列間距離采用歐氏距離，考慮每個區(qū)域上各個點分布是不均勻的，為了使每個區(qū)域上各點的分布整體呈均勻化，LLE中距離的計算改進(jìn)為

(11)

式中M(i)，M(j)分別為與其點距離的均值。

改進(jìn)的LLE算法具體描述如下:

1)計算每個樣本點xi的k個近鄰點xij(j=1,2,…,k)。

(12)

3)計算xi和xj的低維映射值yi和yj，使得ε(Y)最小

(13)

式中ε(Y)為損失函數(shù)，式(14)滿足

(14)

M=(I-W)T(I-W)

(15)

1.3 面向時間序列的核分類方法

由于傳統(tǒng)的SVM[10]距離度量采用歐氏距離，計算時間復(fù)雜度低(O(n2))，但其無量綱且只適用于等長的時間序列[11]。由于所研究的地震信號是時間序列，具有多維、離散的特性，且每一維都有其特定的量綱，為了使SVM更好地適用于地震時間序列數(shù)據(jù)，本文提出面向時間序列的SVM(time series SVM,TS-SVM)分類方法，在TS-SVM方法中，時間序列間相似性度量采用動態(tài)時間彎曲距離[12]。描述如下:

假設(shè)有時間序列X和Y，則其相似性度量測度為

(16)

(17)

式中 < >為空序列,Rest(X)={x2,x3,…，xn},Rest(Y)={y2,y3,…,yn}。

TS-SVM方法在核函數(shù)選取上采用高斯核函數(shù)，用DTW距離代替歐氏距離，核函數(shù)定義為

(18)

得到最優(yōu)分類函數(shù)為

f(x)=sgn{(w*,x)+b*}

(19)

新提出的TS-SVM方法描述如下：

1)由計算時間序列X的動態(tài)彎曲距離DTW。

2)計算時間序列的內(nèi)核函數(shù)，求解二次規(guī)劃問題

得到最優(yōu)解a*={a1,a2,…，ak}。

4)計算f(x)的值，實現(xiàn)對時間序列X的分類。

動態(tài)彎曲距離的時間復(fù)雜度為O(n2)，歐氏距離的時間復(fù)雜度也為O(n2)，TS-SVM并未改變傳統(tǒng)SVM的時間復(fù)雜度。

1.4 基于TS-SVM的 P波拾取方法AMPAT

AMPAT流程如圖1所示。

圖1 基于TS-SVM的P波拾取方法流程

AMPAT方法概括如下：

1)數(shù)據(jù)預(yù)處理：計算15維特征構(gòu)造特征空間集，考慮到特征空間的每一維特征數(shù)量級不同，需要進(jìn)行歸一化處理。對歸一化處理的特征空間用LLE進(jìn)行降維，選取最佳特征降維后的維數(shù)。選取合理數(shù)量的樣本數(shù)據(jù)，50 %樣本數(shù)據(jù)作為訓(xùn)練集，50 %作為測試集。

2)TS-SVM參數(shù)尋優(yōu):在訓(xùn)練樣本集上進(jìn)行有監(jiān)督學(xué)習(xí)的訓(xùn)練，確定最佳參數(shù)。

3)構(gòu)造P波拾取模型：對LLE降維后的特征數(shù)據(jù)，采用合理的TS-SVM尋優(yōu)參數(shù)進(jìn)行測試，拾取P波到時。

2 實驗與結(jié)果分析

2.1 實驗數(shù)據(jù)

本文數(shù)據(jù)由IRIS(Incorporated Research Institutions for Seismology)提供，數(shù)據(jù)格式為SAC。以美國東部海岸 2016年 3月份地震數(shù)據(jù)為樣本，利用MATLAB進(jìn)行仿真實驗來驗證本文提出的方法的準(zhǔn)確性。樣本的采樣頻率均為40Hz，記錄了地震信號在某一時間段內(nèi)的加速度變化情況。

2.2 結(jié)果分析

根據(jù)專家經(jīng)驗對樣本信號進(jìn)行人工拾取，其中P波到時如圖2中直線。拾取得P波到時在[161.65,161.675]s之間，精確P波到時為161.652 8 s。

圖2 人工拾取P波

對樣本信號采用STA/LTA和ML方法處理,P波到時分別為161.618 7 s和161.631 5 s,如圖3。

圖3 2種方法處理結(jié)果

AMPAT方法得到P波到時為161.632 6 s。另取IRIS 2016年3月份美國中部海岸的10組地震數(shù)據(jù)，用這4種方法拾取P波到時如表1。可知，以人工拾取結(jié)果為參考，STA/LTA方法拾取誤差比較大，最大為0.338 8 s，ML拾取結(jié)果比STA/LTA方法好，最大誤差為0.061 4 s,AMPAT拾取結(jié)果比STA/LTA和ML結(jié)果要好，最大誤差為0.055 4 s。進(jìn)一步比較，以人工拾取方法結(jié)果為參照，其他3種方法計算平均誤差為：STA/LTA，0.073 7 s；ML，0.0331 s;AMPAT，0.025 1 s。可知，AMPAT的平均誤差要小于ML方法和STA/LTA方法的，ML拾取誤差較好于STA/LTA方法。

表1 拾取10組地震數(shù)據(jù)的P波初至?xí)r間 s

綜上，AMPAT在最大誤差和平均誤差都優(yōu)于STA/LTA和ML方法，其對P波拾取的精度高于這2種方法。

3 結(jié) 論

針對地震數(shù)據(jù)單個特征受噪聲等其他因素影響準(zhǔn)確拾取P波到時難度之大這一問題，本文提取地震數(shù)據(jù)的15維特征并通過基于LLE降維和SVM拾取P波，實驗對比證明，AMPAT方法獲取的P波到時更加精確。由于AMPAT方法選取了15維特征，導(dǎo)致計算量較多，如何有效降低計算復(fù)雜度是下一階段的研究目標(biāo)。