基于獨立分量分析的工頻干擾消除技術

（中國地質大學 地球物理與信息技術學院，北京 100083）

摘 要：簡要介紹了ICA的基本原理和快速算法，在分析地震信號和工頻干擾特點的基礎上，利用ICA技術來消除地震記錄中的工頻干擾，并與常規方法進行比較。研究結果表明ICA在有效消除工頻干擾的同時，能夠保護有效信號，并且在提高資料的信噪比方面更有優勢，具有良好的應用前景。

關鍵詞：獨立分量分析；盲源分離；快速獨立分量分析算法；工頻干擾

中圖分類號：TP391 文獻標志碼：A

文章編號：10013695(2009)01022703

Power interference removal method based on independent component analysis

WEI Wei，LIU Xuewei

（ School of Geophysics Information Technology， China University of Geosciences， Beijing 100083， China）

Abstract:This paper studied the fundamentals and fast algorithm of ICA， tried to remove the power interference in seismic data on the base of analyzing the features of seismic data and the power interference. Compared with the results of the conventional methods， ICA has good effect on removing the power interference and can protect the available signals effectively. Results show that ICA has better effect on improving the signal to noise ratio and has a good application perspective.

Key words:independent component analysis（ICA）; blind source separation（BSS）; FastICA algorithm; power interference

在地震勘探數據采集時，如果有工業電傳輸線穿過勘探區域，則在采集的地震數據中將包含工頻干擾，這勢必會影響地震勘探的精度。目前壓制工頻干擾主要采用兩種方法：a）陷波濾波法。它是在頻率域中實現的，設計多吸收點的、具有一定阻帶寬度的陷波濾波器組，用于消除工頻干擾的基波和諧波成分。如果信號頻譜與工頻干擾的頻譜有混疊，則陷波濾波器在濾除工頻干擾的同時也會造成目標信號的損失。b）自適應濾波法。該方法的前提是目標信號和干擾波不相關。在此基礎上，利用自適應濾波算法(最小均方或最小二乘法)，自動調整濾波器系數，以跟蹤輸入過程的變化，并實現工頻干擾的自適應抵消。這種方法在時間域進行，其本質與陷波濾波是一致的[1]。

ICA是近期發展起來的一種非常有效的盲源分離技術。其處理的對象是一組相互統計獨立的信源經線性組合而產生的混合信號，最終從混合信號中提取出各獨立的信號分量。目前ICA 在生物醫學信號處理、混合語音信號分離、圖像消噪等方面已取得了較好的應用效果[2~5]。本文在對ICA的基本原理和快速算法的研究基礎上，分析了地震數據和工頻干擾的特點，深入研究了基于ICA的工頻干擾消除技術，并與常規方法進行比較，結果表明ICA技術可以有效消除工頻干擾，并且不傷害有效信號在提高資料的信噪比方面更有優勢，結果令人滿意。

1 ICA基本原理

下面以兩個信號的盲源分離為例介紹ICA的基本原理。在圖1中，s1(t)、s2(t)是兩個非高斯分布的獨立源；x1(t)、x2(t)是兩個觀察信號， 混合矩陣A是未知的模型系數。

x1(t)x2(t)=a11 a12

a21 a22×s1(t)s2(t)(1)

或用矩陣形式描述為

X=A×S（2）

在獨立源信號S和混合矩陣A都是未知的情況下，希望能尋找到一個分離矩陣W=A-1從觀測信號中進行源信號的分離，即

y=A-1×X=W×X=S(3)

其中：X、S分別為觀測信號和源信號向量；A為混合矩陣。

ICA方法可以有效解決上述盲源分離問題，其基本思想就是中心極限定理，即一隨機量如由許多相互獨立的隨機量之和組成，只要各獨立的隨機量具有有限的均值和方差，則不論各獨立隨機量為何種分布，該隨機量必接近高斯分布。由此可以推斷，式(1)中觀測信號xi較之si更接近高斯分布，或si比xi的非高斯性更強。 因此可以在分離過程中測量si的非高斯性，當非高斯性度量達到最大時，表明已完成對各獨立分量的分離。ICA正是通過對分離結果獨立性的度量，使用某種優化算法來尋求分離矩陣W。ICA的關鍵就是一個可度量分離結果獨立性的目標函數及相應的優化算法。人們從不同角度提出多種目標函數及算法，本文所討論的是一種基于負熵的判據和一種非常高效的ICA 算法，即快速ICA(FastICA)算法[7，8]。

2 快速ICA 算法

下面分別從目標函數和優化算法兩個方面對FastICA作介紹。由極限定理可知，分離結果的非高斯性最大與信號統計意義上的獨立是等價的，所以可以通過度量非高斯性的函數作為獨立性的判據。FastICA算法通常用負熵作為度量分離結果獨立性的目標函數

