有界廣義GMM和無(wú)閾值遞歸圖的特征提取方法及應(yīng)用

趙心陽(yáng)，肖 涵，2，呂 勇，2

（1.武漢科技大學(xué)機(jī)械自動(dòng)化學(xué)院，湖北 武漢 430081；2.武漢科技大學(xué)冶金裝備及其控制教育部重點(diǎn)實(shí)驗(yàn)室，湖北 武漢 430081）

1 引言

齒輪傳動(dòng)廣泛應(yīng)用于機(jī)械設(shè)備中，但由于工況的復(fù)雜多變，齒輪自身也容易發(fā)生故障，統(tǒng)計(jì)表明，在齒輪箱的全部零件中，齒輪自身失效引起的設(shè)備故障約占全部故障的60%［1］，因此針對(duì)齒輪故障診斷方法的研究一直是近年來(lái)的熱點(diǎn)。齒輪振動(dòng)信號(hào)更多的表現(xiàn)出非線性、非平穩(wěn)特征是普遍共識(shí)［2］，學(xué)者往往將振動(dòng)信號(hào)的非線性特征，如Lyapunov指數(shù)、分形維數(shù)、熵等用于故障診斷，取得一定的效果［2-3］。但由于振動(dòng)信號(hào)通常信噪比較低，而上述參數(shù)對(duì)噪聲較敏感，導(dǎo)致其作為特征量時(shí)，故障識(shí)別精度不高。

遞歸性是動(dòng)力學(xué)系統(tǒng)的基本屬性之一，1987年Eckman基于相空間重構(gòu)理論提出遞歸圖（Recurrence Plot，RP）分析方法［4］，用以從遞歸圖的內(nèi)部層狀結(jié)構(gòu)揭示時(shí)間序列的復(fù)雜性和可預(yù)測(cè)性，包括傳統(tǒng)閾值遞歸圖（Thresholded Recurrence Plot，RP）和無(wú)閾值遞歸圖（Un-Thresholded Recurrence Plot，URP）。為了解決遞歸圖分析方法無(wú)法定量描述動(dòng)力學(xué)特性等問(wèn)題，文獻(xiàn)［5-6］提出傳統(tǒng)閾值遞歸圖特征提取方法-遞歸定量分析（Recurrence Quantifica?tion Analysis，RQA）。通過(guò)研究遞歸點(diǎn)的分布規(guī)律，研究者已提出的遞歸參量已達(dá)20余種，且已應(yīng)用于設(shè)備故障診斷領(lǐng)域。但遞歸參量的選擇往往由研究者主觀確定，且故障識(shí)別精度有待提高。由于URP中遞歸點(diǎn)間的距離參數(shù)被保留下來(lái)，因此對(duì)于同一時(shí)間序列，URP 中往往包含比RP 更多的信息量，近年來(lái)針對(duì)URP 的研究逐漸成為熱點(diǎn)。2018年，文獻(xiàn)［7］提出利用URP 分析方法檢測(cè)絕緣電力電纜因熱老化引起的劣化程度。文獻(xiàn)［8］根據(jù)URP中的歐式距離分布提出了一種基于無(wú)閾值遞歸矩陣的檢測(cè)器，用于檢測(cè)噪聲環(huán)境中信號(hào)的出現(xiàn)。但上述研究均是從定性角度分析，無(wú)法定量描述URP，限制了其應(yīng)用。這里提出一種采用有界廣義高斯混合模型（BGGMM）進(jìn)行URP特征提取的方法，用于URP的定量分析，并應(yīng)用于齒輪故障分類。

2 算法

2.1 無(wú)閾值遞歸圖

遞歸圖作為時(shí)間序列中重復(fù)模式的圖形表示，是分析時(shí)間序列周期性、混沌性以及非平穩(wěn)性的一個(gè)重要方法［2］。算法流程如下：

（1）對(duì)于一維時(shí)間序列xi(i=1，2，3，…N)，經(jīng)相空間重構(gòu)后，得到高維相空間：

式中：Nd=N-(d-1)τ，τ—延遲時(shí)間；d—嵌入維數(shù)。

（2）對(duì)于相空間中的軌跡x→i(i=1，2，…Nd，x∈Rd)，無(wú)閾值遞歸圖定義為：

式中：‖ ?‖—向量范數(shù)，常用范數(shù)包括1-范數(shù)，2-范數(shù)和∞-范數(shù)，這里采用∞-范數(shù)；與傳統(tǒng)閾值遞歸圖不同，無(wú)閾值遞歸圖由于沒(méi)有閾值ε的限定，其遞歸矩陣Ri，j(unthres)—相空間中兩個(gè)向量x→i和x→j之間的歐式距離。

