余成義，肖 涵，劉安中，李友榮

(武漢科技大學(xué) 冶金裝備及其控制教育部重點(diǎn)實(shí)驗(yàn)室，武漢 430081)

混沌動(dòng)力系統(tǒng)的降噪是20世紀(jì)90年代混沌時(shí)間序列分析的主要研究?jī)?nèi)容之一。由于混沌信號(hào)具有“寬頻”性和“偽隨機(jī)”性，如果用傳統(tǒng)的基于頻譜分析的降噪算法，會(huì)將有用信號(hào)當(dāng)成噪聲去掉。因此在處理混沌信號(hào)時(shí)，必須用非線性的降噪算法取代傳統(tǒng)的線性降噪算法。

目前，Shin等［1－2］提出了基于奇異譜的降噪算法，并將該算法用于軸承、齒輪信號(hào)的降噪中，克服了傳統(tǒng)線性降噪方法無(wú)法處理非線性數(shù)據(jù)的缺點(diǎn)。徐金梧等［3］將局部投影降噪用于微弱特征信號(hào)的檢測(cè)中，提取設(shè)備的早期故障信號(hào)。文獻(xiàn)［4］提出了加權(quán)相空間重構(gòu)投影降噪算法，該算法采用漢寧加權(quán)窗，減少了重構(gòu)誤差，提高算法的效率。但上述算法主要通過(guò)主分量分析去獲取數(shù)據(jù)集內(nèi)總方差最大方向，與主分量分析(PCA)相比，獨(dú)立分量分析(ICA)應(yīng)用了數(shù)據(jù)高階統(tǒng)計(jì)量，且ICA能恢復(fù)數(shù)據(jù)集內(nèi)在變量的各自主要方向。為了克服全局線性ICA不適合處理非線性數(shù)據(jù)的局限性和非線性ICA需要不同約束條件的復(fù)雜性，Karhunen等［5－6］提出了局部線性ICA方法，其基本思想是先對(duì)數(shù)據(jù)進(jìn)行聚類，聚類目的就是找數(shù)據(jù)的局部投影區(qū)間，再對(duì)每一個(gè)聚類進(jìn)行線性 ICA。Gruber等［7］結(jié)合相空間重構(gòu)給出了局部獨(dú)立分量分析(本文中稱其為傳統(tǒng)局部獨(dú)立分量分析)的基本模型框架。文獻(xiàn)［8］指出該降噪算法的一些缺點(diǎn):聚類中心個(gè)數(shù)人為給定和初始聚類中心隨機(jī)選取導(dǎo)致算法的不穩(wěn)定性，從而大大限制了該降噪算法的使用范圍。本文在克服傳統(tǒng)局部獨(dú)立分量分析的不足的基礎(chǔ)上提出了自適應(yīng)獨(dú)立分量分析降噪算法，該算法用聚集模糊K均值聚類和聚類評(píng)價(jià)函數(shù)求取高維數(shù)據(jù)集的聚類個(gè)數(shù)和聚類中心位置，具有自適應(yīng)性和穩(wěn)定性。

1 聚集模糊K均值聚類算法

聚集模糊K均值聚類算法［9］在傳統(tǒng)模糊K均值聚類算法的基礎(chǔ)上引入懲罰項(xiàng)，使得該算法得到的聚類結(jié)果對(duì)初始聚類中心不敏感。該算法和聚類評(píng)價(jià)函數(shù)一起使用能夠確定數(shù)據(jù)集的聚類個(gè)數(shù)和聚類中心位置。

某一維信號(hào)經(jīng)相空間重構(gòu)后為:

X={X1，X2，…，Xn}T

其中相點(diǎn) Xi=［Xi，1，Xi，2，…，Xi，m］，m 為嵌入維數(shù)或代表每個(gè)相點(diǎn)屬性個(gè)數(shù)。用聚集模糊K均值聚類算法將數(shù)據(jù)集X聚成k類，實(shí)質(zhì)求下面目標(biāo)函數(shù)的最小值:

其中λ為懲罰因子，其大小決定式(1)中哪項(xiàng)起主要作用。

約束條件為:

式中:U為n×k的隸屬度矩陣，ui，j表示第i個(gè)相點(diǎn)隸屬于第j個(gè)聚類中心的程度。Z為包含k個(gè)聚類中心矩陣。Di，j用歐氏距離的平方來(lái)度量第i個(gè)相點(diǎn)與第j個(gè)聚類中心的不相似程度:

1.1 目標(biāo)函數(shù)優(yōu)化

可以看出式(1)是帶約束的非線性優(yōu)化問(wèn)題。通常采用如下方法求該目標(biāo)函數(shù)的最小值:先固定U，更新Z使得目標(biāo)函數(shù)P最小;再固定Z，更新U使得目標(biāo)函數(shù)P最小。

固定U，更新的Z由下式求得:

其中1≤j≤k，1≤l≤m，式(3)與懲罰因子 λ 無(wú)關(guān)。

固定Z，用拉格朗日乘數(shù)法將式(1)變?yōu)闊o(wú)約束的優(yōu)化問(wèn)題。更新的U由下式求得:

1.2 聚集模糊K均值聚類算法

先設(shè)置懲罰因子λ，隨機(jī)選取初始聚類中心，通過(guò)交替更新Z和U，使得P(U，Z)不再變化，此時(shí)的U和Z為所求結(jié)果。由(1)式可知，優(yōu)化該算法中的第一項(xiàng)是使數(shù)據(jù)集聚類數(shù)目增加，然而優(yōu)化該算法中的第二項(xiàng)是使數(shù)據(jù)集聚類數(shù)目減少。算法中懲罰因子λ決定哪一項(xiàng)起主要作用，隨著λ增加，聚類過(guò)程向著減少聚類個(gè)數(shù)方向進(jìn)行，當(dāng)λ增大到一定程度時(shí)，該算法就將數(shù)據(jù)集聚成了一類。λ取值合適時(shí)，聚集模糊K均值聚類算法才能準(zhǔn)確得到真實(shí)聚類個(gè)數(shù)和聚類中心位置。下面介紹如何求取合適的懲罰因子λ。

2 聚類評(píng)價(jià)函數(shù)

文獻(xiàn)［10］提出聚類評(píng)價(jià)函數(shù)，該算法同時(shí)考慮了各聚類內(nèi)部的致密度(下式中第一項(xiàng))和各聚類之間的距離(下式中第二項(xiàng))。評(píng)價(jià)函數(shù)值越小，聚類效果越好。該評(píng)價(jià)函數(shù)如下:

求取數(shù)據(jù)集聚類個(gè)數(shù)和聚類中心位置的算法流程如下:

該算法中包含兩個(gè)循環(huán)，第一個(gè)循環(huán)是找出最小的懲罰因子λmin;第二個(gè)循環(huán)是計(jì)算最優(yōu)懲罰因子值。

第一個(gè)循環(huán)的算法流程:

(1)輸入初始聚類個(gè)數(shù)k，要求初始聚類個(gè)數(shù)大于真實(shí)聚類個(gè)數(shù)。

(2)設(shè)置λmin初始值為相點(diǎn)最大距離的0.01倍，并令t=1。

(3)運(yùn)行聚集模糊K均值聚類算法，使位置相同的聚類中心聚集為一個(gè)聚類中心kt=k－kshare

(4)如果kt=k，則執(zhí)行第二個(gè)循環(huán)，否則令λmin=λmin－0.1λmin，跳轉(zhuǎn)到步驟(3)。

第二個(gè)循環(huán)的算法流程:

(2)如果kt=1，則該算法執(zhí)行完畢，否則跳轉(zhuǎn)到步驟(1)。

上述第二個(gè)循環(huán)中最小評(píng)價(jià)函數(shù)值所對(duì)應(yīng)的懲罰因子即為最優(yōu)懲罰因子，將其代入聚類模糊K均值聚類算法中可以得出真實(shí)的聚類個(gè)數(shù)和聚類中心位置。在有無(wú)重疊數(shù)據(jù)集上，該算法都能準(zhǔn)確得到真實(shí)聚類個(gè)數(shù)和精確的聚類中心位置。

