楊 楠,關(guān)佳亮,盛志鵬,王加龍,胥永剛
(北京工業(yè)大學(xué)先進(jìn)制造技術(shù)北京市重點(diǎn)實(shí)驗(yàn)室,北京 100124)
奇異譜分析(Singular Spectrum Analysis,SSA)是一種良好的時(shí)間序列分析方法。目前已被應(yīng)用于故障診斷中,并取得了良好的效果。N Qin 利用奇異譜分析解決了高速列車(chē)轉(zhuǎn)向架故障診斷問(wèn)題[1];Muruganatham,Bubathi 等基于奇異譜分析提出一種簡(jiǎn)單的軸承故障特征提取方法[2]。
奇異譜分析具有的良好的時(shí)間序列分析和預(yù)測(cè)特性[3]也被廣泛應(yīng)用于其他領(lǐng)域。趙炯采用奇異譜分析對(duì)土倉(cāng)壓力時(shí)間序列進(jìn)行預(yù)測(cè)研究[4];周天清采用奇異譜分析對(duì)金融時(shí)間序列進(jìn)行預(yù)測(cè)研究[5]。受此啟發(fā),采用奇異譜分析方法對(duì)軸承時(shí)間序列進(jìn)行分析,對(duì)軸承振動(dòng)趨勢(shì)進(jìn)行預(yù)測(cè)研究。首先,使用奇異譜分析方法對(duì)軸承振動(dòng)信號(hào)進(jìn)行分解,從而將復(fù)雜信號(hào)中的趨勢(shì)成分、周期成分和噪聲分離出來(lái);然后,將分解出的有效成分進(jìn)行重構(gòu);最后,利用數(shù)據(jù)的線性遞歸關(guān)系對(duì)重構(gòu)成分序列進(jìn)行預(yù)測(cè)。
奇異譜分析算法為兩個(gè)階段,第一階段為分解階段,主要包括嵌入和分解兩個(gè)過(guò)程;第二階段為重構(gòu)階段,主要包括分組和對(duì)角平均兩個(gè)過(guò)程。
對(duì)長(zhǎng)度為N 的真值時(shí)間序列xN=(x1,…,xn)。令窗長(zhǎng)L(1<L<N)為整數(shù),且滿(mǎn)足K=N-L+1。
步驟1:嵌入。
對(duì)長(zhǎng)度為N 的真值時(shí)間序列xN=(x1,…,xn)。令窗長(zhǎng)L(1<L<N)為整數(shù),且滿(mǎn)足K=N-L+1。
為實(shí)現(xiàn)嵌入,通過(guò)形成K=N-L+1=個(gè)延遲向量,將原始時(shí)間序列映射到大小為L(zhǎng) 的延遲向量序列,定義延遲向量Xi=(xi,…,xi+L-1)T,i=1,…,K,序列xN的軌跡矩陣為Hankel 矩陣。
式中Qi=XTPi,定義λi=‖Xi‖f2=‖Qi‖2。
步驟3:分組。令d=max{j:λj≠0}。將式(1)求解后,將索引集{1,…,d}分成m 個(gè)不相交的子集I1,…,Im。
定義XI=∑i∈IXi。式(1)可分解為:
選擇集合I1,…,Im的過(guò)程稱(chēng)為特征三元組分組。如果m=d,Ij={j},j=1,…,d,則相應(yīng)的分組稱(chēng)為基本分組。方法A 中幾個(gè)特征三元組的選擇與奇異值分解最優(yōu)性質(zhì)的時(shí)間序列近似相關(guān)。
步驟4:對(duì)角平均。在這一步中,將式(2)中每個(gè)矩陣XIj轉(zhuǎn)變成長(zhǎng)度為N 的新序列。令Y 為一個(gè)L×K 的矩陣,矩陣元素為yij,1≤i≤L,1≤j≤K,為簡(jiǎn)單起見(jiàn),令L≤K。通過(guò)對(duì)角線平均,使用下式將矩陣Y 轉(zhuǎn)換為系列
其中As={(l,k):l+k=s+1,1≤l≤L,1≤k≤K },|As|表示集合As中的元素?cái)?shù)量,即在“反對(duì)角線”上對(duì)矩陣元素求平均值。
將對(duì)角平均式(3)應(yīng)用于合成矩陣XIk,產(chǎn)生重構(gòu)序列χ(k)=。因此,初始序列(x1,…,xN)被分解為m 個(gè)重構(gòu)序列的總和:
由基本分組產(chǎn)生的重構(gòu)序列稱(chēng)為基本重構(gòu)序列。
如果存在a1,…,at,使得式(4)成立,那么時(shí)間序列由線性遞歸關(guān)系控制。
其中,數(shù)字t 是LRR 的順序,a1,…,at是LRR 的系數(shù)。如果t=r 是控制時(shí)間序列SNLRR 的最小階數(shù),則相應(yīng)的LRR 稱(chēng)為最小值。
如果時(shí)間序列的L-軌跡矩陣對(duì)于任何L 具有秩r,使得r<min(L,K),那么時(shí)間序列被稱(chēng)為有限秩r 的時(shí)間序列。注意,如果控制信號(hào)SN的最小LRR 具有r<N/2 的階數(shù)r,則SN具有秩r。
最小LRR 是唯一的。在不大于L-1 的所有非最小LRR 中,LRR 在SSA 周期性預(yù)測(cè)中效果是最好的。
根據(jù)可分離性條件,窗口長(zhǎng)度是L 應(yīng)該在L<N/2 范圍內(nèi)取越大越好,當(dāng)要提取具有已知周期的周期性分量時(shí),選擇可以被周期整除的窗口長(zhǎng)度,雖然窗口長(zhǎng)度的選擇很重要,但是L 的微小變化對(duì)結(jié)果的影響微乎其微。
時(shí)間序列分量可以基于以下原理進(jìn)行識(shí)別:特征向量的形式復(fù)制了生成此特征向量的時(shí)間序列分量的形式。因此,特征向量圖可以用來(lái)輔助鑒定。此外,正弦波可近似地產(chǎn)生具有相同頻率和相移為π/2 的兩個(gè)正弦波分量。因此,一對(duì)特征向量的散點(diǎn)圖可形成近似多邊形,該多邊形可以幫助識(shí)別出周期性的正弦波??蛇x擇幾個(gè)主要的本征因子解決信號(hào)提取中平滑和信號(hào)降噪的問(wèn)題。
W—關(guān)聯(lián)矩陣是一種包含重構(gòu)時(shí)間序列分量之間加權(quán)相關(guān)的矩陣,他可以用來(lái)識(shí)別有效的分離信息。分離良好的組件具有很小的相關(guān)性,而分離效果差的組件有很大的相關(guān)性。因此,通過(guò)查看W—關(guān)聯(lián)矩陣,可以找到識(shí)別重構(gòu)序列分量的相關(guān)性,并且可將該信息用于隨后的分組。
使用奇異譜分析方法對(duì)張家口宣龍三高強(qiáng)度棒線材廠28#精軋機(jī)2010 年9 月6 日到2010 年12 月30 日共116 d 數(shù)據(jù)的每日峰值進(jìn)行分析。
為便于分析,將該信號(hào)分為兩部分,取該信號(hào)前98 d 的數(shù)據(jù)作為原始測(cè)試數(shù)據(jù),對(duì)該信號(hào)進(jìn)行奇異譜分析,將該信號(hào)的周期成分、趨勢(shì)成分和噪聲進(jìn)行分離,取余下的18 d 數(shù)據(jù)作為驗(yàn)證數(shù)據(jù),與奇異譜分析得到的的預(yù)測(cè)信號(hào)進(jìn)行對(duì)比,驗(yàn)證預(yù)測(cè)效果。
對(duì)該信號(hào)前98 d 的數(shù)據(jù),取窗口長(zhǎng)度為L(zhǎng)=48,對(duì)其進(jìn)行奇異譜分析,得到長(zhǎng)度為48 的特征值序列,如圖1 所示,由圖中可知,前10 個(gè)特征值包含了信號(hào)的絕大部分信息,因此取前10 個(gè)特征值對(duì)應(yīng)的特征向量進(jìn)行分析。經(jīng)過(guò)計(jì)算,第1 個(gè)特征三元組在原始信號(hào)中占據(jù)97.91%的比重,包含了原始信號(hào)的絕大部分信息,其他特征三元組在原始信號(hào)占據(jù)的比重較低。
前10 個(gè)特征值對(duì)應(yīng)的特征向量如圖2 所示,從特征向量圖可以大體看出,第1 個(gè)特征向量與原始信號(hào)的趨勢(shì)一致,說(shuō)明第一個(gè)特征三元組包含了信號(hào)的趨勢(shì)成分;第2 個(gè)和第3 個(gè)特征向量是具有相同頻率和相移為π/2 的兩個(gè)正弦波分量,根據(jù)2.2的分析可初步判定第2 個(gè)和第3 個(gè)特征三元組包含了信號(hào)的周期成分,其他特征三元組包含了信號(hào)的噪聲成分。
為進(jìn)一步驗(yàn)證通過(guò)奇異譜分析方法是否可以將周期性成分分離出來(lái),對(duì)前10 個(gè)特征向量做散點(diǎn)圖,如圖3 所示。由散點(diǎn)圖可知,第2 個(gè)和第3 個(gè)特征三元組對(duì)應(yīng)的特征向量對(duì)是一個(gè)多邊形,根據(jù)2.2 的分析可知,第2 個(gè)和第3 個(gè)特征三元組包含了信號(hào)的周期成分,這說(shuō)明奇異譜分析方法可以將周期成分分離出來(lái)。
為評(píng)估不同特征三元組之間的獨(dú)立性,由2.3 分析可知,W—關(guān)聯(lián)矩陣圖可有效反映分離信息,因此對(duì)前10 個(gè)重構(gòu)成分做W—關(guān)聯(lián)矩陣圖,W—關(guān)聯(lián)矩陣圖如圖4 所示,顏色越淺表示相關(guān)性越小。由W—關(guān)聯(lián)矩陣圖可知,序列C1 與其他序列完全獨(dú)立;C2 和C3 屬于周期性成分,有較高的相關(guān)性,;C3 和C6、C7、C8 有較小的相關(guān)性;C4 與其它序列幾乎不相關(guān),C5 與C6、C7、C8 之間有較小的相關(guān)性;但是其他重構(gòu)成分之間的相關(guān)性較強(qiáng),甚至未完全分離。因此,選擇前5 個(gè)特征三元組進(jìn)行重構(gòu),并在前5 個(gè)特征三元組基礎(chǔ)上對(duì)信號(hào)進(jìn)行預(yù)測(cè)。
為驗(yàn)證所取重構(gòu)成分序列的合理性,分別取第1個(gè)特征三元組的重構(gòu)成分序列和原始信號(hào)進(jìn)行對(duì)比,第1-3 個(gè)特征三元組的重構(gòu)成分序列和原始信號(hào)進(jìn)行對(duì)比,第1-5 個(gè)特征三元組的重構(gòu)成分序列和原始信號(hào)進(jìn)行對(duì)比。由圖5 可知,第1 個(gè)特征3 元組的重構(gòu)成分序列和第1-3 個(gè)特征三元組的重構(gòu)組件都不能夠很好的反映信號(hào)的趨勢(shì)和周期性,只有第1-5 個(gè)特征三元組可以良好的反映信號(hào)的趨勢(shì)和周期特性,因此選取第1-5 個(gè)特征三元組的重構(gòu)成分序列對(duì)信號(hào)進(jìn)行預(yù)測(cè)。
圖1 奇異值分解后的前48 個(gè)特征值
取前5 個(gè)重構(gòu)組件C1-C5 對(duì)后18 d 數(shù)據(jù)進(jìn)行線性遞歸預(yù)測(cè),預(yù)測(cè)結(jié)果如圖6 所示,由圖6 可知,因此,對(duì)由奇異譜分析得到的重構(gòu)成分序列進(jìn)行線性遞歸預(yù)測(cè),可以實(shí)現(xiàn)良好的預(yù)測(cè)效果。將該方法應(yīng)用于軸承振動(dòng)信號(hào)中,可以實(shí)現(xiàn)對(duì)軸承振動(dòng)趨勢(shì)的預(yù)測(cè)。
圖2 前10 個(gè)特征值對(duì)應(yīng)的特征向量及前10 個(gè)特征值能量占比
圖3 奇異譜分析處理后的特征向量對(duì)散點(diǎn)
圖4 奇異譜分析處理后的-關(guān)聯(lián)矩陣
圖5 重構(gòu)成分序列與原始信號(hào)對(duì)比
圖6 預(yù)測(cè)結(jié)果對(duì)比
奇異譜分析對(duì)信號(hào)進(jìn)行處理可得到如下結(jié)論:
(1)奇異譜分析方法可以把信號(hào)的趨勢(shì)、周期成分和噪聲有效地分離出來(lái)。
(2)對(duì)奇異譜分析方法處理得到的特征三元組進(jìn)行重構(gòu),并對(duì)其重構(gòu)組件進(jìn)行預(yù)測(cè),可以實(shí)現(xiàn)良好的預(yù)測(cè)效果。
(3)奇異譜分析方法的窗口長(zhǎng)度選擇主要靠經(jīng)驗(yàn)選取,如果可以自適應(yīng)選擇窗口長(zhǎng)度,將提高該方法的準(zhǔn)確性和使用效率。