2.2 有界廣義高斯混合模型

基于高斯分布的混合模型一直都是研究的熱點(diǎn)，廣泛應(yīng)用于模式識(shí)別、圖像處理和故障診斷等領(lǐng)域。其中高斯混合模型（GMM）［9］是應(yīng)用較多的混合模型之一，但由于其具有奇異性，即樣本數(shù)不足時(shí)參數(shù)估計(jì)會(huì)更加困難，從而導(dǎo)致算法發(fā)散。為解決GMM 的不足，學(xué)者提出更具魯棒性的廣義高斯混合模型（GGMM）［10］。相較于傳統(tǒng)GMM，GGMM 嵌入額外參數(shù)λ用以控制擬合曲線尾部凸起與否，因此可適應(yīng)更多類型數(shù)據(jù)。但是上述混合模型都未考慮數(shù)據(jù)邊界問(wèn)題，當(dāng)數(shù)據(jù)存在支撐邊界時(shí)，將導(dǎo)致參數(shù)估計(jì)存在較大的偏差。因此文獻(xiàn)［11］2014年提出有界廣義高斯混合模型（BGGMM）用以擬合高斯分布和非高斯分布數(shù)據(jù)，并考慮數(shù)據(jù)的邊界問(wèn)題，近年來(lái)廣泛用于模式識(shí)別等領(lǐng)域。文獻(xiàn)［12-13］成功將BGGMM 應(yīng)用于視頻和圖像降噪以及圖像分割領(lǐng)域。文獻(xiàn)［14］結(jié)合獨(dú)立分量分析方法（ICA）和BGGMM成功實(shí)現(xiàn)對(duì)關(guān)鍵字進(jìn)行無(wú)監(jiān)督識(shí)別。

2.2.1 算法簡(jiǎn)介

高斯混合模型（GMM）是單高斯概率密度函數(shù)的延伸，廣泛應(yīng)用于模式識(shí)別、圖像處理、故障診斷［6］等領(lǐng)域。

K個(gè)分量的GMM，其定義為：

式中：Θj—模型參數(shù)，ωj作為混合權(quán)重系數(shù)，滿足ωj≥0。

式中：μj和σj—混合模型中的均值和方差，參數(shù)λj控制該分布的尾部凸起與否；Γ(?)—伽馬函數(shù)。值得提出的是，GMM 和GGMM是BGGMM的特殊情況：當(dāng)λj=2，H(xi|Ωj)=1時(shí)，式（4）可表示為GMM的概率密度函數(shù)；當(dāng)H(xi|Ωj)=1時(shí)，表示為GGMM的概率密度函數(shù)。

2.2.2 參數(shù)估計(jì)

期望最大化（Expectation-Maximum，EM）算法在統(tǒng)計(jì)學(xué)中廣泛被用于參數(shù)的最大似然估計(jì)，同時(shí)也是基于高斯分布的混合模型中常用的參數(shù)估計(jì)方法［15］。EM算法以迭代的形式實(shí)現(xiàn)最大似然函數(shù)的最大化，包括E-Step 和M-Step，在這里所應(yīng)用的BGGMM中，最大似然函數(shù)定義如下：

為了最大化式（8），采用其負(fù)數(shù)方程作為誤差方程，此時(shí)最大化似然函數(shù)即等效為最小化誤差式（9），誤差方程定義如下：

在E-Step中，后驗(yàn)概率定義為：

在M-Step 中，根據(jù)式（10）得出的后驗(yàn)概率更新模型參數(shù)如下：

其中：

式中：γ—比例因子—從混合模型特定組分j的廣義高斯分布中提取的隨機(jī)變量；M—大正整數(shù)；當(dāng)x≥0時(shí)，sign(x)=1，否則sign(x)=0。

最后，對(duì)式（9）中的誤差函數(shù)進(jìn)行評(píng)估，檢查收斂準(zhǔn)則。若不滿足收斂條件，則從E-Step重新開(kāi)始下一次迭代，直到滿足收斂準(zhǔn)則為止。

3 URP統(tǒng)計(jì)定量分析

