基于特征提取與聚類識別的人工源電磁偽隨機信號處理方法

張 賢 李帝銓* 李 晉 胡艷芳

(①中南大學有色金屬成礦預測與地質環境監測教育部重點實驗室，湖南長沙 410083；②中南大學地球科學與信息物理學院，湖南長沙 410083；③湖南師范大學信息科學與工程學院，湖南長沙 410081)

0 引言

電磁法是一種重要的地球物理勘探方法，按照源的種類可分為天然源電磁法(MT、AMT等)和人工源電磁法(CSEM、TEM等)[1-3]。相比天然源電磁法，人工源電磁法克服了天然場源信號弱且隨機性強的缺點[4-5]，電磁信號具有更高的信噪比和分辨率?？煽卦匆灾芷谛缘姆讲ê蛡坞S機信號為主。廣域電磁法(Wide Field Electromagnetic Method, WFEM)擁有完備的理論體系及成熟的儀器設備，中國對其具有完全自主知識產權[3,6]。該方法通過一次激發—接收可獲得多個頻率的地電信息，提高了野外工作效率和抗干擾能力[7-8]。與可控源音頻大地電磁法(CSAMT)相比，WFEM摒棄了只能在“遠區”觀測而導致信號微弱的劣勢，將“過渡區”與“遠區”有機統一，顯著擴大了觀測適用范圍以及探測深度[9]；與磁偶源頻率測深法(MELOS)相比，WFEM保留了計算公式中的高次項，定義了適合全域的廣域視電阻率，不需要對測量結果進行“校正”，且只需要觀測一個分量就可以獲得視電阻率[3,6,10-11]。目前，WFEM已經在頁巖氣勘探[12-14]、金屬礦勘探[15]、工程物探及城市物探等領域得到越來越廣泛的應用。

隨著現代工業與科技的發展及人類活動范圍的不斷擴大，電磁干擾越來越強，噪聲壓制[16]一直是困擾廣大電磁工作者的關鍵問題，在一定程度上制約了電磁技術方法的發展。WFEM采用人工場源，信號發射機功率很高(最高可達250kW)，相較于天然場源電磁法，采集到的電磁信號信噪比明顯提高。在實際勘探時，WFEM同樣受到噪聲干擾的問題，去噪仍然是數據處理的關鍵環節。此外，為了提高WFEM方法的縱向分辨率，須加強WFEM數據去噪方法的研究。因此，如何利用新方法對WFEM數據進行有效去噪是目前亟需解決的關鍵問題之一。

針對人工源電磁信號受噪聲影響的問題，張必明等[17]提出了自適應雙向均方差閾值法，先對原始數據樣本進行排序，采用迭代或遞歸的方式，以中點為界計算前、后兩部分數據的均方差；再利用單一的閾值控制異常點，在少量頻點信息被噪聲淹沒時，可有效剔除頻譜數據中的粗大誤差。但是，僅采用閾值控制可能會造成有效數據的過處理及部分噪聲不明顯的數據未被處理，同時該方法也沒有考慮到整體曲線的圓潤性。針對此問題，Mo等[18]提出了一種灰色系統理論和穩健估計方法，通過灰色建模求解測量數據的標準差，結合閾值法識別、剔除頻譜數據中的異常值，將最終保留數據的M估計值視為真實值。陳超健等[19]提出了一種灰色判別準則和基于有理函數濾波的電磁數據去噪方法，通過灰色判別準則剔除各個頻點頻譜數據的明顯異常值，再利用有理函數濾波充分平滑數據曲線，進而約束數據的二次優化處理，以剔除殘余噪聲。

上述頻率域處理方法主要針對干擾小的電磁數據，可顯著提高數據質量。若電磁采集工作受到持續性強人文噪聲的影響，大部分頻點會發生畸變，依賴于頻率域功率譜挑選的去噪方法就會失效。為此，楊洋等[20]提出頻率域基于小波變換和Hilbert解析包絡的CSEM噪聲評價方法，在不增加野外工作量的同時可提取更豐富的頻率信息，針對高階諧波信號，實現了頻率加密，提高了CSEM數據的縱向分辨率。針對時間域CSEM數據處理，Yang等[21]利用不含噪聲或僅含隨機噪聲的時序片段，求解出噪聲頻譜，獲得了高質量的數據，但該方法的前提是選出一段高質量的時間序列片段，而精準地挑選出高質量數據的位置較難。Li等[22]提出了基于快速傅里葉變換、互補總體經驗模態分解以及與移不變稀疏編碼結合的CSEM去噪方法，可精準去除工頻干擾，校正基線漂移，并基于稀疏表示和字典學習剔除信號中的強干擾，重構獲取高質量的偽隨機有用信號，獲得的視電阻率曲線和電場曲線更穩定、連續。上述時間域廣域電磁信號處理方法均是對采集的電磁數據進行整體處理，雖在一定程度上提高了數據質量，但缺乏對信號和噪聲的辨識環節[23-24]。若能從CSEM信號和噪聲類型的特征入手，剖析CSEM信號與噪聲之間的定性及定量辨識關系，這會對后續電磁法高精度的信噪分離起到關鍵作用，并能避免時間域處理技術在噪聲壓制過程中的過處理，最大限度地保留CSEM有效信號。

