劉 巖 ，王金東 ，李 穎

（1.東北石油大學 機械科學與工程學院，黑龍江 大慶 163318；2.齊齊哈爾大學 機電工程學院，黑龍江 齊齊哈爾 161006）

1 引言

對于天然氣與石化行業(yè)的大型往復式壓縮機故障診斷而言，軸承類故障因其監(jiān)測的困難性及后果的突發(fā)性與危害性，決定了對其開展有效的故障診斷的意義重大[1，6]。故障狀態(tài)及壽命預測作為故障診斷技術不可或缺的環(huán)節(jié)，在往復機軸承類故障診斷中其重要性尤為突出。

大型往復機作為典型的非線性系統(tǒng)，其預測不準確的主要原因在于：系統(tǒng)動力模型與實際的差異、初始條件的影響2方面[1]。從預測模型角度分析，確定性的參數類回歸模型由于模型影響參數過多、影響因素考慮不足以及模型修正過分依賴主觀性等不足，造成可靠度和精度不高而限制其應用；非參數回歸方法是一種近年來在預測學中較熱門建模方法，特別適合于不確定性、非線性動態(tài)系統(tǒng)。另外，對于初始條件的影響，帕科德（Packard）在1980年提出用時間序列相空間重構理論，從靜態(tài)點與時間關系結合的角度再現了系統(tǒng)動態(tài)演化軌跡；為復雜非線性系統(tǒng)演化趨勢分析提供有效途徑[1-2]。學者們基于KNN及與相空間重構方法結合的改進算法方面開展了大量研究，如文獻[3]在股票價格預測中提出一種基于徑向基網絡的KNN預測模型；文獻[4]提出EMD分解與KNN結合的股票價格改進預測模型；文獻[5]將改進KNN算法和相空間重構結合在風電功率預測中取得較好的效果。

然而，往復機振動信號表現出的類周期震蕩和復雜的多重分形特征—即長期不可預測性，使得對該類振動信號預測而言，局限在振動時序離散點值的預測，而忽略因初始條件敏感性決定的預測可信時間的方法顯然不可取，同時預測結果受熵值復雜度和類周期振動的制約，有限個振動幅值因無法展現系統(tǒng)整周期包含的細節(jié)特征而失去意義。

基于以上分析，筆者將最大預測可信時間引入預測算法，對時序信號等間隔周期采樣，利用VMD算法提取主分量并提取各采樣時段的多重分形奇異譜特征參數構造預測序列，計算各模態(tài)分量的最大預測可信時間以確定預測可信時段，基于相空間重構的改進KNN創(chuàng)建初始建模域，以C_C算法變參數尋優(yōu)還原系統(tǒng)演化趨勢；算法將度量預測時間的熵值和基于模態(tài)分解的多重分形譜特征有效結合，將系統(tǒng)動力學模型和初始條件影響因素作為整體統(tǒng)一在算法中。

2 系統(tǒng)預測可信時間

在香農（Shannon）信息論中，對于M維時間序列的相空間點對，若將其劃分為k個尺寸為rM的單元，令單元序列si的生成概率為pi；其信息值Ik它正比于以精度r確定系統(tǒng)單元區(qū)域序列si（或稱軌道）所需信息；即Ik+1-Ik度量了系統(tǒng)從時間（k-1）τ到kτ的信息損失[1，6]。

系統(tǒng)的最大預測可信時間可用K（柯爾莫哥洛夫，Kolmogorov）熵表示，設t時刻信息量為I（t），經過Δt后的信息量為I（t+Δt），則有：I（t+Δt）=I（t）-KΔt；

令I（t）=1，當I（t+Δt）=0時，即可得到系統(tǒng)的最大預測可信時間T=1/K。由于Renyi熵具有性質Kq≤Kq+1，可用Renyi熵的下界K2作為K熵的近似估計[1]。式（1）可近似用式（2）表示：

式中：Cm（r）—m維相空間距離小于r的點對占的比例，即C_C算

法的相似容限，故有：

式（4）中當嵌入維數m飽和時，可得T≈1/K2。具體算法流程，如圖1所示。