無(wú)閾值遞歸圖（URP）可以檢測(cè)遞歸點(diǎn)之間距離的微小變化，而如果這些變化小于預(yù)定義的截止距離，則在RP中無(wú)法檢測(cè)到。因此，URP中的模式可以比RP提供更多的信息。為了從統(tǒng)計(jì)學(xué)角度量化URP，這里提出采用BGGMM進(jìn)行URP統(tǒng)計(jì)定量分析方法，算法流程如下：

（1）輸入時(shí)間序列x(t)；

（2）通過(guò)式（1）和式（2）經(jīng)相空間重構(gòu)獲得無(wú)閾值遞歸矩陣Ri，j(unthres)，然后將遞歸矩陣中值的直方圖進(jìn)行歸一化處理；

（3）初始化模型參數(shù)Θ=(ωj，μj，σj，λj)；

（4）通過(guò)式（10）式（14）迭代更新模型參數(shù)Θ；

（5）檢測(cè)迭代收斂準(zhǔn)則。若不滿足，返回步驟（4）繼續(xù)迭代，直至滿足條件，迭代停止；

（6）輸出均值和方差，組成特征向量(μ，σ)，對(duì)URP進(jìn)行定量描述。

本節(jié)中，為了驗(yàn)證該方法的有效性，利用兩組仿真信號(hào)進(jìn)行分析。（1）周期信號(hào)：x(t)=sin(2π10t)；（2）Lorenz信號(hào)：

周期信號(hào)和Lorenz信號(hào)的時(shí)域波形圖、對(duì)應(yīng)的URPs以及歸一化后的歐式距離分布直方圖，如圖1所示。從圖1可看出，正弦周期信號(hào)URP的遞歸點(diǎn)沿主對(duì)角線平行或垂直分布，而混沌信號(hào)URP中出現(xiàn)大量平行與主對(duì)角線的線段。其對(duì)應(yīng)的分布直方圖近似服從于高斯分布，所不同的是，圖1（a）中的直方圖明顯存在數(shù)據(jù)邊界，數(shù)據(jù)也更加分散；圖1（b）中，數(shù)據(jù)更接近于高斯分布，且相對(duì)集中，這表明原始信號(hào)存在差異。經(jīng)BGGMM擬合（K=1），提取均值分別為1.45和1.92，方差分別為7.04和1.91，提取的均值和方差均存在差異，其中方差存在較大的差異也解釋了直方圖數(shù)據(jù)的集中程度的不同，這表明混合模型的參數(shù)可作為URP的定量分析指標(biāo)。

圖1 仿真信號(hào)時(shí)域圖及其對(duì)應(yīng)的URPs和歐式距離分布直方圖Fig.1 Time-Domain Plot and their URPs as well as the Euclidean Distance Distribution Histogram of Simulate Signals

4 齒輪故障分類識(shí)別實(shí)驗(yàn)

這里將該方法應(yīng)用于齒輪故障的分類，提出的齒輪故障分類算法的流程圖，如圖2所示。

圖2 齒輪故障分類算法的流程圖Fig.2 The Flowchart of Fault Classification Algorithm

4.1 實(shí)驗(yàn)裝置

為了驗(yàn)證所提方法的有效性，通過(guò)齒輪故障實(shí)驗(yàn)臺(tái)采集不同齒輪狀態(tài)下的振動(dòng)信號(hào)，實(shí)驗(yàn)臺(tái)及其結(jié)構(gòu)簡(jiǎn)圖，如圖3 所示。實(shí)驗(yàn)臺(tái)包括550w（220V-50Hz）電機(jī)，帶有一對(duì)嚙合齒輪的齒輪箱，聯(lián)軸器和磁粉制動(dòng)器。其中，電機(jī)的轉(zhuǎn)速和負(fù)載在（300～1217）r/min和（0～20）N·m范圍內(nèi)變動(dòng)；小齒輪的齒數(shù)為20，大齒輪37，模數(shù)均為3；振動(dòng)信號(hào)由安裝在輸入軸承座的垂直方向上的加速度傳感器采集。

圖3 齒輪故障實(shí)驗(yàn)臺(tái)及其結(jié)構(gòu)簡(jiǎn)圖Fig.3 Gear Failure Test Table and its Structure Diagram

4.2 時(shí)域信號(hào)及無(wú)閾值遞歸圖

采用上述實(shí)驗(yàn)臺(tái)，分別采集四種齒輪狀態(tài)的振動(dòng)信號(hào)各100組，采樣頻率10kHz，采樣時(shí)間0.1s。四種齒輪故障狀態(tài)的時(shí)域信號(hào)波形圖及其對(duì)應(yīng)的URP，如圖4、圖5所示。原始信號(hào)不進(jìn)行降噪處理，采用虛假鄰域法確定嵌入維數(shù)d=4，平均互信息算法［16］確定延遲時(shí)間τ=4。

