基于多重分形與BGSA-PNN的水電機組振動信號狀態識別

安宇晨，鄭 陽，陳啟卷，席 慧，閆懂林，游仕豪

（武漢大學動力與機械學院，武漢430072）

0 引 言

水電機組是電站的重要設備，其運行的穩定性對電站安全運行影響重大，振動則是影響水電機組穩定性的重要因素。由于水電機組受到機械、水力和電氣因素的共同作用，其振動信號具有明顯的非平穩、非線性特點［1］。對于水電機組的振動信號，通常從時域、頻域或時頻域的角度，通過幅值域分析［2］、傅里葉變換、小波變換、Hilbert-Huang 變換等方法進行處理［3］。但是這些方法各有不足之處：傅里葉變換不能分析非平穩信號的時域特征，小波變換雖然能夠分析信號的時頻特征，但是基函數的選擇對分析結果的影響較大，且缺乏統一理論依據，Hil?bert-Huang 變換雖然能夠實現非平穩信號的自適應分析，但存在端點效應和模態混疊現象，可能導致分解結果失真，甚至產生虛假的模態分量等情形，嚴重影響信號處理的效果，多位學者針對上述問題進行了深入分析［2-5］。

研究表明，水電機組振動信號呈現出分形特征，是典型的分形信號［6］，因此可使用多重分形去趨勢波動分析（multifractal detrend fluctuation analysis，MFDFA）處理。該方法以分形幾何為基礎，研究整體與局部間自相似性與不規則性，具有不依賴系統模型即可實現系統狀態信息提取的特點［6，7］，從而避免了常見的時頻域信號分析方法所存在的不足，基于此即可實現信號的特征提取和分類識別。陳龍［8］分析了不同類型軸承的時間序列信號并進行研究，證明了軸承振動信號具有多重分形特性；林近山［9］將多重分形分析方法引入至齒輪箱的故障診斷中，證明了基于該方法提取的特征對齒輪箱狀態變化十分敏感，并根據多重分形譜特征的均值與均方差選擇特征，最終實現了對齒輪箱正常、輕度磨損、中度磨損和斷齒4 種狀態的識別；裴峻峰［10］基于聲信號的多重分形譜特征和改進的SVM 實現了往復泵的故障診斷，提高了往復泵故障診斷準確率；唐利平［11］提取刀具磨損實驗聲發射信號的多重分形譜參數，對比不同狀態下多重分形譜特征的分布情況選擇特征，實現了刀具磨損狀態的識別。

目前雖已有學者從多重分形的角度實現了設備的故障診斷，但由于多重分形分析方法提取出的特征種類較為豐富，通常根據不同故障下各種特征數值的分布情況進行特征的選擇［8-11］，或者直接將譜圖曲線輸入至分類器中實現分類［7］，因此亟須在故障特征選擇方面開展進一步研究。

1 基于多重分形分析與BGSA-PNN 的狀態識別算法

本文提出的狀態識別算法主要由特征提取、特征選擇與分類器參數的聯合尋優、狀態識別模型構建三部分組成，涉及的理論和方法如下。

1.1 多重分形去趨勢波動分析

多重分形去趨勢波動分析方法（multifractal detrend fluctua?tion analysis，MFDFA）由Kantelhardt于2002年提出［12］，可以精確地計算出振動信號的多重譜。多重分形譜圖的特征對振動信號的微弱變化較為敏感，因此其可作為振動信號的特征，用于設備的狀態監測與故障診斷中［9，13-15］。多重分形去趨勢波動分析的計算過程如下。

對于時間序列xk(k= 1，2，…，n)，計算其累積離差Y(i)：

將序列Y劃分成Ns個長度為s的子序列，由于序列Y的長度可能不是子序列長度的整數倍，從反方向重復這一過程，共得到2Ns個子序列。對于每個子序列，使用最小二乘法擬合并計算其方差：

其中yv(i)為各個子序列的擬合函數。求取所有子序列方差的平均值以計算q階波動函數：

式中：q可取任意整數。當q取各個整數時，分析Fq(s)與s的雙對數圖，確定Fq(s)與s間的冪律關系：