圖1 預測可信時間計算流程Fig.1 The Flowchart of Maximum Dependable Prediction Time

3 變分模態(tài)分解

變分模態(tài)分解（VMD）是由Dragomiretskiy和Zosso于2014年提出的一種新的自適應信號分解方法。通過求解變分框架中的約束變分模型，自適應地將復雜的多分量振動信號分解成一組有限帶寬的模態(tài)函數，實現各信號分量的頻率分離[7]。

對任意信號x（t），分解步驟可分為兩部分：

3.1 建立變分模型：

（1）定義有限帶寬的內稟模態(tài)函數（Band-Limited Intrinsic Mode Function，BLIMF）；

（2）構建解析信號獲取單邊譜，然后混頻得：

（3）估計帶寬：通過高斯平滑（即L2范數梯度的平方根）對信號解調獲得各模態(tài)函數帶寬；

（4）建立優(yōu)化模型：引入約束條件，構造優(yōu)化變分模型如下：

式中：K—BLIMF 分量｛uk｝=｛u1，u2，L，uK｝的數目；

｛ωk｝=｛ω1，ω2，L，ωK｝—uk（t）中心頻率。

3.2 變分模型求解

（1）為獲得上述變分模型的最優(yōu)解，VMD通過引入二次罰因子 α 和拉格朗日乘數 λ（t）構建增廣拉格朗日函數 L（｛uk｝，｛ωk｝，λ），把約束問題轉化為無約束問題。

其中，二次懲罰因子在高斯噪聲的存在下保證了信號重構的準確性，而拉格朗日乘子可以保證模型約束的剛度；

（2）用乘數交替方向法更新｛uk｝，｛ωk｝和 λ，并尋求增廣拉格朗日函數的鞍點，變分模型的最優(yōu)值顯示在方程（8）里；由此，輸入信號x（t）被分解成K個BLIMF分量。

4 基于PSR的KNN預測模型

傳統(tǒng)的KNN預測算法原理是對原n個數據構造m個序列，找到其最鄰近的k個值，然后依據k個點的均值作為第（n+1）個數據的預測值[8]。算法認為，系統(tǒng)的內在復雜的映射關系均包含在歷史數據中，它既沒有針對歷史數據建立訓練模型，也沒有將歷史數據作處理，特別適合于非線性問題的快速分類[9]。

基于C_C算法對序列相空間重構（Phasespacereconstruction，PSR），改進了KNN算法對模型的決策性，由于空間重構與原動力系統(tǒng)方程微分同胚，能使相點分布表現出最優(yōu)空間支撐；將構造的空間矩陣行向量作為KNN的建模域，每次計算得到的k個鄰近點均值作為狀態(tài)空間的下一個演化軌跡點；以預測值滾動更新建模域，并重新計算嵌入維度m和延遲時間τ實現變參數預測，依此循環(huán)直至序列最大可信時間尺度。

具體預測步驟如下：（1）等間隔抽取測點的等長短時周期振動信號，構造長度為L的振動時序X0；（2）計算X0分解尺度并VMD分解，計算各BLIMF的最大預測可信時間尺度向量Ti=［T1，T2…Tn］；（3）提取各主模態(tài)不同時段的多重分形奇異譜特征值，形成特征序列；（4）創(chuàng)建不同模態(tài)特征序列初始建模域，通過C_C算法計算嵌入維度m和時延；（5）以該重構相空間點列為初始待分類序列，KNN算法計算該相型分布的k個鄰近；（6）以K個鄰近的加權均值為該奇異譜特征預測值，建模域滾動前移再預測；（8）重復（2）～（7）步；直至各模態(tài)的min｛Ti｝為截止時間尺度，獲得特征值預測序列。因此對于多重分形奇異譜特征向量而言，基于VMD和PSR（相空間重構）的改進KNN預測算法（稱為VMD＆PSR_KNN法）計算框圖，如圖2所示。

圖2 VMD＆PSR_KNN預測流程Fig.2 The Flowchart of VMD＆PSR_KNN Method