圖4 四種齒輪狀態(tài)時(shí)域波形圖Fig.4 Time-Domain Signals of Four Gear States

圖5 四種齒輪狀態(tài)的無(wú)閾值遞歸圖Fig.5 Un-Thresholded Recurrence Plots of Four Gear States

從時(shí)域波形圖看，原始信號(hào)摻雜噪聲較多，無(wú)法直接從原始時(shí)域信號(hào)區(qū)分四種齒輪故障狀態(tài)；從URP中看出，遞歸圖中的遞歸點(diǎn)主要呈垂直分布，動(dòng)力學(xué)系統(tǒng)所決定。

4.3 故障識(shí)別結(jié)果及對(duì)比

無(wú)閾值遞歸圖中，向量之間的歐式距離大小由顏色深淺不一的遞歸點(diǎn)映射在圖像中，如圖5 所示。四種齒輪故障模式下的URP 顏色深度不同，表明遞歸點(diǎn)的歐式距離分布存在差異。應(yīng)用BGGMM 擬合URP 中遞歸點(diǎn)歐式距離分布直方圖，并提取均值和方差組成特征向量(μ，σ)用以齒輪的故障分類。分別以均值(μ)和方差(σ）為橫縱坐標(biāo)繪制，如圖6 所示。由圖6 可以看出，四種齒輪故障狀態(tài)下的振動(dòng)信號(hào)特征向量的分布存在較大的差異，除了斷齒和周節(jié)兩種狀態(tài)存在少量重疊外，四種故障齒輪的特征向量基本完全被分開(kāi)。

圖6 BGGMM提取的特征向量分布情況Fig.6 Distribution of Feature Vectors Extracted by BGGMM

為了驗(yàn)證該方法的有效性，將所采集的各齒輪故障狀態(tài)數(shù)據(jù)按1：1的比例分為訓(xùn)練數(shù)據(jù)和測(cè)試數(shù)據(jù)，即前50組數(shù)據(jù)用于訓(xùn)練，后50組數(shù)據(jù)用于測(cè)試，建立高斯混合模型，應(yīng)用貝葉斯最大似然分類器進(jìn)行分類。并與傳統(tǒng)RQA的10種遞歸定量［16］（RR、DET、L、Lmax、ENTR、LAM、TT、Vmax、T1/T2）作為特征向量時(shí)的識(shí)別率進(jìn)行比較，結(jié)果如表1所示。

表1 URP定量分析和傳統(tǒng)RQA方法的識(shí)別率Tab.1 Identification Rates of URP Quantitative Analysis and Traditional RQA Methods

從表1可以看出，基于BGGMM的URP定量分析方法的識(shí)別率明顯高于傳統(tǒng)RQA方法，四種齒輪故障狀態(tài)的識(shí)別率均高于90%，其中對(duì)齒輪正常狀態(tài)的識(shí)別率較高，達(dá)到100%。對(duì)于齒輪磨損和正常兩種狀態(tài)，較傳統(tǒng)RQA 方法分別提高10%和26%。表明這里所提出的齒輪故障分類方法具有一定的優(yōu)勢(shì)。

5 結(jié)論

近年來(lái)，遞歸定量分析廣泛應(yīng)用于故障診斷等領(lǐng)域，但其診斷精度有待提高，此外，目前缺乏針對(duì)無(wú)閾值遞歸圖定量分析的有效方法。這里提出一種采用BGGMM進(jìn)行URP特征提取的方法，用于定量分析URP，并成功應(yīng)用于齒輪的故障分類中。實(shí)驗(yàn)結(jié)果顯示，與傳統(tǒng)的RQA方法相比，該方法在齒輪故障識(shí)別方面具有更好的性能和更高的分類精度。

所提方法是從模式識(shí)別的角度，對(duì)無(wú)閾值遞歸圖的一種統(tǒng)計(jì)定量分析方法，將遞歸矩陣中歐式距離的分布直方圖的統(tǒng)計(jì)參數(shù)作為定量分析指標(biāo)。但提出的定量指標(biāo)所對(duì)應(yīng)的動(dòng)力學(xué)特征還需進(jìn)一步研究，此外，研究其它的統(tǒng)計(jì)學(xué)參數(shù)作為定量分析指標(biāo)也是后續(xù)的研究?jī)?nèi)容。