自適應奇異值分解局放信號降噪方法

孫傳銘， 魏 隆， 張夢楠， 劉 凱， 潘貴翔， 高國強

(1. 中車青島四方機車車輛股份有限公司， 山東 青島 266111；2. 西南交通大學電氣工程學院， 四川 成都 610031)

1 引言

高壓電氣設備局部放電(Partial Discharge，PD)檢測過程中往往會受到各種背景噪聲的干擾，尤其在工程現(xiàn)場的電磁干擾對檢測結(jié)果的影響更加嚴重，有時甚至會出現(xiàn)局放信號完全被背景噪聲湮沒的情況，對后續(xù)電氣設備絕緣狀態(tài)判斷及檢修帶來一定的困難[1]。局放測試中的干擾主要分為以下三種[2-5]：隨機性脈沖干擾、周期性窄帶干擾和白噪聲干擾。其中，隨機脈沖干擾通常強度大、頻率低，易于識別和濾除。周期性窄帶干擾主要來源于電網(wǎng)絡內(nèi)部及環(huán)境中的無線電廣播等信號，其出現(xiàn)頻率高，幅值大，常在時域中將局放信號湮沒，且在頻域范圍內(nèi)經(jīng)常與局放信號發(fā)生混疊，對局放信號的檢測影響很大。白噪聲干擾主要是電氣設備的熱噪聲引起的寬帶干擾隨機信號，在頻域上與局放信號具有相似特征。因此，如何有效濾除周期性窄帶干擾和白噪聲成為局放信號研究的一大難點。

針對周期性窄帶干擾和白噪聲干擾混合噪聲的抑制方法，國內(nèi)外學者進行了大量的研究。文獻[6]提出結(jié)合廣義S時頻變換和奇異值分解去噪方法抑制信號中的混合噪聲，該方法可有效地抑制混合噪聲，但去除窄帶干擾信號時需要人為判斷窄帶干擾區(qū)域，受人為因素影響存在一定誤差，不具備自適應性。文獻[7]針對奇異值分解耗時較長且有效奇異值數(shù)量難以選擇問題，提出基于滑動短時數(shù)據(jù)能量窗的奇異值分解降噪方法，該方法無需預先假設信號中含有周期性窄帶干擾，可自動實現(xiàn)周期性窄帶干擾的甄別和混合噪聲的抑制。但該方法在去除窄帶干擾時往往受計算精度的影響而存在誤差，最終影響降噪結(jié)果。文獻[8]提出基于經(jīng)驗模態(tài)分解(Empirical Mode Decomposition，EMD)和獨立成分分析(Independent Component Analysis，ICA)的高壓電纜局放信號降噪方法，利用快速ICA方法進一步濾除含噪IMF分量中的噪聲信號，但該方法并未明確提及是否具備同時濾除兩種噪聲的能力。文獻[9]提出一種基于總體經(jīng)驗模態(tài)分解和補充總體經(jīng)驗模態(tài)分解的局部放電閾值去噪新方法，該方法通過對總體經(jīng)驗模態(tài)分解的IMF分量進一步通過補充總體經(jīng)驗模態(tài)分解提高降噪能力，但該方法耗時較長，應用受限。文獻[10]提出基于經(jīng)驗小波和小波變換的局放信號降噪方法，通過兩種方法結(jié)合實現(xiàn)降噪優(yōu)化，但對于小波變換依舊存在基函數(shù)和分解層數(shù)的選擇問題，自適應性能較差。

針對上述局放信號混合噪聲干擾抑制存在的問題，本文提出一種自適應奇異值分解降噪方法。該方法首先對測試信號構建Hankel矩陣，以此作為軌跡矩陣進行奇異值分解[11]。通過提取前兩個奇異值進行重構并結(jié)合功率譜熵自適應判斷染噪信號中是否存在窄帶干擾；隨后確定奇異值系列特征量，結(jié)合K類均值聚類[12]對窄帶干擾、有效PD信號和白噪聲所對應的奇異值進行劃分，對有效PD信號對應的奇異值進行重構進而還原PD信號。該方法可自適應地判斷是否存在窄帶干擾，從而決定聚類區(qū)間；通過1次K類均值聚類分類即可獲取有效奇異值數(shù)據(jù)，自適應性能良好。