式中：h(q)為廣義Hurst 指數，對于具有多重分形特征的序列，h(q)將隨q的變化而變化。對于多重分形體，標度指數τ(q)可通過如下關系計算：

根據Legendre 變換可以計算奇異性指數α和多重分形譜f(α)：

對于多重分形體，其多重分形f(α) -α譜圖為上凸單峰函數，多重分形譜圖特征的含義如下［8，11，13，14］：多重分形譜寬度Δα=αmax-αmin反映信號多重分形特征的強弱，即信號波動的劇烈程度，Δα越大，信號的波動越劇烈，信號的能量越強，多重分形特征越強，反之則信號波動越平均，多重分形特征越弱；Δf=f(αmin)-f(αmax)反映了振動信號中大、小峰值所占有的比例，Δf> 0 則信號的奇異性較強，Δf< 0 則信號的奇異性較弱，其中f(αmin)反映信號的大波動特征，f(αmax)反映了信號的小波動特征；α0=α(f=fmax)描述了振動信號的不規則性，該數值越大，則信號隨機性越強，信號越趨向于隨機振動。

將多重分形譜圖的特征組成信號的特征向量輸入至分類器中，即可實現信號的模式識別，由于提取到的特征較多，可基于二進制引力搜索算法對特征進行降維。

1.2 二進制引力搜索算法

二進制引力搜索算法［16］（binary gravitational search algo?rithm，BGSA）由Reshedi 提出，適合解決離散和二進制變量的優化問題。在粒子屬性的計算方式上與常規的引力搜索算法類似，主要區別在于粒子每個維度下坐標的取值范圍為0 或1［17，18］。該算法主要流程如下：

定義搜索空間中第i個粒子的質量Mi(t)為：

式中：i= 1，2，…，N，N為粒子的個數為粒子在第l維下的位置坐標，為第i個粒子在時間為t時相對于最佳和最差粒子的適應度的占比，用于計算粒子的質量；fiti(t)為時間為t時第i個粒子的適應度，worst(t)和best(t)分別代表該時間下，所有粒子中最差和最佳適應度的值。在時間為t時，第j個粒子對第i個粒子在第l維空間內的引力為：

式中：Rij(t)為第j個粒子和第i個粒子間的距離；ε為保證分母不為0 的極小常量為引力常數，G0為初始引力值；α為衰減系數；T為最大時間即最大迭代次數。因此第i個粒子受到的合力粒子在合力的作用下將產生加速度迭代過程中粒子的位置按照以下公式更新：

在引力搜索算法中，粒子的位置坐標即為待尋優的參數信息，因此需要將待尋優的參數轉換為粒子的位置坐標。根據本文研究內容，將粒子的尋優空間維數確定為17 維，每個維下的坐標只有0 或1，其中第1～10 維對應MFDFA 提取出的10 個特征，將坐標為1 的特征組合成為降維后的特征集，將第11～17 維坐標組合成的7位二進制數轉換為十進制后作為分類器參數。

1.3 概率神經網絡

概率神經網絡（Probabilistic Neural Network，PNN）由D.F.Specht 于1989年提出［19］，與BP 神經網絡和徑向基神經網絡（Radial Basis Function Neural Network）相比，具有網絡結構簡單、收斂速度快、分類準確等優點［20］，更加適用于分類問題，由輸入層、隱層、求和層和輸出層4層組成，每層的作用如下：

輸入層：接收輸入數據并傳至隱層，該層節點個數即為樣本的特征個數。

隱層：該層神經元個數即為輸入樣本個數，該層通過連接權值與輸入層連接并接收數據，計算輸入數據與隱層中心的距離，該層輸入輸出之間的關系如下：

式中：i= 1，2，…，M，Φij為第i類神經元的第j個神經元對應的隱層輸出；σ為平滑參數；d為樣本的特征維數；i為輸入樣本的類別；xij為第i類的第j個中心。

求和層：該層神經元與隱層中對應的神經元連接，接收隱層數據，對隱層同一類神經元的輸出結果求加權平均，并將計算結果輸出至輸出層。

輸出層：根據求和層計算結果，判斷該輸入所對應的輸出。

1.4 基于多重分形分析與BGSA-PNN 的特征提取及狀態識別