J(y)=[E{G（y）}-E{G(yGauss)}]2（4）

其中:E(#8226;)為均值運算，G(#8226;)可取

G1（u）=(1/a1)×lg cos (a1u)G2(u)=(1/4)×u4（5）

ICA的統計特性取決于所選的非高斯性度量目標函數的性質，而算法的收斂性、穩定性和有效性取決于所構造的最優化算法。FastICA算法大致分為兩個步驟：a）對觀測數據x進行去白化處理，假設白化后的信號為，則E[×T]=I。用類似主分量分析方法(PCA)即可實現信號的白化處理，其目的是去除各觀測信號間的相關性，簡化了后續獨立分量提取算法[9]。b)對白化后的信號進一步處理，即尋找分離矩陣，以實現獨立分量提取。這一過程是一迭代逼近過程，變量n表示迭代步數，令i(n)是源信號中S中的某一分量，ωi(n)為分離矩陣W中與i(n)對應的某一行向量，即

s^i(n)=i(n)× (n=1，2，3，…)（6）

分離過程中，用式(4)所定義的目標函數對分離結果i(n)的非高斯性進行度量，并對ωi(n)進行調整，FastICA 算法的調整式為[8]

i(n+1)=E{G′(Ti(n))}-E{G″(Ti(n))}i(n)（7）

當相鄰兩次的ωi(n)無變化或變化很小時，即可認為i(n)≈si，迭代過程結束。式(7)中的均值計算可通過時間平均獲得，所要注意的是：每次迭代后， 都要對ωi(n)進行歸一化處理:

ni=ni/‖ni‖(8)

以確保式(6)的分離結果具有單位能量， 這是ICA 問題的一個基本限制。對于多個的獨立分量，可重復使用上述過程進行分離[2]。 但每提取出一個獨立分量后，要從觀測信號中減去這一獨立分量，如此重復，直至所有獨立分量完全分離。

3 ICA在地震信號工頻干擾消除的應用

在野外采集地震數據時，探區內存在嚴重的工頻干擾，嚴重影響地震資料的品質。含有工頻干擾的地震信號表示為

x(t)=a11×s(t)+a12×n(t)=a11×s(t)+a12×A×sin(2πf0t+θ0)(9)

其中：s為有效信號；n(t)為工頻干擾，參數f0、A、θ0分別為工頻干擾的頻率、幅度和相位。a11與a12是混合矩陣的常系數。若設地震記錄中的有效信號和工頻干擾滿足ICA的應用條件，即源信號中最多只允許一個高斯信號；源信號之間相互統計獨立，則可嘗試用ICA 進行數據處理。地震記錄當中的有效信號可以看做是地震子波與地層反射系數序列的卷積，事實已經證明地層反射系數是滿足非高斯分布的，所以地震信號也是非高斯分布的。又因為有效信號和工頻干擾信號分別由不同的信源產生，因此彼此相互獨立，可見地震信號和工頻干擾的特點滿足ICA的應用條件。

根據前面對ICA 問題的描述可知，如果要分離出s(t)與n(t)，還需包含工頻干擾和有效信號的另一道觀測信號x2(t)=a21×s(t)+a22×n(t)作為ICA的另一個輸入。但在實際中，不可能找到這樣一組觀測信號，這就需要構造一個觀測信號，由于s(t)是未知的，而工頻干擾信號是規則的，且滿足n(t)=A sin(2πf0t+θ0)，在不考慮振幅的影響下，若已知工頻干擾的頻率、初始相位，很容易構造出符合要求的另一路觀測信號x2(t)=0×s(t)+a22×(sin 2πf0t+θ0)。工頻干擾的頻率是已知的（一般為50 Hz）或可以精確估計出來[2]，為了避免估計初始相位θ0，文獻[2]給出了一種人工構造觀測信號的方法。含有工頻干擾的觀測信號進一步表示為

x(t)=a11×s(t)+a12×n(t)=

a11×s(t)+a12×A sin (2π×f0t+θ0)=a11×s(t)+a12×Acos θ0 sin 2π×f0t+a12×A sin θ0 cos 2π×f0t（10）

或表示為

x(t)=a×s(t)+b×sin 2πf0t+c×cos 2π×f0t(b=a12×A cos θ0，c=a12×A sin θ0)（11）

可見，式（11）將含有未知相位的工頻干擾信號用兩個初始相位為零的，相互正交的工頻信號來代替，因此，可通過構造另外兩組參考信號p1=sin 2πf0t和p2=cos 2πf0t，與含有工頻干擾的地震信號共同作為ICA的輸入，分離出有效地震信號s(t)，這樣就避免了對θ0的直接估計，其數學模型表示為

xp1p2=a b c

0 1 0

0 0 1×s(t)sin 2πf0tcos 2π f0t

或X=A×S（12）