2 奇異值分解

2.1 奇異值分解原理

2.1.1 軌跡矩陣的構建

本文選取Hankel矩陣作為奇異值分解的軌跡矩陣。Hankel矩陣具體構建方式如下：

設染噪信號X為：

X=[x(1),x(2),…，x(N)]

(1)

對采樣系列X構造Hankel矩陣：

(2)

式中，N=m+n-1，本文中n取為N/2。

2.1.2 奇異值獲取

矩陣A是一個m×n的矩陣，其秩為r，則必存在m×m的正交矩陣U和n×n的正交矩陣V，使得：

A=UΛVT

(3)

其中

(4)

對角矩陣Λ除了前r階對角元素外，其他元素均為零。對角元素?i即為矩陣A的奇異值，且數(shù)值由大到小排列，奇異值的大小反映了能量的集中情況。通過對窄帶干擾、有效PD信號和白噪聲對應的奇異值規(guī)律進行剖析，進而選取合適的奇異值進行重構，即可還原真實的局放信號。

2.2 信號仿真

結(jié)合以往數(shù)據(jù)仿真經(jīng)驗，本文選取單指數(shù)振蕩衰減模型和雙指數(shù)振蕩衰減模型來模擬理想局放信號[13]。具體表達式如下：

(5)

(6)

式中，f為振蕩頻率；τ為衰減系數(shù)；A1、A2為脈沖幅值。

本文仿真了四種局部放電脈沖，其中脈沖模型1和模型2根據(jù)式(5)得出，脈沖模型3和模型4根據(jù)式(6)得出。脈沖仿真參數(shù)見表1。

表1 局放仿真信號參數(shù)Tab.1 PD simulation signal parameters

周期性窄帶干擾通常呈正弦或余弦波形[14]，且PD信號實際檢測中往往會存在窄帶干擾與PD信號混疊的問題。故本文選取窄帶干擾頻率分別0.5 MHz、1 MHz、2 MHz、5 MHz、7 MHz。周期性窄帶干擾的具體表達式如下：

(7)

式中，Ai對應各窄帶干擾信號幅值；fi為頻率。模擬窄帶干擾仿真信號參數(shù)見表2。

表2 窄帶干擾仿真信號參數(shù)Tab.2 Parameters of narrow-band interference simulation signals

實際運行環(huán)境中除受窄帶干擾影響外，往往還會受到白噪聲的干擾，白噪聲利用高斯白噪聲模擬產(chǎn)生。仿真獲取理想PD仿真信號如圖1(a)所示，添加周期性窄帶干擾和白噪聲后的信號如圖1(b)所示，圖1(c)為染噪信號頻域譜圖。

圖1 仿真波形Fig.1 Simulation waveforms

2.3 窄帶干擾奇異值特征分析

針對窄帶干擾與奇異值分解間存在的關系問題，文獻[15]發(fā)現(xiàn)對含有單個頻率窄帶干擾的PD信號進行SVD分解時，提取前兩個奇異值可有效提取窄帶干擾信號。文獻[16]發(fā)現(xiàn)每個頻率的窄帶干擾都對應兩個非0奇異值，通過提取窄帶干擾頻率個數(shù)n對應的前2n個奇異值即可提取窄帶干擾并濾除。文獻[17]進一步指出隨著采樣數(shù)據(jù)長度增加，窄帶干擾對應奇異值幅值越來越大，而局放信號對應奇異值變化較小。因此，通過增加數(shù)據(jù)長度，可保證窄帶干擾被全部濾。染噪信號奇異值隨數(shù)據(jù)長度變化情況如圖2所示。

圖2 奇異值與數(shù)據(jù)長度的關系Fig.2 Relation between singular values and data length

3 局放混合噪聲抑制方法

3.1 窄帶干擾判別