模糊C均值(Fuzzy C-means, FCM)聚類屬于無監督機器學習方法，主要通過優化目標函數計算每個樣本點對所有樣本類中心的隸屬度，從而決定樣本點的類屬以自動進行數據分類，主要應用于圖像分割、地震屬性分析等[25-26]。本文主要以廣域電磁法7頻波2頻組(7-2頻組)和7頻波3頻組(7-3頻組)的數據為例進行分析與處理，其中7-2頻組的頻率為1、2、4、8、16、32、64Hz，7-3頻組的頻率為0.75、1.5、3、6、12、24、48Hz。由于廣域電磁法偽隨機信號屬于周期信號，具備一定的特征及規律性，因而根據時域統計特征能對廣域電磁7頻波數據進行信號與噪聲的定量識別。本文通過對比電場曲線和視電阻率曲線，分析、評價方法的處理效果。通過提取時域統計特征并結合FCM聚類方法進行廣域電磁偽隨機信號處理，提出了一種特征分析與聚類識別相結合的人工源電磁偽隨機信號處理方法。首先，建立偽隨機信號和典型噪聲的樣本庫，分析樣本庫信號的時、頻域特征；然后，利用時域統計特征結合FCM聚類進行信、噪識別，重構有用信號，并利用數字相干技術[27]提取有效頻點的頻譜；最后，通過本文方法對合成數據及實測數據進行處理，驗證方法的可靠性與有效性。

1 偽隨機信號的時域特征

分析采集的廣域電磁7頻波時域信號可知，時間域序列中的偽隨機信號是有用信號[3]，但時域信號難免會受到各種異常波形的影響，導致信號發生跳變或衰減，從而在頻譜上出現頻點信息的淹沒或丟失。因此，識別異常波形并壓制其影響能有效提升數據質量。

本文通過提取簡單的時域統計特征參數進行信、噪的定量識別，通過探討量綱特征參數和無量綱特征參數識別廣域電磁信號與噪聲。信號中的無量綱特征參數主要包括峰值因子、脈沖因子、裕度因子、波形因子、峭度和偏度，量綱特征參數主要有最大值、最小值、峰峰值、平均值、方根幅值、方差、標準差和有效值等。本文主要研究最大值、峰峰值、峰值因子、脈沖因子和裕度因子等特征參數，這些參數有助于快速、直觀地表征廣域電磁信號與噪聲之間的關系，為后續聚類算法提供有效的特征參數。本文重點介紹無量綱特征參數。

1.1 無量綱特征參數

1.1.1 峰值因子

峰值因子FFZ是信號峰值Xp與有效值Xrms的比值，表征峰值在波形中的極端程度。其計算公式為

(1)

1.1.2 脈沖因子

脈沖因子FMC是檢測信號中是否存在沖擊的統計指標，其計算公式為

(2)

由式(1)和式(2)可知，FMC和FFZ的區別在分母，由于同一組數據的絕對值的平均值小于有效值，因而FMC必定大于FFZ，根據這兩個特征參數可檢測信號中的突變成分和沖擊因素。

1.1.3 裕度因子

裕度因子FYD的物理意義與峰值因子和脈沖因子相似，常用于檢測設備與信號中的損傷和變化程度。其定義為

(3)

式中Xr表示信號幅值的方根。

1.2 樣本庫特征分析

為了剖析實測廣域電磁7頻波數據中的偽隨機信號與異常干擾波形之間的定量辨識關系，基于實測WFEM數據的信號和噪聲特點，利用Matlab構建一個由模擬典型干擾和模擬偽隨機信號形成的樣本庫，如圖1所示。從樣本庫中任意選取三類信號的時域波形及其頻譜，其中樣本庫包含了30個偽隨機信號、30個脈沖干擾信號和30個衰減干擾信號，每個樣本信號的采樣長度均為1200，采樣率為300Hz。

通過觀測噪聲樣本庫(圖1)可知，時域脈沖干擾和時域衰減干擾導致原始偽隨機信號(圖1a)發生異常突變(圖1b、圖1c)，幅值增大，相應的頻譜及主頻信息也出現不同程度的失真。由圖1a可見，未受干擾的偽隨機信號呈周期性，幅值穩定，頻譜也相對穩定，其頻點信息被完整地保留。分析圖1b、圖1c可見，若時域信號中存在著干擾或異常波形，信號在時、頻域出現嚴重畸變，無法反映原始偽隨機7頻波信號的固有特征。

