基于SARIMA 和SVR 組合模型的轉(zhuǎn)向架系統(tǒng)壽命評(píng)估

師 蔚，范 喬，楊 洋，胡定玉，廖愛(ài)華

（1 上海工程技術(shù)大學(xué) 城市軌道交通學(xué)院，上海 201620；2 上海地鐵維護(hù)有限公司 車輛分公司，上海 200031）

隨著我國(guó)經(jīng)濟(jì)的快速發(fā)展，地鐵成為人們出行的主要交通工具，轉(zhuǎn)向架作為地鐵車輛的重要組成部分，其性能對(duì)車輛運(yùn)行的穩(wěn)定性和安全性有重要影響?；诎踩缘目紤]，對(duì)地鐵車輛轉(zhuǎn)向架系統(tǒng)進(jìn)行健康狀態(tài)和剩余壽命評(píng)估，有助于針對(duì)性地制定檢修、退役或延壽方案，對(duì)城市軌道列車運(yùn)營(yíng)的安全性與經(jīng)濟(jì)性均具有重要意義。

車輛轉(zhuǎn)向架屬于機(jī)械系統(tǒng)，目前機(jī)械系統(tǒng)的壽命評(píng)估預(yù)測(cè)方法主要有基于模型的預(yù)測(cè)和數(shù)據(jù)驅(qū)動(dòng)的預(yù)測(cè)方式兩種［1］。其中基于模型的方法是指利用模型對(duì)系統(tǒng)剩余壽命進(jìn)行預(yù)測(cè)，目前轉(zhuǎn)向架系統(tǒng)的壽命預(yù)測(cè)大多是基于動(dòng)力學(xué)模型［2-3］，如有限元分析模型對(duì)轉(zhuǎn)向架零部件的壽命預(yù)測(cè)分析，但很少有基于故障數(shù)據(jù)對(duì)轉(zhuǎn)向架整體進(jìn)行壽命評(píng)估分析。數(shù)據(jù)驅(qū)動(dòng)主要指基于故障數(shù)據(jù)采用機(jī)器學(xué)習(xí)方法進(jìn)行壽命預(yù)測(cè)，主要方法包括神經(jīng)網(wǎng)絡(luò)［4］、支持向量機(jī)等。但單個(gè)預(yù)測(cè)模型難以滿足數(shù)據(jù)多樣復(fù)雜化的條件，一般將支持向量機(jī)和其他模型或算法結(jié)合，以提高預(yù)測(cè)精度。Chen 將相似性理論和SVM 模型結(jié)合，對(duì)發(fā)動(dòng)機(jī)的壽命進(jìn)行預(yù)測(cè)，驗(yàn)證了模型的有效性［5］。劉珍翔提出了改進(jìn)增量SVR 模型對(duì)齒輪壽命進(jìn)行預(yù)測(cè)，并驗(yàn)證了該模型更可靠穩(wěn)定［6］。NIETO 等為提高預(yù)測(cè)精度，采用粒子群優(yōu)化SVM 模型對(duì)航天發(fā)動(dòng)機(jī)進(jìn)行壽命預(yù)測(cè)，解決了訓(xùn)練時(shí)超參數(shù)的優(yōu)化問(wèn)題［7］。但車輛轉(zhuǎn)向架系統(tǒng)屬于可修復(fù)系統(tǒng)，其故障數(shù)據(jù)因維修計(jì)劃存在一定的周期性，因此需結(jié)合具有周期性的預(yù)測(cè)模型，以適應(yīng)轉(zhuǎn)向架故障數(shù)據(jù)特征的需要，建立更優(yōu)的組合預(yù)測(cè)模型。