5 實例分析

結合以上分析與算法流程，課組對大慶天然氣公司壓氣站2D12型往復式壓縮機為研究對象；其中，壓縮機主要參數為：軸功率500kW，排氣量70m3/min，活塞行程240mm，電機轉速496r/min。采集設備大修前30天的大頭軸瓦敏感測點（缸體外曲軸箱處）振動信號；數據采集系統(tǒng)為湖北優(yōu)泰公司的uT3416型16通道數據采集儀；壓縮機結構與測點位置，如圖3所示。其中，黑色三角形表示二級軸瓦敏感測點位置。

5.1 算法參數確定

由于振動信號內部的分形本質限定了預測只能是有限時間長度的可信值；而基于故障識別的要求，預測得到的信號要展現至少一個周期內的振動特征，以滿足對其進一步分析和處理的需求；同時，過分強調增加采樣間隔，似乎能增到預測可信時間，但必然會與復雜系統(tǒng)動力學演化的長期不可預測性特性相悖。

為保證采樣的有效長度，對長期采樣序列每間隔1h取得2周期約2410個采樣點（采樣頻率10kHz），獲得30h的間隔采樣數據，為便于程序化計算，采用MRMR（最大相關最小冗余）法計算各個采樣時段的VMD預分解個數均為4個，計算各模態(tài)序列的最大預測可信時間間隔向量為：T=［43Vt 33Vt 35Vt 37Vt］，可知在計算點后33h內的預測值是可信的；KNN算法中取鄰近值K=3。

圖3 2D12往復式壓縮機結構如圖Fig.3 The Structure of 2D12 Reciprocating Compressor

多重分形奇異譜特征量中，譜中心值α0（fmax＝f（α0））代表信號長程相關性，其值越大相關性越強，其穩(wěn)定性雖好而可分性并不好[5-6]；通過構造包含對稱性和奇異性的譜差比△αdr=（αmax-α0）/（α0-αmin）；以物理意義明確的α0和△αdr形成穩(wěn)定、可分的特征向量作為預測指標，擬將其用于相空間重構改進的KNN預測模型加以預測。多重分形奇異譜形態(tài)參數示意，如圖4所示。

圖4 多重分形奇異譜形態(tài)參數示意Fig.4 Schematic Diagram of Shape Parameter for MSS

5.2 特征值預測結果

基于上述算法參數計算與分析，將該狀態(tài)時間序列經VMD分解的4個BLIMF的譜參數α0和Δαdr計算結果與實測結果，如圖5、圖6所示。

圖5 各主模態(tài)奇異譜中心值Fig.5 Central Valueof MSS for Each BLIMF

圖6 各主模態(tài)奇異譜譜差比Fig.6 Difference Ratioof MSS for Each BLIMF

從各個模態(tài)分量的特征參數預測趨勢可見，譜寬和中心值兩個參數其幅值和波動趨勢均與實測數據偏差較小，特別是前20個預測可信時間間隔以內，預測值的吻合度較高。

5.3 方法評估

為證明算法有效性，課組采用線性擬合回歸的方式[4]；構造方程為Ypredict=R×Xactual+Err，Ypredict和Xactual分別表示預測值和實測值，Err為擬和誤差，R為曲線斜率，其值越接近1則曲線預測效果越好，用相同數據分別對EMD_KNN、VMD_KNN和VMD&PSR_KNN三種預測方法做對比分析（其中在EMD分解中取前4個主模態(tài)，其信號能量比率已達到80%，已包含了高頻的特征信息[10]），兩特征參數回歸比較結果（圖略）。

回歸分析可見，EMD_KNN算法由于模態(tài)混疊等因素造成特征值預測結果不佳，而VMD法因模態(tài)分解方面由于EMD而使預測精度明顯提高，特別是結合空間重構變參數建模后，算法考慮了時序熵值的變化的影響，預測回歸結果最優(yōu)。

表1 各方法a0預測值誤差比較Tab.1 Comparison of Predicted Errors for a0by Different Methods