圖1 樣本庫信號(上)及頻譜(下)

圖2為樣本庫信號的特征參數統計結果。分析圖2可知，含干擾信號的特征參數明顯大于偽隨機信號的特征參數，偽隨機信號的特征參數值穩定且幅值小?；谶@些時域統計特征參數(最大值、峰峰值、FFZ、FMC和FYD)，可較好地識別偽隨機信號和含干擾的信號。因此，結合這些特征參數可提高后續聚類分析的信噪區分度和信噪識別效果。

圖2 圖1樣本庫信號的特征參數統計結果

2 FCM聚類分析

FCM聚類是一種典型的無監督機器學習方法[28]，該方法按照“物以類聚，人以群分”的思想將樣本點按照某種規律進行劃分，這些規律是通過樣本點的某些特征確定的，無需事先給定或約束。FCM聚類算法通過定義下面的最小化目標函數得到最優解[29]

(4)

式中：m表示聚類的簇數；N表示樣本數；C表示聚類中心數；uij表示第j個待分類對象屬于第i個聚類中心的隸屬度；xi表示第i個樣本；cj表示第j個聚類中心；‖·‖表示數據相似性的度量。

通過選用0～1的隨機數確定初始化隸屬度矩陣UC×N，使其元素uij滿足以下約束條件

(5)

針對上述約束問題，利用拉格朗日乘數法對目標函數進行求導，可得

(6)

計算每組的聚類中心Vi

(7)

使得目標函數最小，即歐式距離最短，相似度最高。迭代更新uij和Vi，直到前后兩次隸屬度的最大變化值低于預設的誤差閾值，即目標函數達到最小值。

圖3所示為對圖1樣本庫中的90個樣本提取的峰峰值、峰值因子、脈沖因子和裕度因子分別進行K均值聚類和FCM聚類的結果。

圖3 圖1樣本庫信號的K均值聚類(a)和FCM聚類(b)分析

分析圖3可知：雖然K均值聚類分析能將樣本庫中的90個樣本分為兩類，但對于海量的電磁數據特征值較亂時，K均值聚類無法選取合適的聚類中心，導致聚類效果不夠理想；FCM聚類分析通過計算每個樣本點到聚類中心的歐式距離，自動將兩類含干擾的信號與偽隨機信號劃分為不同的類別，聚類效果明顯。因此，FCM聚類有利于后續的干擾信號剔除，保留偽隨機信號。

3 仿真實驗分析

為了驗證方法的識別效果，對模擬樣本庫信號中的典型干擾類型進行仿真實驗。圖4所示為模擬圖1所示樣本庫中的7-2頻組數據的信、噪識別結果及其頻譜，以及作為對比的小波去噪結果及其頻譜。

從圖4可知，含噪信號(圖4a)經時域特征提取，根據聚類分析結果能準確地識別信號中的異常干擾，保留偽隨機信號(圖4b)。結合頻譜分析(圖4b)可知，主頻信息受干擾影響，主頻值不穩定，且異常程度不同。對比時頻分析小波去噪結果(圖4c)可知，針對時域中的噪聲進行去噪處理的同時，偽隨機信號也被濾除，過處理導致無法還原偽隨機信號的形態；同時，對應的主頻信息因殘余少量噪聲導致主頻幅值降低。經本文方法處理后(圖4d)，典型的干擾信號被有效識別，并被準確去除，噪聲頻譜也同時被濾除，偽隨機信號的有效波形及其頻譜原始特征被有效重構、還原。

圖4 模擬7-2頻組信號的信噪識別(左)及頻譜分析(右)

為了定量分析本文方法的有效性，統計了模擬7-2頻組信號經不同方法去噪處理后的電場幅值及相對誤差(表1)。

分析表1可知，由于噪聲添加在不同的時刻，導致時間域和頻率域信號出現突變或混亂現象，且含噪信號在1～16Hz的5個頻點的電場幅值均大于真實值，最大達1.301mV，對應頻率為16Hz，相對誤差為34.95%。因而，噪聲在時域會嚴重影響偽隨機信號的特征，在頻域會影響不同的主頻信號及其諧波。根據小波去噪方法處理結果(圖4c)可見，7-2頻點的主頻幅值整體更低，相對誤差平均達14.10%，小波去噪方法雖能壓制異常波形，但嚴重過處理導致了該方法失效。因而，經本文方法處理后的電場幅值更接近于真實電場值，相對誤差也較低，最大誤差僅0.52%。

表1 模擬7-2頻組信號經不同方法處理后的電場幅值及相對誤差