為實現(xiàn)對局放信號混合噪聲的自適應抑制，首先需要對染噪信號中是否存在窄帶干擾進行判別。根據(jù)前人對窄帶干擾與PD信號奇異值規(guī)律的剖析，本文提取信號前兩個奇異值進行重構，根據(jù)重構后的信號是否符合正(余)弦規(guī)律即可判斷是否存在窄帶干擾。對于功率譜熵而言，信號混亂程度越高，其功率譜熵越大，混亂程度越低，功率譜熵越小。正弦信號混亂程度較局放信號較小，因此，本文引入功率譜熵[18,19]的概念對正(余)弦信號進行檢測。利用正弦信號與局放信號間功率譜熵大小的差異，判斷是否存在窄帶干擾?；诠β首V熵檢測的具體步驟為：

(1)將信號x(t)經(jīng)FFT變換得到功率譜為：

(8)

(9)

式中，N為數(shù)據(jù)點個數(shù)。

(2)求取信號功率譜熵H為：

(10)

(11)

(3)H作為檢測統(tǒng)計量為：

(12)

當檢測統(tǒng)計量H小于檢測閾值T時，即可確定重構信號是窄帶干擾。

圖3為局放信號和窄帶干擾信號幅值和頻率改變時分別對應的功率譜熵幅值。從圖3中可以看出，無論窄帶干擾信號幅值和頻率如何變化，對應功率譜熵幅值均小于1；而局放信號功率譜熵幅值始終大于1。經(jīng)過多組數(shù)據(jù)分析，本文最終設定檢測閾值T=1，若重構信號對應檢測閾值T<1，判定染噪信號中含有窄帶干擾。

圖3 功率譜熵求取結(jié)果Fig.3 Result of power spectrum entropy obtained

3.2 有效奇異值選取

奇異值有效個數(shù)選取問題一直是奇異值分解降噪的關鍵[20]。如果有效奇異值個數(shù)選取過少，將會損失局放信號部分有用信息；如果有效奇異值個數(shù)選取過多，則降噪效果不明顯。此外，如何實現(xiàn)有效奇異值個數(shù)的自適應選取也是一大研究重點。根據(jù)以往的研究得知，窄帶干擾信號奇異值遠大于局放混合白噪聲信號對應奇異值，且該數(shù)值位于奇異值分解前列。同時，局放信號相對白噪聲而言，其奇異值相對較大且數(shù)據(jù)較分散。因此，本文引入K類均值聚類算法[21]對局放混合噪聲進行抑制。為凸顯奇異值大小，本文選取奇異值本身F1作為奇異值系列特征參量，同時，為表征數(shù)據(jù)離散程度，同時選取奇異值子集標準差F2作為另一奇異值系列特征量，并與文獻[17]中選取以奇異熵增量F3及其能量F4作為奇異值系列特征量進行比較，對混合染噪信號進行分類，結(jié)果如圖4和圖5所示。各類特征參量計算公式及表達式如下。

圖4 F1、F2特征參量K類均值聚類結(jié)果Fig.4 Results of K-means clustering of characteristic parameters F1 and F2

圖5 F3、F4特征參量K類均值聚類結(jié)果Fig.5 Results of K-means clustering of characteristic parameters F3 and F4

奇異值子集ssi構建及子集標準偏差σi計算公式為：

ssi=[?r,…,?i]

(13)

(14)

式中，?i對應第i個奇異值；μi為ssi數(shù)據(jù)均值。

奇異熵增量計算公式：

(15)

式中，r為奇異值總數(shù)。

奇異熵增量能量計算公式：

ei2=(Δei)2

(16)

特征量F1、F2、F3、F4表達式如下：

F1=[?r,?r-1,…,?1]

(17)

F2=[σ1,σ2,…,σr]

(18)

F3=[Δe1,Δe2,…,Δer]

(19)

F4=[Δe12,Δe22,…,Δer2]

(20)