文中建立基于地鐵車輛轉(zhuǎn)向架系統(tǒng)歷史故障數(shù)據(jù)的健康指數(shù)，將具有周期性的SARIMA 模型和SVR 模型進(jìn)行組合，同時(shí)采用協(xié)方差優(yōu)選法將二者的權(quán)重進(jìn)行分配，預(yù)測(cè)轉(zhuǎn)向架系統(tǒng)的健康指數(shù)，最后將歷史和預(yù)測(cè)的健康指數(shù)結(jié)合后，建立與運(yùn)行時(shí)間的數(shù)學(xué)模型，分析得到轉(zhuǎn)向架系統(tǒng)的剩余壽命。

1 轉(zhuǎn)向架系統(tǒng)健康指數(shù)及剩余壽命評(píng)估

文中依據(jù)轉(zhuǎn)向架系統(tǒng)未來(lái)的健康狀態(tài)估計(jì)其剩余壽命，選用健康指數(shù)HI作為評(píng)價(jià)其健康狀態(tài)的指標(biāo)，HI越大表示轉(zhuǎn)向架系統(tǒng)的健康程度越高，反之健康程度越低，計(jì)算公式為式（1）［8-9］：

式中：K為比例系數(shù)；C為曲率系數(shù)；λ（t）為轉(zhuǎn)向架系統(tǒng)的故障率，故障率為式（2）：

式中：dN（t）為［t，t+dt］內(nèi)發(fā)生故障的次數(shù)；當(dāng)dt趨于0 時(shí)，λ（t）為瞬時(shí)故障率；當(dāng)dt為一段時(shí)間時(shí)，λ（t）為這段時(shí)間內(nèi)的平均故障率。

當(dāng)健康指數(shù)HI下降到轉(zhuǎn)向架允許的最低值時(shí)，表明轉(zhuǎn)向架已達(dá)到最大使用年限，不能再繼續(xù)使用，其所對(duì)應(yīng)的時(shí)間為轉(zhuǎn)向架的總壽命T0為式（3）：

式中：HIL為轉(zhuǎn)向架所允許的最低健康指數(shù)，文中設(shè)定為70。由此可得轉(zhuǎn)向架系統(tǒng)的剩余壽命ΔT為式（4）：

式中：T1為轉(zhuǎn)向架的已使用年限。

因此，根據(jù)轉(zhuǎn)向架的歷史健康指數(shù)對(duì)其未來(lái)健康指數(shù)進(jìn)行預(yù)測(cè)，并建立歷史健康指數(shù)和預(yù)測(cè)的健康指數(shù)HI與運(yùn)行時(shí)間T的數(shù)學(xué)模型，利用式（3）和式（4）估計(jì)轉(zhuǎn)向架剩余壽命。

2 基于SARIMA 和SVR 組合預(yù)測(cè)模型

文中基于轉(zhuǎn)向架歷史的健康指數(shù)預(yù)測(cè)其未來(lái)值時(shí)，轉(zhuǎn)向架系統(tǒng)屬于可修復(fù)系統(tǒng)，其故障數(shù)據(jù)因維修計(jì)劃存在一定的周期性，因此需結(jié)合具有周期性的SARIMA 預(yù)測(cè)模型，以適應(yīng)轉(zhuǎn)向架故障數(shù)據(jù)特征的需要，由于單一預(yù)測(cè)模型的精度和穩(wěn)定性不高，提出基于SARIMA 和SVR 的組合預(yù)測(cè)模型，采用協(xié)方差優(yōu)選法確定2 個(gè)模型的權(quán)重再進(jìn)行組合，建立更優(yōu)的組合預(yù)測(cè)模型。

2.1 SARIMA 預(yù)測(cè)模型

SARIMA 是常見(jiàn)的時(shí)間序列預(yù)測(cè)方法，由差分自回歸移動(dòng)平均模型演變而來(lái)，主要用于分析具有周期性變化的時(shí)間序列［10］。

SARIMA（p，d，q）（P，D，Q）S季節(jié)性回歸移動(dòng)平均模型可表示為式（5）：