4 實測數據分析

4.1 實測信/噪識別與處理

實測數據來自四川盆地某頁巖氣廣域電磁法勘探項目。以某一段含典型干擾的實測CSEM信號(采樣率1200Hz)為例(圖5a)，采用本文方法進行特征提取與信噪識別，結果見圖5b和圖5c。

分析圖5可知，實測信號中包含了脈沖干擾和衰減干擾，通過本文方法提取時域統計特征參數，結合FCM聚類，處理結果有效地識別了噪聲(圖5b)和偽隨機信號(圖5c)得到有效識別并得以保留。

圖5 實測電場時域信號及本文方法識別結果

CSEM的發射信號和接收信號都是偽隨機信號，若接收到的時域信號為明顯的偽隨機信號且無任何異常波形的干擾，則頻譜也會穩定，無主頻點出現紊亂現象。對實測的無強干擾波形的7-2頻組信號(圖6a、圖6b)，對其添加脈沖干擾和衰減噪聲，得到加噪數據(圖6c、圖6d)，經小波去噪和本文方法處理，提取的噪聲及其頻譜見圖6e和圖6f，重構數據及頻譜見圖6g和圖6h，經小波去噪后的數據及頻譜見圖6i和圖6j。對應的廣域電磁法Ex分量頻域電場曲線見圖7。

從圖6a的實測數據(詳見前4s的放大圖)可知，數據整體上屬于周期性偽隨機信號，電場曲線相對平穩且無異常跳變。因此，該數據不存在強干擾的影響，數據質量較高[20]。通過在此信號中添加典型干擾，信號發生突變，其頻譜也相應地發生變化，體現在電場曲線上主頻點的跳變及頻譜混亂(圖6b和圖6c)。分析小波方法處理結果(圖6i、圖6j)可知，由于小波基函數、分解層數和閾值較難選取，小波去噪效果也不甚理想。由于時間域信號過處理現象和大量有用信息的丟失，頻率域中噪聲殘余使得電場曲線不夠穩定、不光滑。經過本文方法處理后(圖6g和圖6h)，基于電場幅值和頻譜信息可有效識別干擾信號段并將其剔除，重構信號基本消除了人為添加的典型干擾，成功地還原出原始偽隨機信號。進一步地采用數字相干提取技術(圖7)，對比加噪前、后的WFEM電場曲線可知，噪聲干擾是主頻點波動的主要原因，將其從數據中有效識別并剔除后，利用去噪后的數據可以基本恢復原始WFEM電場曲線的趨勢及形態特征。

圖6 加噪實測信號及處理結果

圖7 圖6數據對應的WFEM電場幅值對比

4.2 實測點分析

廣域電磁7頻波數據的頻段范圍為0.0117～

8192Hz，實測點僅涉及7-2、7-3頻組數據，因而本文分析的電場曲線和視電阻率曲線僅涉及這兩個頻組信號。這里以測點S1和S2為例，實測數據經本文方法處理后的結果見圖8和圖9。

分析圖8、圖9可知，測點S1和測點S2觀測電場曲線和視電阻率曲線在頻率低于24Hz時均呈鋸

圖8 測點S1本文方法處理結果

圖9 測點S2本文方法處理結果

齒狀，發生跳變，說明觀測數據的時間域波形中包含了脈沖干擾和衰減干擾；采用本文方法提取特征參數并進行FCM聚類，處理后的數據可有效識別噪聲干擾，并進行剔除，利用數字相干技術提取的歸一化電場曲線和計算得出的視電阻率曲線平穩、連續，可見本文方法基本消除了噪聲干擾引起的頻點跳變現象，重構了高質量的偽隨機信號。

5 結論

針對CSEM數據受噪聲干擾的問題，提出了基于特征提取和聚類識別的廣域電磁數據偽隨機信號處理方法。基于仿真實驗及實測數據處理結果得到如下結論。

(1)所提方法針對時間域數據進行信、噪辨識處理，能快速、準確地識別信號中的干擾或異常波形，重構并還原出偽隨機信號，提升了數據質量，彌補了傳統CSEM數據僅在頻率域內開展去噪處理的短板。

(2)模擬實驗和實測數據處理結果表明，經本文方法處理得到了穩定、連續的歸一化電場曲線和視電阻率曲線，能有效應用于電磁數據等周期信號的特征提取及信號與噪聲辨識處理，為廣域電磁法偽隨機信號數據處理提供了技術支持。

本文方法簡單且易于操作，能快速而有效地識別噪聲，可提高廣域電磁數據質量，為后續反演解釋提供可靠數據。但本文方法僅在時間域進行信號與噪聲識別處理，未分析如何去噪，因此，如何去除識別的噪聲干擾是后續的研究重點。