3 自適應(yīng)局部獨(dú)立分量分析算法

本文提出的改進(jìn)算法主要思想是在進(jìn)行傳統(tǒng)的K均值聚類前，先使用聚集模糊K均值聚類和聚類評(píng)價(jià)函數(shù)，求出數(shù)據(jù)集的真實(shí)聚類個(gè)數(shù)和聚類中心位置，從而減少傳統(tǒng)局部獨(dú)立分量分析的參數(shù)輸入，消除初始聚類中心的隨機(jī)選取對(duì)該降噪算法穩(wěn)定性的影響，同時(shí)提高了降噪的效果。

自適應(yīng)局部獨(dú)立分量分析算法流程:

(1)對(duì)一維含噪信號(hào)進(jìn)行相空間重構(gòu)得到軌跡矩陣X。

(2)對(duì)相空間重構(gòu)后的數(shù)據(jù)去均值處理:X'=XX。對(duì)X'用聚集模糊K均值聚類和評(píng)價(jià)函數(shù)求出其真實(shí)的聚類個(gè)數(shù)k和聚類中心。

(3)對(duì)X'運(yùn)用傳統(tǒng)的K均值聚類得到k個(gè)聚類Yi。其中聚類個(gè)數(shù)和初始聚類中心為步驟(2)中求得的聚類個(gè)數(shù)和聚類中心。

(4)對(duì)每一個(gè)聚類 Yi去均值處理:Y'i=Yi－i。對(duì)每一個(gè)聚類Y'i進(jìn)行獨(dú)立分量分析，并由高維空間向低維空間投影，達(dá)到信噪分離的目的。其中投影階數(shù)由最小描述長(zhǎng)度準(zhǔn)則來(lái)確定，得到低維數(shù)據(jù)為L(zhǎng)i。

(5)將所得的低維數(shù)據(jù)加上第(4)步中去除的均值Y'i，并按數(shù)據(jù)原來(lái)的順序?qū)⑦@k個(gè)聚類排列成完整的數(shù)據(jù)L'。

(6)將第(5)步所得到的數(shù)據(jù)L'加上第(2)步中去除的均值重構(gòu)成一維數(shù)據(jù)，即為降噪后的數(shù)據(jù)。

如果噪聲較大，可以重復(fù)步驟(1)～(6)，直到達(dá)到理想的效果為止。

4 仿真信號(hào)實(shí)驗(yàn)

Lorenz系統(tǒng)是典型的非線性動(dòng)力系統(tǒng)，Lorenz信號(hào)由下式產(chǎn)生:

圖1 降噪前信噪比為5.76 dB的Lorenz信號(hào)Fig.1 Noisy Lorenz signal waveform before noise reduction(SNR 5.76 dB)

圖2 降噪后信噪比為17.01 dB的Lorenz信號(hào)Fig.2 Noisy Lorenz signal waveform after noise reduction(SNR 5.76 dB)

在Lorenz信號(hào)中添加不同水平的白噪聲，比較本文提出的降噪算法與傳統(tǒng)局部獨(dú)立分量分析(LICA)和全局投影(GP)降噪算法的降噪效果，如表1所示。從表1可以看出:在高信噪比的情況下，全局投影降噪的降噪效果優(yōu)于LICA和本文提出的降噪算法;當(dāng)信噪比較低時(shí)，LICA和本文提出的改進(jìn)降噪算法的降噪效果優(yōu)于GP，并且本文提出的降噪算法的降噪效果要優(yōu)于LICA，同時(shí)改進(jìn)算法的每次降噪效果都不變，具有很好的穩(wěn)定性。

表1 幾種不同的降噪算法對(duì)不同噪聲水平的Lorenz信號(hào)降噪效果比較Tab.1 Comparison among the effects of several methods on Lorenz signal of different noisy level