式中：p、q、P、Q分別為非季節(jié)和季節(jié)性的自回歸階數(shù)、移動(dòng)平均階數(shù)；φ(B)為非季節(jié)自回歸特征多項(xiàng)式；Ф（BS）為季節(jié)自回歸特征多項(xiàng)式；θ（B）為非季節(jié)移動(dòng)平均特征多項(xiàng)式；Θ（BS）為季節(jié)移動(dòng)平均特征多項(xiàng)式；（1-B）d為d階差分；（1-BS）D為D階季節(jié)差分，二者可將時(shí)間序列Yt變成平穩(wěn)時(shí)間序列。εt為白色噪聲，服從均值為0 且方差為常數(shù)的正態(tài)分布；BS為季節(jié)后移算 子；S為周期；c為 常數(shù)；下標(biāo)t為時(shí)刻；γ1，γ2，…，γP為自回歸系數(shù)；η1，η2，…，ηQ為移動(dòng)平均系數(shù)。

2.2 SVR 預(yù)測(cè)模型

支持向量機(jī)是通過(guò)尋找一個(gè)最優(yōu)超平面，使樣本所有個(gè)體到該平面的幾何間隔最大，反之其間隔最小就可尋找一個(gè)回歸面，即支持向量回歸［11］。

SVR 回歸估計(jì)函數(shù)表達(dá)式為式（6）：

式中：w為特征空間的維數(shù)；φ(x)為從低維到高維的映射；b為偏置項(xiàng)。最優(yōu)化問(wèn)題表述為式（7）：

式中：ε為不敏感因子；ξi，為松弛變量；C′為懲罰參數(shù)。為求解式（7），引入拉格朗日乘數(shù)法得式（8）：

式中：K(xi，xj)=φ(xi)φ(xj)為核函數(shù)。文中選用徑向基（Radial Basis Function，RBF）核函數(shù)為式（9）：

式中：σ0為核寬度。

2.3 SARIMA 和SVR 組合預(yù)測(cè)

組合預(yù)測(cè)是利用不同的單個(gè)模型預(yù)測(cè)，并選取適當(dāng)?shù)臋?quán)重進(jìn)行加權(quán)，可消除單個(gè)預(yù)測(cè)模型可能存在的偏差，提高預(yù)測(cè)精度［12］，其模型組合預(yù)測(cè)值為式（10），組合預(yù)測(cè)流程如圖1 所示。

圖1 組合預(yù)測(cè)模型流程圖

式 中：YSARIMA和YSVR分別為SARIMA 預(yù)測(cè)模型和SVR 預(yù)測(cè)模型的預(yù)測(cè)結(jié)果；ω1、ω2分別為SARIMA模型和SVR 模型在組合預(yù)測(cè)中所占權(quán)重。

合理有效的權(quán)重分配，可提高組合預(yù)測(cè)的準(zhǔn)確性。文中選用協(xié)方差優(yōu)選法，根據(jù)誤差的方差來(lái)確定單個(gè)預(yù)測(cè)模型的權(quán)重，當(dāng)誤差的方差大時(shí)，其所占權(quán)重小，反之其所占權(quán)重大。令SARIMA、SVR、組合模型預(yù)測(cè)的誤差分別為e1、e2、e3。SARIMA 和SVR 預(yù)測(cè)誤差的方差分別為σ1、σ2。則組合模型預(yù)測(cè)誤差的方差為式（11）：

用拉格朗日乘數(shù)法求var（e3）極小值，在ω1+ω2=1 的條件下，得到SARIMA 模型和SVR 模型的權(quán)重為式（12）：

3 應(yīng)用分析

