周海軍，李 磊

(滁州學院 電氣學院，安徽 滁州 239000)

一、引言

在科學技術(shù)迅速發(fā)展的今天，人類的活動越來越頻繁，同時，地震記錄儀響應的帶寬顯著擴大，能檢測到頻率更低的或更高的、能量更加微弱的地下震動信號，而實際記錄到的信號往往同時包含多種不同來源的環(huán)境噪聲信號[1]。不確定的噪聲信號對實際測得的地震信號的分析會產(chǎn)生很大的影響，因此，對強噪聲下的實測地震信號的去噪研究對地震定位，地震預測與識別，地殼速度結(jié)構(gòu)反演等具有十分重要的現(xiàn)實意義。

二、小波閾值去噪原理

(一)理論依據(jù)

假設原始地震信號為s(t)，受到噪聲干擾后為地震信號為x(t)，其基本公式模型表示為：

x(t)=s(t)+σe(t)

(1)

式中σ為噪聲強度，e(t)為干擾噪聲。

對帶噪聲的地震信號進行小波變換，分解后的各層信號去噪使用的閾值是根據(jù)原信號的信號和噪聲比值來確定的。小波變換過程可以分三個階段[1]14-22：

(1)小波分解：系統(tǒng)會選定一種小波基，對帶噪聲的地震信號進行N層小波分解。用小波分解函數(shù)對信號進行分解的過程就是依次把信號的各種頻率成分從高到低逐步分離為不同頻帶的過程，對于db4小波而言，其頻帶近似分為均等的幾份[2]。

(2)計算分層的閾值：對小波分解得到的各層小波系數(shù)選取一個閾值，并用閾值函數(shù)對小波系數(shù)進行處理。在高斯模型中，維數(shù)趨向無窮時，噪聲的幅值幾乎小于閾值thr，見公式2。

(2)

其中，σ為噪聲標準差，N為去噪信號的長度值。

(3)重建信號：將去噪處理后的小波分解系數(shù)通過小波重建方式恢復原始地震信號。

(二)閾值函數(shù)的選擇

系統(tǒng)選取適合地震信號去噪的閾值對分解后的小波系數(shù)進行運算是小波閾值去噪的最主要問題。硬閾值法和軟閾值法是兩種比較常用的閾值處理方法。硬閾值法是把信號經(jīng)過小波分解后的小波系數(shù)的絕對值與理論閾值進行比較，大于閾值的小波系數(shù)保持不變，小于或等于閾值的小波系數(shù)設置為零，其函數(shù)如公式(3)所示[1]14-22。軟閾值法是把信號經(jīng)過小波分解后的小波系數(shù)中大于閾值的點變?yōu)樵擖c與閾值的差值，小波系數(shù)等于或小于閾值的點設置為零,其函數(shù)如式(4)所示[3]。它們的基本思想就是去除較小的小波幅值系數(shù)，對幅值較大的系數(shù)進行收縮或者保留[4]。

硬閾值函數(shù)表示為：

(3)

其中，λ為閾值，ωj,k為小波系數(shù)。

軟閾值函數(shù)表示為：

(4)

其中，λ為閾值，sgn(ωj,k)是符號函數(shù)，ωj,k為小波系數(shù)。

(三)改進的小波閾值選擇

去噪閾值的選取方法是小波閾值去噪的關(guān)鍵。為了使小波閾值去噪法更適合強干擾的地震信號去噪。本文主要對閾值的運算進行了改進。分別計算出小波分解后的每一層系數(shù)對應的閾值[1]14-22。如公式(5)所示：

(5)

其中，thrm表示第m層的閾值，m是表示小波分解后的第m層，n是天然地震信號的長度，nm是小波分解m層信號的長度，σm是m層天然地震信號的噪聲強度，σ是天然地震信號的噪聲強度[1]14-22。

強干擾下的地震信號采用普通小波閾值去噪的步驟跟采用改進小波閾值去噪的步驟不同之處在第二步求解分層閾值的方法，改進的小波閾值去噪法的閾值計算是根據(jù)公式(5)求得。

三、強干擾的地震信號去噪

對采集的原始地震信號預處理之后，用db4小波進行45層分解[5]。如圖1所示，地震原始波形信號、分解的低頻部分和分解的高頻部分的系數(shù)圖。從圖中可以得知地震信號的有效信息主要集中在分解的低頻部分、分解的第一層波形和分解的第二層波形，分解的第三層波形和分解的第四層波形實際有用信息很少，可以忽略不計。

圖2中最上面的一幅圖是某臺站記錄的某次天然地震事件垂直方向的原始地震波形圖，中間一幅圖是用默認的小波閾值去噪后的地震波形圖，最下面一幅圖是用改進的小波閾值去噪的地震波形圖。對比兩種去噪方法，改進的小波閾值去噪法更加強化天然地震信號的特征信號，在識別的時候提取的波形特征信息更加明顯。

四、特征提取