從圖4中可以看出，本文選取的奇異值特征量實現(xiàn)了窄帶干擾、有效PD信號和白噪聲信號對應奇異值的有效分類，重構信號在保證濾除混合噪聲的同時保留了局放信號的完整性，本文方法更適用于同時實現(xiàn)3種不同信號的有效分離。此外，在局放信號中混疊窄帶干擾的情況下，本文分別比較了2次K類均值聚類分類次數(shù)為2的奇異值分類和1次K類均值聚類分類次數(shù)為3的分類。結(jié)果表明，采用本文方法選取的特征參量進行兩種分類算法獲取的有效奇異值基本吻合。為節(jié)省計算時間，本文最終選取僅作1次K類均值聚類分類次數(shù)為3的計算。

3.3 混合噪聲抑制步驟

要實現(xiàn)混合噪聲的自適應抑制，窄帶干擾信號的判別至關重要。本文首先對信號進行奇異值分解，通過提取前兩個奇異值重構判斷信號中是否存在窄帶干擾；然后采用K類均值聚類最終提取有效奇異值，進而獲取局放信號?；旌显肼曇种撇襟E具體如下：

(1)對染噪信號構建Hankel矩陣，以此作為軌跡矩陣進行奇異值分解；

(2)重構前兩個奇異值獲取重構信號，判斷信號對應功率譜熵T是否大于1。若T<1，則存在窄帶干擾，設置K類均值聚類分類個數(shù)n=3；若T>1，則不存在窄帶干擾，設置n=2；

(3)分別以奇異值本身F1和奇異值子集標準差F2為奇異值系列特征量，通過K類均值聚類將奇異值系列分為n類；

(4)若n=3，則選取第二類奇異值數(shù)據(jù)進行重構還原PD信號；若n=2，則選取第一類奇異值數(shù)據(jù)進行重構還原PD信號。

綜上所述，本文還原PD信號的流程圖如圖6所示。

圖6 本文降噪方法流程圖Fig.6 Flow chart of noise reduction method in this paper

4 去噪效果對比

為分析本文的自適應奇異值分解降噪方法對局放信號的降噪效果，對原始PD仿真信號加入周期性窄帶干擾并疊加分布為(0，10)的高斯白噪聲進行降噪處理。通過引入FFT-小波變換降噪、S時頻變換-EEMD聯(lián)合去噪方法與本文方法進行對比。各方法降噪結(jié)果如圖7所示。經(jīng)對比得出，F(xiàn)FT-小波變換降噪可有效濾除噪聲，但局放信號也被部分濾除，導致信號減??；S變換-EEMD聯(lián)合降噪既不能保證噪聲的高精度濾除，同時PD信號存在部分衰減。本文所選方法能同時滿足噪聲信號的高精度濾除和PD信號的高度還原。

圖7 3種方法降噪結(jié)果對比Fig.7 Denoising results of three method

本文引入去噪評價參數(shù)信噪比(Signal to Noise Ratio，SNR)、均方誤差(Mean Square Error，MSE)和波形相似參數(shù)(Normalized Correlation Coefficient，NCC)[22]進一步對降噪效果進行評估。

去噪評價參數(shù)計算結(jié)果如表3所示。從表3中可以看出，本文所采用的方法具有明顯的優(yōu)勢，無論從信噪比、均方誤差還是波形相似參數(shù)上都顯示出非常好的效果，對于混合噪聲干擾的抑制效果最好，且信號還原度最高。

表3 去噪評價參數(shù)計算結(jié)果Tab.3 Calculation results of denoising evaluation parameters

5 實測信號去噪分析

為檢驗本文方法對于實測PD信號濾除混合噪聲的能力。基于實驗室條件下搭建電纜終端刀痕缺陷測試模型如圖8所示，測試采用的高頻脈沖電流傳感器-6 dB帶寬為80 kHz～40 MHz，采樣率為50 MSa/s。測試得局放波形如圖9(a)所示。實測PD信號基于理想試驗條件下測得，故而PD信號明顯，而環(huán)境噪聲干擾很小。為測試本文降噪方法對實測信號的去噪效果，通過對實測PD信號施加3個幅值為0.5 mV，頻率分別為0.5 MHz、2 MHz和8 MHz的周期性窄帶干擾信號，并疊加分布為(0，10)的高斯白噪聲。染噪信號如圖9(b)所示。