基于多重分形分析與BGSA-PNN 的特征提取方法流程如圖1所示，其步驟為：

圖1 基于多重分形分析與BGSO-PNN的特征提取與識別流程Fig.1 Flowchart of state identification of vibration signal based on multifractal and BGSA-PNN

（1）特征提取：使用多重分形方法分析水電機振動信號數據，繪制其多重分形譜圖并提取多重分形特征。

（2）初始化BGSA 算法粒子：運行BGSA 算法，根據設定參數初始化粒子。

（3）解釋粒子信息：對于BGSA 算法中的某個粒子，根據粒子的信息進行特征降維，并選擇PNN的參數。

（4）訓練分類器：將降維后的數據按照80%和20%的比例劃分為訓練集和測試集，使用訓練集訓練PNN網絡。

（5）輸出PNN 分類準確度：將測試集數據輸入至經過訓練的PNN網絡中，輸出得到分類器分類精度，即BGSA算法中該粒子的適應度函數值。

（6）重復步驟3至步驟5，得到每個粒子的適應度函數值。

（7）判斷是否達到BGSA 算法停止的條件，若未達到停止條件，則更新粒子位置，并返回至步驟3，若已達到停止條件，則停止BGSO算法并結束。

2 軸承數據驗證

軸承數據來源于美國凱斯西儲大學軸承數據中心網站，由正常、軸承內圈故障、外圈故障、滾動體故障4 種狀態下的數據組成。本文選擇4 種軸承狀態下的數據各50 條，基于多重分形分析與BGSA-PNN 方法提取并選擇數據的特征，并使用PNN進行狀態識別，以驗證本文提出的方法。

2.1 軸承信號特征提取

提取4 種狀態下的軸承信號，計算并繪制軸承信號的多重分形譜圖如圖2 所示。對比發現軸承信號在4 種狀態下的多重分形譜圖差異明顯，可使用多重分形譜的圖形特征進行軸承狀態的識別。本文提取多重分形譜圖中的αmin、αmax、f(αmin)、f(αmax)、Δf=f(αmin)-f(αmax)、Δα=αmax-αmin、α0=α(fmax)、αL=α0-αmin、αR=αmax-α0以及曲線的二階矩μ作為軸承信號的MFDFA 特征，并基于此進行特征降維工作。部分軸承信號的MFDFA特征數據見表1。

表1 部分軸承信號特征數據Tab.1 Feature data of the partial bearing signal

圖2 4種狀態下的軸承多重分形譜圖Fig.2 Multifractal spectra of bearings in four states

考慮到提取的特征量較多，構成的特征向量維度較高，其中可能存在一些冗余特征，將增加算力的消耗，延長運算時間，甚至可能會影響分類器的識別效果，因此有必要進行降維處理。

2.2 基于BGSA-PNN 聯合尋優的特征降維與狀態識別

本文選擇PNN 作為分類器，使用二進制引力搜索算法進行特征與分類器參數的聯合優化。以特征向量和PNN 網絡參數作為BGSA 算法的優化對象，PNN 網絡模型作為適應度函數，PNN 網絡輸出的準確度作為其適應度函數值。尋優算法的參數為：群體規模N= 20，最大迭代次數max_it= 50。分類準確度隨迭代次數的變化曲線如圖3所示。

圖3 分類準確度隨迭代次數的變化曲線Fig.3 Accuracy curve of the iterations

為避免隨機性，重復進行50 次實驗并記錄，取實驗結果的平均值作為PNN 網絡參數的最終確定值，選擇出現概率大于50%的特征作為最終降維結果。PNN 參數值最終確定為0.15，各特征出現頻次如圖4 所示，其中特征序號1～10 分別代表特征：αmin、αmax、Δα、Δf、α0、αL、αR、f(αmin)、f(αmax)、μ，根據圖4 選擇αmin、α0、αL、f(αmin)、μ作為特征降維結果，構成特征向量。

圖4 特征出現頻次統計Fig.4 Feature appearance statistics

將降維后的特征向量和未經過降維的MFDFA 特征輸入至同參數PNN 進行識別。各重復50 次實驗，統計其準確度與耗時，結果如表2 所示。從表2 可以看出，對于軸承故障數據，降維前后的MFDFA 特征識別的平均準確度均接近99%，證明了MFDFA 特征可應用于信號的模式識別中，且經過降維后，計算時間得到縮短。