文中以某地鐵運(yùn)營(yíng)公司提供的轉(zhuǎn)向架系統(tǒng)故障數(shù)據(jù)為基礎(chǔ)，對(duì)健康狀態(tài)和壽命進(jìn)行評(píng)估，首先統(tǒng)計(jì)轉(zhuǎn)向架系統(tǒng)的歷史故障率，并轉(zhuǎn)化為對(duì)應(yīng)的健康指數(shù)。然后依據(jù)歷史健康指數(shù)，基于SARIMA和SVR 模型分別對(duì)健康指數(shù)進(jìn)行預(yù)測(cè)，并采用協(xié)方差優(yōu)選法得到組合模型中各自的權(quán)重，基于權(quán)重進(jìn)行組合得到組合模型預(yù)測(cè)的健康指數(shù)。最后將實(shí)際健康指數(shù)值和組合預(yù)測(cè)結(jié)果結(jié)合，建立健康指數(shù)和運(yùn)行時(shí)間的數(shù)學(xué)模型，當(dāng)健康指數(shù)下降到轉(zhuǎn)向架允許的最低值時(shí)，所對(duì)應(yīng)的時(shí)間為轉(zhuǎn)向架系統(tǒng)的總壽命，以此得到其剩余壽命。

3.1 轉(zhuǎn)向架系統(tǒng)健康指數(shù)

（1）轉(zhuǎn)向架系統(tǒng)故障率仿真

文中統(tǒng)計(jì)了轉(zhuǎn)向架系統(tǒng)從開(kāi)始運(yùn)營(yíng)到2021 年7 月的故障數(shù)據(jù)，以2 個(gè)月為時(shí)間單位，通過(guò)式（2）得到該系統(tǒng)故障率如圖2 所示，圖中5 個(gè)虛線部分分別表示第1 次架修、第1 次大修、第2 次架修、第2 次大修和第3 次架修。

圖2 轉(zhuǎn)向架系統(tǒng)故障率

（2）轉(zhuǎn)向架系統(tǒng)健康指數(shù)

通過(guò)轉(zhuǎn)向架最小故障率來(lái)獲得故障率和健康指數(shù)的轉(zhuǎn)換關(guān)系［13］，具體為：當(dāng)故障率為最小故障率λmin時(shí)，健康指數(shù)HI=100。由圖2 可知，165 個(gè)數(shù)據(jù)點(diǎn)中最小的故障率為0.006 9。采用威布爾分布分段擬合圖2 轉(zhuǎn)向架系統(tǒng)的故障率，根據(jù)轉(zhuǎn)向架30年的理論壽命，得出健康指數(shù)為70 時(shí)，故障率為0.192。即λmin=0.006 9，HI=100；λ0=0.192，HI=70；代入式（1），得比例系數(shù)K=0.006 9，曲率系數(shù)C=0.11。將求得K和C代入式（1），將圖2 的故障率轉(zhuǎn)化為健康指數(shù)，如圖3 所示。

圖3 轉(zhuǎn)向架系統(tǒng)健康指數(shù)

從圖3 看出轉(zhuǎn)向架系統(tǒng)經(jīng)過(guò)架大修后健康指數(shù)逐漸降低，健康指數(shù)隨維修計(jì)劃存在一定的周期性，因此預(yù)測(cè)其未來(lái)健康指數(shù)時(shí)，文中采用具有周期性的SARIMA 預(yù)測(cè)模型，但單一預(yù)測(cè)模型的預(yù)測(cè)精度和穩(wěn)定性不高，文中將SARIMA 和SVR模型進(jìn)行組合，同時(shí)采用協(xié)方差優(yōu)選法確定兩個(gè)模型的權(quán)重，得到更精確的組合預(yù)測(cè)模型。

3.2 轉(zhuǎn)向架健康指數(shù)預(yù)測(cè)

3.2.1 SARIMA 預(yù)測(cè)

（1）ADF 平穩(wěn)性檢驗(yàn)

SARIMA 預(yù)測(cè)的前提是要求數(shù)據(jù)為平穩(wěn)序列，通過(guò)ADF 檢驗(yàn)判定原序列的平穩(wěn)性，見(jiàn)表1。

表1 ADF 平穩(wěn)性檢驗(yàn)