表2 各方法Δadr預測值誤差比較Tab.2 Comparison of Predicted Errors for Δadrby Different Methods

如表1、表2所示對三種方法得到的特征參數預測誤差用標準差（STD）和均方差（RMSE）匯總并比較，兩種基于VMD分解的預測算法得到的預測特征值與實測序列計算的特征值誤差較?。豢梢宰C實，文中多重分型譜參數KNN預測算法的有效性；更為重要的是，從基于多重分形譜形參數衰減指標程度判定，該故障已達到重度磨損閾值，這與停機大修檢查的實際結果是一致。

6 結論

（1）從熵值變化反映復雜系統(tǒng)初始敏感性角度將最大預測可信時間算法引入預測程序，優(yōu)化預測時間尺度；（2）從相空間重構與原系統(tǒng)狀態(tài)方程微分同胚原理出發(fā)，利用VMD算法改進KNN非參數預測模型缺陷，結合相空間重構變尺度構造建模域，再現系統(tǒng)動力學演化趨勢；（3）基于多重分形理論，選擇可分性和魯棒性較好的奇異譜特征參數，結合VMD分解方法提取表征精細分形的特征向量，提高非參數預測算法理論的可解釋性。通過特征提取與實測分析證實，相比于EMD和傳統(tǒng)KNN算法，基于VMD多重分形奇異譜的變參數相空間重構預測算法為設備健康狀態(tài)評估提供新的途徑和理論支撐。

［1］石博強，申焱華．機械故障診斷的分形方法：理論與實踐［M］．北京：冶金工業(yè)出版社，2001：288-300．（Shi Bo-qiang，Shen Yan-hua．Fractal Methods for Mechanical Fault Diagnosis：Theory and Practice［M］.BeiJing：Metallurgical Industry Press，2001：288-300.）

［2］劉巖，王金東，郭建華．基于參量優(yōu)化的多變量相空間重構方法［J］．機械傳動，2013（6）：5-7．（Liu Yan，Wang Jin-dong，Guo Jian-hua．Method of Phase Space Reconstruction of Multivariate Data based on Parameters Optimization［J］.Journal of Mechanical Transmission，2013（6）：5-7．）

［3］周廣旭．一種新的時間序列分析算法及其在股票預測中的應用［J］．計算機應用，2005（9）：2179-2181．（Zhou Guang-xu．RBF-based tmi e-series forecasting［J］.Journal of Computer Applications，2005（9）：2179-2181．）

［4］Lin Ai-jing，Shang Peng-jian，Feng Guo-chen.Application of EMD combined with KNN approach in financial time series forecasting［J］.Fluctuation Noise Letters，2012（11）：1250018.

［5］楊茂，賈云彭，穆鋼.基于改進KNN算法的風電功率實時預測研究［J］.電測與儀表，2014（24）：38-43.（Yang Mao，Jia Yun-peng，Mu Gang.Wind Power Real-time Prediction Research Based on the Improved KNN Algorithm［J］.Electrical Measurement&Instrumentation，2014（24）：38-43.）

［6］J.Wang，W.Wei，H.Li.Multifractal features extraction of bearing faults for reciprocating compressor［J］.Journal of Vibration and Shock，2008，27（S）：313-315.

［7］Dragomiretskiy K and Zosso D．Variational mode decomposition IEEE Trans［C］.Signal Process，2014，62（3）：531-544.

［8］Bagui S C，Bagui S，Pal K.Breast Cancer Detection using Nearest Neighbor Classification Rules［J］.Pattern Recognition，2003（36）：25-34.

［9］Zeng Y，Yang Y，Zhou L.Pseudo Nearest Neighbor Rule for Pattern Recognition ［J］.Expert Systems with Applications，2009 （36）：3587-3595.

［10］張維強．基于MSKPCA和SVM的轉子故障診斷模型及應用［J］.機械設計與制造，2015（10）：4-8.（Zhang Wei-qiang.Rotor fault diagnosis model and application based on multi-scale kernel principal component analysis and support vector machine［J］.Machinery Design&Manufacture，2015（10）：4-8.）