表2 基于軸承數據的不同方法的識別準確度與計算時間統計Tab.2 Statistics of different Recognition accuracy and calculation time of different methods of bearing data

3 目標電站實測數據分析

實測數據采集自國內某電站機組，該電站水輪機型號為HLA835i-LJ-305，發電機型號為SF55-30/6400，機組設計水頭69 m，最高水頭85 m，額定有功功率為55 MW，額定轉速200 r/min。分別選擇該電站機組在穩定運行區、小負荷工況區和渦帶運行區3種工況下的數據各50條進行分析。

機組在3種工況下的水導軸承X向擺度的波形如圖5所示，其多重分形譜圖如圖6 所示，3 種運行工況下部分數據的特征值如表3所示。

表3 部分水導軸承X向擺度的特征數據Tab.3 Feature data of the partial water guide bearing X-direction swing waveform

圖5 水導軸承X向擺度波形Fig.5 X-direction swing waveform of water guide bearing

圖6 3個工況下水導軸承X向擺度的多重分形譜圖Fig.6 Multifractal spectra of X-direction swing of water guide bearing under three working conditions

使用BGSA-PNN 進行特征降維和PNN參數的聯合尋優，其中二進制引力搜索算法的參數設置為：群體規模N= 50，最大迭代次數max_it= 50，優化算法中的適應度函數值即分類準確度隨迭代次數變化的曲線如圖7 所示。重復實驗50 次，根據實驗結果確定PNN 參數值為0.013，繪制各個特征的出現頻次如圖8 所示，由于αmin、Δα、αL、f(αmin)、μ五維特征出現的頻次明顯高于其他特征，故將其選擇為分類器特征降維后的特征集。

圖7 電站實測數據的分類準確度隨迭代次數的變化曲線Fig.7 Classification accuracy curve of the measured data of the power station with iterations

圖8 電站實測數據特征出現頻次統計Fig.8 Feature appearance statistics of the measured data of the power station

經過降維，水電機組振動信號的MFDFA 特征維數從10 降至5。使用經過降維的MFDFA 特征進行一次識別的結果如圖9所示，其中類別1 表示小負荷運行工況、類別2 表示渦帶工況、類別3表示穩定運行工況，由圖9可以看出，基于BGSA-PNN 降維的MFDFA 特征的識別結果與實際情況十分接近，幾乎能夠完全識別渦帶工況和穩定運行工況，對于小負荷運行工況則能夠完全識別。

圖9 電站實測數據特征降維后的分類結果Fig.9 Classification results after dimension reduction of measured data features of power plants

作為對比，使用降維前后的MFDFA 特征和機組實測信號的EMD 模糊熵作為特征向量，輸入至同參數的PNN 進行識別分類，重復進行50 次實驗并將其統計結果記錄于表4。對比可見，水電機組振動信號的MFDFA 特征在識別準確度與時間上優于EMD 模糊熵方法，且經過降維后，神經網絡分類器對于機組狀態的識別效率得到了提升，對于水電機組的實測振動信號識別效果較好。

表4 基于實測數據的不同方法的識別準確度與計算時間統計Tab.4 Statistics of Recognition accuracy and calculation time of different methods of the measured data

4 結 論

本文提出了一種基于多重分形分析與BSGA-PNN 的水電機組信號特征提取方法。其特點在于通過多重分形分析繪制水電機組振動信號的多重分形譜圖并提取其特征，并采用BG?SA-PNN 進行聯合尋優，對PNN 參數進行尋優的同時開展多維特征向量降維操作。根據4種故障狀態下的軸承數據和三種運行工況下的水電機組實測數據進行對比驗證，證明了MFDFA特征能夠有效表征振動信號中包含的信息，且經過BGSA-PNN聯合尋優后，識別效率得到了提高。經與基于EMD-PNN 的傳統狀態識別方法進行對比，發現本文方法在識別準確度和計算時間上均優于傳統方法。本文方法有望推廣應用于更多類型的水電機組振動信號狀態識別中，提高電站安全、穩定運行的水平。 □