ADF 檢驗(yàn)結(jié)果表明，統(tǒng)計(jì)量-0.435，大于3 個(gè)臨界值，說(shuō)明原序列不平穩(wěn)，需進(jìn)行差分運(yùn)算，一階差分后統(tǒng)計(jì)量為-20.201 8，小于3 個(gè)臨界值，則一階差分后序列平穩(wěn)，d值為1。地鐵車輛架大修周期約為4～5 年，即約為48～60 個(gè)月進(jìn)行一次架大修，上述數(shù)據(jù)均為每?jī)蓚€(gè)月統(tǒng)計(jì)1 次，因此周期S取25。同理對(duì)序列進(jìn)行周期為25 的一階季節(jié)差分后該序列是平穩(wěn)的，因此d值取1。

（2）參數(shù)值確定

AIC信息準(zhǔn)則用來(lái)確定模型參數(shù)，選取最小AIC值作為最優(yōu)模型，見(jiàn)表2，AIC=976.533 最小，因此自回歸AR 的參數(shù)p取2，移動(dòng)平均模型MA 的參數(shù)q取2 時(shí)。本數(shù)據(jù)的自相關(guān)與偏自相關(guān)系數(shù)顯著不為0，p、q取值為1。

表2 不同參數(shù)對(duì)應(yīng)AIC 值

綜上SARIMA 模型SARIMA（p，d，q）（P，D，Q）S可確定為SARIMA（2，1，2）（1，1，1）25。

（3）模型殘差檢驗(yàn)

利用SARIMA 模型預(yù)測(cè)前，需對(duì)模型殘差進(jìn)行檢驗(yàn)，如圖4 所示，圖中殘差的有序分布基本遵循正態(tài)分布采樣的線性趨勢(shì)，認(rèn)為上述模型可行。

圖4 殘差檢驗(yàn)圖

（4）預(yù)測(cè)結(jié)果

文中訓(xùn)練集為1993 年06 月～2016 年10 月轉(zhuǎn)向架的健康指數(shù)，測(cè)試集為2016 年11 月～2021 年07月數(shù)據(jù)，測(cè)試集預(yù)測(cè)結(jié)果如圖5 所示，誤差的方差為9.397。并往后進(jìn)行預(yù)測(cè)，預(yù)測(cè)結(jié)果如圖6所示。

圖5 SARIMA 測(cè)試集預(yù)測(cè)值和真實(shí)值對(duì)比圖

圖6 SARIMA 往后預(yù)測(cè)結(jié)果

3.2.2 SVR 預(yù)測(cè)

（1）數(shù)據(jù)劃分

文中采用交叉和漸消記憶原理將165 個(gè)數(shù)據(jù)按每4 個(gè)數(shù)據(jù)為一組分成162 組，其中前137 組（1993 年06 月～2016 年10 月）作為訓(xùn)練數(shù)據(jù)，后25組（2016 年11 月～2021 年07 月）作為測(cè)試數(shù)據(jù)。每組數(shù)據(jù)中前3 個(gè)為輸入樣本，最后一個(gè)為輸出樣本。

（2）核函數(shù)及參數(shù)的確定

核函數(shù)的選取決定SVR 回歸模型的精度，文中采用徑向基核函數(shù)，具有較強(qiáng)的回歸能力。同時(shí)參數(shù)的選取對(duì)SVR 預(yù)測(cè)精度也具有較大的影響，設(shè)定懲罰參數(shù)C′=12.5，不敏感因子ε=1，徑向基核函數(shù)σ0=0.03。

（3）預(yù)測(cè)結(jié)果

測(cè)試集的真實(shí)值和預(yù)測(cè)值對(duì)比如圖7 所示，測(cè)試集誤差的方差為10.029，往后預(yù)測(cè)結(jié)果如圖8所示。

圖7 SVR 測(cè)試集預(yù)測(cè)值和真實(shí)值對(duì)比圖

圖8 SVR 往后預(yù)測(cè)結(jié)果

3.2.3 組合預(yù)測(cè)

