楊 楠，關(guān)佳亮，盛志鵬，王加龍，胥永剛

（北京工業(yè)大學(xué)先進(jìn)制造技術(shù)北京市重點(diǎn)實(shí)驗(yàn)室，北京 100124)

0 引言

奇異譜分析（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è)。

1 奇異譜分析算法

奇異譜分析算法為兩個(gè)階段，第一階段為分解階段，主要包括嵌入和分解兩個(gè)過(guò)程；第二階段為重構(gòu)階段，主要包括分組和對(duì)角平均兩個(gè)過(guò)程。

對(duì)長(zhǎng)度為N 的真值時(shí)間序列xN=（x1，…，xn）。令窗長(zhǎng)L（1＜L＜N）為整數(shù)，且滿足K=N-L+1。

1.1 第1 階段：分解

步驟1：嵌入。

對(duì)長(zhǎng)度為N 的真值時(shí)間序列xN=（x1，…，xn）。令窗長(zhǎng)L（1＜L＜N）為整數(shù)，且滿足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。

1.2 第2 階段：重構(gòu)

步驟3：分組。令d=max{j：λj≠0}。將式（1）求解后，將索引集{1，…，d}分成m 個(gè)不相交的子集I1，…，Im。

定義XI=∑i∈IXi。式（1）可分解為：

選擇集合I1，…，Im的過(guò)程稱為特征三元組分組。如果m=d，Ij={j}，j=1，…，d，則相應(yī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)序列稱為基本重構(gòu)序列。

1.3 線性遞歸關(guān)系

如果存在a1，…，at，使得式（4）成立，那么時(shí)間序列由線性遞歸關(guān)系控制。

其中，數(shù)字t 是LRR 的順序，a1，…，at是LRR 的系數(shù)。如果t=r 是控制時(shí)間序列SNLRR 的最小階數(shù)，則相應(yīng)的LRR 稱為最小值。

如果時(shí)間序列的L-軌跡矩陣對(duì)于任何L 具有秩r，使得r＜min（L，K），那么時(shí)間序列被稱為有限秩r 的時(shí)間序列。注意，如果控制信號(hào)SN的最小LRR 具有r＜N/2 的階數(shù)r，則SN具有秩r。

最小LRR 是唯一的。在不大于L-1 的所有非最小LRR 中，LRR 在SSA 周期性預(yù)測(cè)中效果是最好的。

2 可分離性和參數(shù)選擇

2.1 窗口長(zhǎng)度選擇

根據(jù)可分離性條件，窗口長(zhǎng)度是L 應(yīng)該在L＜N/2 范圍內(nèi)取越大越好，當(dāng)要提取具有已知周期的周期性分量時(shí)，選擇可以被周期整除的窗口長(zhǎng)度，雖然窗口長(zhǎng)度的選擇很重要，但是L 的微小變化對(duì)結(jié)果的影響微乎其微。

2.2 周期性成分識(shí)別

時(shí)間序列分量可以基于以下原理進(jìn)行識(shí)別：特征向量的形式復(fù)制了生成此特征向量的時(shí)間序列分量的形式。因此，特征向量圖可以用來(lái)輔助鑒定。此外，正弦波可近似地產(chǎn)生具有相同頻率和相移為π/2 的兩個(gè)正弦波分量。因此，一對(duì)特征向量的散點(diǎn)圖可形成近似多邊形，該多邊形可以幫助識(shí)別出周期性的正弦波。可選擇幾個(gè)主要的本征因子解決信號(hào)提取中平滑和信號(hào)降噪的問(wèn)題。

2.3 分離效果識(shí)別

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)性，并且可將該信息用于隨后的分組。

3 工程信號(hào)分析

3.1 信號(hào)來(lái)源

使用奇異譜分析方法對(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è)效果。

3.2 奇異譜分析結(jié)果

對(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è)特征值

3.3 奇異譜分析預(yù)測(cè)

取前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ì)比

4 結(jié)束語(yǔ)

奇異譜分析對(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)確性和使用效率。