5 自適應(yīng)局部獨(dú)立分量分析在齒輪箱故障診斷中的應(yīng)用

試驗(yàn)在傳動(dòng)系統(tǒng)故障試驗(yàn)臺(tái)上進(jìn)行，減速箱為一級(jí)減速。主動(dòng)齒輪為20齒，從動(dòng)齒輪為37齒，輸入轉(zhuǎn)速為14.23 Hz，采樣頻率為2 000 Hz，取1 024 個(gè)數(shù)據(jù)點(diǎn)。主動(dòng)齒輪有一個(gè)齒斷裂，則主動(dòng)齒輪產(chǎn)生的振動(dòng)沖擊時(shí)間間隔約為0.07 s。

圖3 原始信號(hào)及降噪信號(hào)Fig.3 The original signal and the signal after noise reduction

圖3(a)為減速箱振動(dòng)信號(hào)的時(shí)域波形，從時(shí)域圖中很難看出周期性的沖擊信號(hào);圖3(b)為采用本文降噪算法處理后的信號(hào)。從降噪后的信號(hào)時(shí)域圖中，可以明顯的看出周期沖擊，且沖擊間隔為0.07 s。說(shuō)明本文提出的降噪算法能有效地提取斷齒故障信號(hào)中沖擊成分。

6 結(jié)論

(1)與傳統(tǒng)的局部獨(dú)立分量分析相比，本文提出自適應(yīng)局部獨(dú)立分量分析降噪算法，減少了輸入?yún)?shù)(聚類中心數(shù)目)，克服了由于傳統(tǒng)K均值聚類導(dǎo)致的降噪效果的不穩(wěn)定性。將該算法應(yīng)用到含噪Lorenz信號(hào)和含有噪聲的齒輪沖擊信號(hào)中，都能得到比較好的降噪效果。

(2)噪聲較大時(shí)，可以采取多次循環(huán)的方法來(lái)逐步消除噪聲，在剛開始循環(huán)中，嵌入維數(shù)應(yīng)取大一些，噪聲變小時(shí)，嵌入維數(shù)取小些，但必須滿足嵌入定理。

［1］Shin K，Hammond J K，White P R.Iterative SVD method for noise reduction of low-dime-nsional chaotic time series［J］.Mechanical Sys-tem and Signal processing，1999，13(1):115－124.

［2］呂志民，張武軍，徐金梧.基于奇異譜的降噪方法及其在故障診斷技術(shù)中的應(yīng)用［J］.機(jī)械工程學(xué)報(bào)，1999，35(3):95－99.

［3］徐金梧，呂 勇，王海峰.局部投影降噪算法及其在非線性時(shí)間序列分析中的應(yīng)用［J］.機(jī)械工程學(xué)報(bào)，2003，39(9):146－150.

［4］呂 勇，李友榮，肖 涵，等.基于加權(quán)相空間重構(gòu)降噪及樣本熵的齒輪故障分類［J］.振動(dòng)工程學(xué)報(bào)，2009，22(5):462－466.

［5］Karhunen J， Malaroiu S. Locally linearind-ependent component analysis［C］.Proceedin-gs of the International Joint Conference on Neural Networks，1999，2:882－887.

［6］Karhunen J， Malaroiu S. Local linear inde-pendent component analysis using clustering［J］.Int.Journal of Neural Systems，2000，10(6):439 －451.

［7］Gruber P，Stadlthanner K.Denoising using local projective subspace methods［J］.Neuroc-omputing，2006，69(13 －15):1485－1501.

［8］肖 涵.基于高斯混合模型與子空間技術(shù)的故障識(shí)別研究［D］.武漢:武漢科技大學(xué)，2007.

［9］Li M J，Ng M K.Agglomerative fuzzy K-mea-ns clustering algorithm with selection of number of clusters［J］.IEEE Transactions on Knowledge and Data Engineering，2008，20(11):1519－1534.

［10］Sun H，Wang S.FCM-based model selection algorithms for determining the number of clusters［J］.Pattern Recognition，2004，37(10):2027－2037.