文中采用協(xié)方差優(yōu)選法確定組合模型中各自的權(quán)重，上述SARIMA 和SVR 模型預(yù)測(cè)誤差的方差σ1=9.397，σ2=10.029，代入式（12）中，得SARIMA 預(yù)測(cè)模型的權(quán)重ω1=0.52，SVR 預(yù)測(cè)模型的權(quán)重ω2=0.48。

將上述2 個(gè)模型權(quán)重代入式（10）中，得組合預(yù)測(cè)值Y=0.52YSARIMA+0.48YSVR，測(cè)試集3 種預(yù)測(cè)模型的預(yù)測(cè)值和真實(shí)值對(duì)比結(jié)果如圖9 所示。

圖9 測(cè)試集3 種模型預(yù)測(cè)和真實(shí)值對(duì)比圖

選用均方誤差（Mean Square Error，MSE）評(píng)價(jià)指標(biāo)對(duì)3 種預(yù)測(cè)模型進(jìn)行評(píng)價(jià)，公式為式（13）：

式中：n為測(cè)試集的樣本數(shù)；Ym為測(cè)試集真實(shí)值；′為測(cè)試集預(yù)測(cè)值。

3 種預(yù)測(cè)模型的均方誤差MSE比較見(jiàn)表3，SARIMA-SVR 組合預(yù)測(cè)模型相比于單個(gè)模型能有效減少誤差，提高預(yù)測(cè)精度。因此采用該組合模型往后進(jìn)行預(yù)測(cè)，結(jié)果如圖10 所示。

表3 3 種預(yù)測(cè)模型誤差分析結(jié)果

圖10 組合預(yù)測(cè)模型的預(yù)測(cè)結(jié)果

3.3 轉(zhuǎn)向架系統(tǒng)剩余壽命

由表3 可知，SARIMA 和SVR 組合預(yù)測(cè)模型的預(yù)測(cè)效果最好，因此從第3 次架修后將實(shí)際的健康指數(shù)和組合預(yù)測(cè)的結(jié)果結(jié)合，采用K鄰近法將雜點(diǎn)用前后兩個(gè)均值代替，再進(jìn)行擬合，如圖11 所示，擬合度R2=0.887 6，得到轉(zhuǎn)向架健康指數(shù)和運(yùn)行年份的函數(shù)關(guān)系為式（14）：

圖11 健康指數(shù)擬合結(jié)果

由式（14）轉(zhuǎn)化得運(yùn)行時(shí)間T為式（15）：

由上述可知HIL=70，即將y=70 代入式（15），得該車型轉(zhuǎn)向架系統(tǒng)總壽命T0為33.0 年，該車型從開(kāi)始運(yùn)行到2021 年7 月，轉(zhuǎn)向架系統(tǒng)已使用壽命T1為28.1 年，則代入式（4）得該車型轉(zhuǎn)向架系統(tǒng)剩余壽命ΔT約為4.9 年。

4 結(jié)論

（1）根據(jù)地鐵車輛轉(zhuǎn)向架系統(tǒng)歷史的故障率，通過(guò)轉(zhuǎn)化得到其對(duì)應(yīng)的健康指數(shù)，提出SARIMA和SVR 模型結(jié)合的組合模型，并采用協(xié)方差優(yōu)選法得到SARIMA 和SVR 模型各自的權(quán)重，基于組合模型對(duì)其健康指數(shù)進(jìn)行預(yù)測(cè)，結(jié)果表明組合模型的均方誤差MSE均低于單個(gè)預(yù)測(cè)模型，提高了預(yù)測(cè)精度。

（2）根據(jù)建立的歷史和預(yù)測(cè)健康指數(shù)與運(yùn)行時(shí)間的數(shù)學(xué)模型，分析其健康狀態(tài)隨運(yùn)行年限的變化規(guī)律，以此得到剩余壽命。以某地鐵車輛轉(zhuǎn)向架系統(tǒng)為例進(jìn)行算例分析，驗(yàn)證了該方法的可行性，為后期維護(hù)、維修、退役及延壽決策提供理論數(shù)據(jù)分析支撐。