用FastICA對式（12）中的混合矩陣A和分離矩陣W=A-1進行估計，并利用式（3）求出不含工頻干擾的地震信號s(t)，即可實現工頻干擾的消除。

為了考察方法的有效性，先提取一道不含工頻干擾的地震信號（圖2），加入工頻干擾模擬信號n1=10-4×sin(2π×50×t+π/3)，得到包含50 Hz工頻干擾的地震信號如圖3所示。

根據上述所提出的方法，首先要構造參考信號源sin(2π×50×t)和

cos(2π×50×t)，并與含工頻干擾的地震信號一起構成ICA 算法的三路輸入信號，如圖4 所示。利用MATLAB編制了前面介紹的FastICA算法程序，并對圖4 所示輸入進行盲源分離，ICA 輸出波形如圖5 所示。可見，ICA成功提取出了有效地震信號。

為了與常規方法相比較， 分別采用幅頻響應和相頻響應如圖6 所示的50 Hz陷波濾波器和基于NLMS算法的自適應噪聲抵消技術對圖3所示信號進行處理，結果如圖7所示。可以看出，ICA方法的結果最接近實際有效地震信號，在消除工頻干擾方面更有效。從三種方法處理結果與實際地震信號的頻譜分析圖中，如圖8所示，可以看出自適應噪聲抵消和頻率域濾波方法的結果不僅消除了工頻干擾，而且也損失了相同頻率的有效信號，而ICA的結果更接近實際地震信號的頻譜。由此可知，ICA方法不僅成功消除了工頻干擾，并且不傷害有效信號。ICA方法之所以較其他兩種方法有如此的優勢，是因為獨立分量分析是在統計獨立意義下，來對混合信號進行分離的，所以不存在傷害有效信號問題，而頻率域濾波方法則會造成相同頻率的有效地震信號的缺失，自適應濾波方法，則會由于算法收斂需要一定的時間，導致在開始一段時間內工頻干擾的消除效果較差，提取的有效信號存在一定的失真。

另外通過計算，處理前信號的信噪比為0.580 3 dB，自適應噪聲抵消后信號的信噪比11.324 3 dB，頻率域濾波后信噪比變為9.177 3 dB，而ICA的輸出結果的信噪比為32.599 3 dB，信噪比提高了大約32 dB，較頻率域濾波方法、自適應噪聲抵消算法分別提高了23 dB、21 dB，可見ICA方法在提高資料的信噪比方面要優于常規方法。 

4 結束語

在野外進行地震數據采集時，經常受到工頻干擾，通常的濾波方法在消除工頻干擾的同時會濾除一些有效信號，而自適應噪聲抵消技術達到穩態需要一定的時間，因此結果也難以令人滿意，而ICA可以在統計獨立意義下對混合信號進行分離，以達到消除噪聲的目的。本文簡要介紹了ICA的基本原理和快速算法，通過對地震信號和工頻干擾的分析得知，地震記錄中的有效信號和工頻干擾信號在統計上相互獨立，且有效地震信號服從非高斯分布，滿足ICA的適用條件。在仿真實驗中，為了利用FastICA消除地震記錄中的工頻干擾，構造了兩兩相互正交的參考信號作為ICA輸入，并將ICA與傳統方法的結果進行比較，從波形上觀察ICA方法的結果最接近有效地震信號，頻譜分析則進一步說明ICA在有效消除工頻干擾的同時，能夠很好地保護有效信號，而且ICA方法在提高資料的信噪比方面更有優勢。總之，ICA在消除地震信號中的工頻干擾方面具有良好的應用前景。

參考文獻：

［1］劉洋.強工頻干擾波的提取與消除方法[J].石油物探，2003，42(2):154159.

[2]吳小培，詹長安，周荷琴，等.采用獨立分量分析的方法消除信號中的工頻干擾[J].中國科技大學學報，2000，30(6):671676.

[3]甘武，孫云蓮，蒲曉羽.用獨立分量分析消除工頻通信中的諧波干擾[J].電力系統通信，2005，26(148):2528.

[4]劉長生，唐艷，湯井田.基于獨立分量分析的腦電中眼電偽跡消除[J].計算機工程與應用，2007，43(17):230232.

[5]李鴻燕，郝潤芳，馬建芬，等.基于維納濾波和快速獨立分量分析的有噪混合圖像盲分離[J].計算機應用研究，2007，24(10): 161162，165.

[6]COMMON P.Independent component analysis， a new concept[J].Signal Processing，1994，36(3):287314.

[7]HYVARINEN A，OJA E.A fast fixedpoint algorithm for independent component analysis[J].Neural Computation，1997，9(7):14831492.

[8]HYVARINEN A，OJA E.Independent component analysis:algorithm and application[J].Neural Network，2000，13(4):411430.

[9]張發啟，張斌，張喜斌.盲信號處理及應用[M].西安：西安電子科技大學出版社，2006:1743.