對強干擾下的地震信號進行小波去燥后，進行EMD分解，再取第一個IMF分量，對這個分量進行Hilbert黃變換，變換后的數(shù)據(jù)進行歸一化處理，再提取特征。提取本文所用到的兩個特征，即HHT譜幅比特征和自相關(guān)函數(shù)的均值Meanxcorr。

圖3 特征提取流程圖

(一)HHT譜幅比特征

(6)

其中A表示IMF分量經(jīng)過Hilbert變換后的值，Amaxs是max(A(1:n))，Amaxp是max(A(n:2000))，n為S波的起始點。

(二)自相關(guān)函數(shù)的均值Mc特征

(7)

(8)

其中A表示IMF分量經(jīng)過Hilbert變換后的值。

五、線性不可分支持向量機分類

地震數(shù)據(jù)從長時看是非線性的，需要將低維非線性數(shù)據(jù)映射到高維線性數(shù)據(jù)，然后對天然地震和人工爆破進行分類，因此，考慮假設集是這種類型的函數(shù)：

(9)

這里的關(guān)系φ:x→f是從輸入空間到特征空間的映射[6-8]。天然地震和人工爆破數(shù)據(jù)通過非線性支持向量機會得到一個分類器，過程如圖4所示。

圖4 建立分類器

實驗采用的數(shù)據(jù)來自人工挑選的首都圈地震帶的300個天然地震事件數(shù)據(jù)集和300個人工爆破數(shù)據(jù)集，合計600份地震波形數(shù)據(jù)。從天然地震數(shù)據(jù)和人工爆破數(shù)據(jù)集合中隨機的各選取200份數(shù)據(jù)作為訓練模型數(shù)據(jù)，余下100份天然地震數(shù)據(jù)和100份人工爆破數(shù)據(jù)作為識別源。數(shù)據(jù)分好類之后，將每份數(shù)據(jù)進行EMD分解，取第一個IMF分量，經(jīng)過希爾伯特變換后歸一化處理，提取HHT譜幅比、自相關(guān)函數(shù)均值作為特征，利用非線性支持向量機建立天然地震和人工爆破的模型，再根據(jù)模型識別未知的震源類型。

(1)最大識別率

ri,max=max(r(i,1),…r(i,200)),i=1，2

(10)

(2)最小識別率

ri,min=min(r(i,1),…r(i,200)),i=1,2

(11)

(3)平均識別率

(12)

(4)識別率標準差

σri=std(r(i,1),…r(i,200)),i=1,2

(13)

其中i表示第幾個IMF分量。在比較采用普通小波閾值方法和改進小波閾值方法去噪的識別中，上述四個公式中i都取1。經(jīng)過200次反復隨機抽樣后，采用普通小波閾值方法和改進小波閾值方法去噪的識別實驗的結(jié)果如表1所示。

表1 基于兩種去噪方法的第一個IMF分量提取的特征識別效果的統(tǒng)計量

從表1可以看出，其他條件一致，文中提出的采用改進的小波閾值去噪法處理數(shù)據(jù)后，對天然地震數(shù)據(jù)和人工爆破數(shù)據(jù)識別，在最大識別率、最小識別率、平均識別率、識別率標準差等四個方面都比未改進的去噪法處理數(shù)據(jù)效果好。

六、結(jié)論

通過上面實驗和地震波形理論可得出，對強干擾的地震信號去噪，小波變換可以保留地震信號中的有效信息，用改進的小波閾值去噪法對采集的地震信號去噪，識別效果更好，圖形更符合理論要求。同時，在小波域閾值去噪中，小波基選擇和閾值的選取對去噪效果影響很大，因此對于不同的地震帶要選擇合適的閾值和小波基進行去噪，才能在后期的地震震源類型的識別中取得更加準確的值。

參考文獻：

[1] 周海軍．基于地震波形的震源類型自動識別研究[D]．桂林：廣西師范大學，2012：14-22.

[2] 肖本夫,萬永革,祁玉萍．強干擾環(huán)境下有效數(shù)字地震信號的提取[J]．華北地震科學,2011，29(1)：6-9.

[3] DONOHO D L．De-Noising by Soft-Thresholding [J].IEEE Transaction Information Theory,2002，41(3):613-627.

[4] 石巖峻．地震資料的小波去噪方法[J]．內(nèi)蒙古石油化工,2010(6):47-48.

[5] 潘泉，張磊，孟晉麗，等．小波濾波方法及應用[M]．北京：清華大學出版社，2005：123-134.

[6] 鄧乃揚，田英杰．支持向量機：理論、算法與拓展[M]．北京：科學出版社，2009：87-112.

[7] 邊肇祺，張學工．模式識別[M]．北京：清華大學出版社，2000：231-242.

[8] 王宏強，尚春陽，高瑞鵬，等．基于小波系數(shù)變換的小波閾值去噪算法改進[J]．振動與沖擊,2011,30(10):165-168.