圖8 局放檢測平臺原理Fig.8 Schematic of PD detection circuit

圖9 實驗室實測PD信號Fig.9 Measured PD signals in laboratory

分別采用FFT-小波變換降噪、S時頻變換-EEMD聯(lián)合去噪方法和本文方法對添加周期性窄帶干擾和白噪聲的實驗室環(huán)境下實測PD信號進行降噪處理，各方法降噪結(jié)果如圖10所示。從圖10中可以很明顯地看出，本文降噪方法能夠高度還原實測PD信號，且抑制噪聲效果最好。FFT-小波變換降噪明顯改變了PD信號的特征。同時，PD信號明顯減??；S時頻變換-EEMD聯(lián)合降噪方法較FFT-小波變換降噪效果更佳，但同樣存在信號衰減問題。此外，S時頻變換-EEMD降噪去除噪聲效果相對較差。

圖10 實驗室實測PD信號降噪結(jié)果Fig.10 Noise reduction results of PD signal measured in laboratory

由于無法測得完全不含噪聲的PD信號，故無法使用上述去噪評價參數(shù)對各方法降噪效果進行定量分析。因此，本文引入噪聲抑制比μ1和幅值衰減比μ2對降噪效果進行評價[11]。其中，μ1反映了降噪后信號的凸顯程度。μ1越大，說明降噪方法去噪效果越好。μ2反映了降噪前后PD信號的衰減程度。μ2越大，說明降噪后PD信號衰減越嚴重。μ1、μ2具體定義如式(21)、式(22)所示。

μ1=10(lgδ12-lgδ22)

(21)

(22)

式中，δ1、δ2分別為降噪前后信號的標準偏差;Am1、Am2分別為降噪前后信號的最大幅值。各降噪方法降噪評價參數(shù)計算結(jié)果如表4所示。

從表4中可以看出，文本降噪方法降噪效果最好，信號衰減程度最小。故而選取本文的方法在降噪上占有很大的優(yōu)勢。

表4 降噪評價參數(shù)計算結(jié)果Fig.4 Results of evaluation parameters of noise reduction

6 結(jié)論

本文基于奇異值分解自適應降噪，提出了一種有效濾除白噪聲和周期性窄帶干擾的降噪方法。通過與FFT-小波變換降噪和S時頻變換-EEMD聯(lián)合降噪方法進行對比，分析結(jié)果發(fā)現(xiàn)本文方法具有更優(yōu)的降噪效果，且PD信號衰減最小，還原度最高。具體結(jié)論如下：

(1)信號降噪前首先對染噪信號提取前兩個奇異值重構，自適應判斷信號中是否含有窄帶干擾，避免了因窄帶干擾存在與否問題導致的降噪失誤，從而為后續(xù)利用SVD實現(xiàn)降噪奠定基礎。

(2)利用Hankel矩陣作為軌跡矩陣進行奇異值分解時，窄帶干擾頻率個數(shù)對應兩倍奇異值個數(shù)，且采樣數(shù)據(jù)足夠長時，窄帶干擾對應奇異值數(shù)值遠大于PD信號。

(3)采用奇異值本身和奇異值子集標準差作為奇異值系列特征量時，利用K類均值聚類分類方法可一次實現(xiàn)對窄帶干擾、有效PD信號和白噪聲對應奇異值的有效分離，且分類效果與分別進行兩次分類的效果一致，分類所需用時有效縮短。

(4)根據(jù)窄帶干擾判斷結(jié)果可自適應確定K類均值分類個數(shù)并進行分類，進而對有效奇異值進行重構得到PD信號，自適應性能良好。

(5)通過與FFT-小波變換降噪和S時頻變換-EEMD聯(lián)合降噪方法進行對比，發(fā)現(xiàn)本文方法能更好地抑制噪聲，同時保證PD信號的高度還原。