結(jié)合小波去噪的PCA工業(yè)過(guò)程傳感器故障檢測(cè)

（武漢大學(xué) 動(dòng)力與機(jī)械學(xué)院，武漢 430072）

主元分析PCA是一種多變量分析技術(shù)，在數(shù)據(jù)壓縮、模式識(shí)別、信號(hào)處理和過(guò)程監(jiān)控等領(lǐng)域均獲得廣泛運(yùn)用[1-3]。PCA將原始的高維度數(shù)據(jù)投影到低維空間，以消除各個(gè)變量之間的弱相關(guān)性，降低數(shù)據(jù)的復(fù)雜度；從而能夠更容易地從海量歷史數(shù)據(jù)中提取有用信息[4]。PCA將數(shù)據(jù)空間分解為2個(gè)正交子空間：主元空間PCS（principal component sub space）和殘差空間 RS（residual subspace），利用主元模型求取SPE、T2統(tǒng)計(jì)量[5]。由于工業(yè)現(xiàn)場(chǎng)存在噪聲且噪聲不滿足高斯分布，使得PCA在工業(yè)現(xiàn)場(chǎng)傳感器故障檢測(cè)中存在一定程度的故障誤檢測(cè)問(wèn)題。為此，如何削弱工業(yè)現(xiàn)場(chǎng)噪聲的影響，降低故障誤檢測(cè)率值得進(jìn)行研究。

小波分析是信號(hào)處理的有力工具，可保證在信號(hào)品質(zhì)不受影響的前提下濾去信號(hào)中的噪聲[6-8]。1984 年，文獻(xiàn)[9]提出了“小波”（wavelet）的概念，提出用一個(gè)函數(shù)的時(shí)移性和尺度組合表示信號(hào)的新思想。1989年，文獻(xiàn)[10]提出信號(hào)的多分辨率分析與快速重構(gòu)算法，隨后小波變換快速走向?qū)嵱没?。小波去噪在工業(yè)過(guò)程也獲得了廣泛運(yùn)用，文獻(xiàn)[11]通過(guò)滑動(dòng)窗將小波分析與PCA結(jié)合，用于TE化工過(guò)程故障檢測(cè)，論證了該方法的可行性。

本文結(jié)合小波去噪和PCA，通過(guò)滑動(dòng)窗對(duì)PCA的SPE、T2統(tǒng)計(jì)量進(jìn)行小波去噪，削弱噪聲對(duì)SPE、T2統(tǒng)計(jì)量的影響。文中提出的方法在電廠水處理過(guò)程中進(jìn)行了仿真，結(jié)果表明該方法可有效減少故障誤檢測(cè)。

1 PCA故障誤檢測(cè)

設(shè)工業(yè)生產(chǎn)過(guò)程有m個(gè)傳感器，PCA故障檢測(cè)包含以下3個(gè)步驟：

步驟1數(shù)據(jù)收集和標(biāo)準(zhǔn)化。構(gòu)造測(cè)量數(shù)據(jù)矩陣 X1=［x1，x2，…，xn］T∈Rn×m，其中 X1的每一列表示一個(gè)變量，每一行表示一個(gè)樣本。對(duì)X1的每一個(gè)樣本 xi（i∈［1，n］）進(jìn)行標(biāo)準(zhǔn)化處理，得到新樣本集 X。

步驟2特征值分解及SPE和Hotelling T2控制限求取。X按照式（1）分解：

步驟3在線計(jì)算SPE和Hotelling T2統(tǒng)計(jì)量并檢測(cè)故障。對(duì)實(shí)測(cè)樣本x∈R1×m（已經(jīng)標(biāo)準(zhǔn)化處理），計(jì)算其SPE和Hotelling T2統(tǒng)計(jì)量。按照故障檢測(cè)邏輯診斷當(dāng)前系統(tǒng)運(yùn)行狀態(tài)：當(dāng)SPE＞SPElim或時(shí)，則判定故障產(chǎn)生，否則判定無(wú)故障。有關(guān)PCA故障檢測(cè)可參考文獻(xiàn)[12-14]。

PCA故障檢測(cè)在工業(yè)現(xiàn)場(chǎng)運(yùn)用時(shí)，SPE、T2統(tǒng)計(jì)量波動(dòng)劇烈，容易超過(guò)控制限導(dǎo)致故障誤檢測(cè)。為了降低故障誤檢測(cè)率，本文利用小波方法濾除SPE、T2統(tǒng)計(jì)量的噪聲。

2 小波PCA故障檢測(cè)

2.1 小波變換與去噪

連續(xù)小波變換CWT（continuous wavelet transform）的數(shù)學(xué)表達(dá)式為

式中：a，τ分別為尺度因子和平移因子；ψ（t）為母函數(shù)；ψ*（t）為其共軛函數(shù)。尺度因子a控制函數(shù)的伸縮：當(dāng)a增大時(shí)，函數(shù)時(shí)窗變寬；當(dāng)a減小時(shí)，函數(shù)時(shí)窗變窄。這樣小波變換就可以對(duì)信號(hào)進(jìn)行多尺度細(xì)化，同時(shí)得到良好的時(shí)間、頻率分辨率。

式（2）中尺度因子a連續(xù)變換，導(dǎo)致計(jì)算量大，實(shí)際工程應(yīng)用中，常用離散小波變換DWT（discrete wavelet transform）。能在保證信號(hào)不失真的同時(shí)，減小計(jì)算量。文獻(xiàn)[10]中提出了正交小波快速算法，得到了廣泛運(yùn)用。本文也采用這種算法進(jìn)行小波分解，得到各尺度的小波系數(shù)和近似系數(shù)，然后通過(guò)離散小波逆變換IDWT（inverse discrete wavelet transform）重構(gòu)信號(hào)。

以信號(hào)f（t）的小波三層分解為例，小波分解如圖1所示。

圖1 信號(hào)的小波三層分解Fig.1 Three level decomposition of signal

圖1中：H0，H1分別為低通、高通濾波器；↓2代表下采樣過(guò)程。信號(hào)f（t）三層分解后得：

文獻(xiàn)[6-8]提出了非線性小波變換閾值去噪法，得到了廣泛運(yùn)用。小波變換可將信號(hào)的“能量”集中于少數(shù)小波系數(shù)；而白噪聲由于其無(wú)序性，其“能量”分散于各小波系數(shù)。隨著分解的進(jìn)行，白噪聲的小波系數(shù)變小，信號(hào)的小波系數(shù)遠(yuǎn)大于白噪聲的小波系數(shù)。因此，可找到合適的閾值δ，將小于δ的小波系數(shù)濾去，得到純凈的小波系數(shù)。最后將去噪后的小波系數(shù)進(jìn)行重構(gòu)即可獲得濾波后的信號(hào)。

Matlab小波工具箱提供了小波去噪函數(shù)，本文即采用Matlab工具箱實(shí)現(xiàn)信號(hào)去噪[15]。

2.2 結(jié)合小波去噪的PCA故障檢測(cè)

2.2.1 小波PCA故障檢測(cè)流程

傳統(tǒng)PCA故障檢測(cè)方法在實(shí)際工程應(yīng)用中SPE和T2統(tǒng)計(jì)量波動(dòng)劇烈，易導(dǎo)致故障誤檢測(cè)。本文采用小波去噪濾除SPE和T2統(tǒng)計(jì)量的噪聲，以減少PCA故障誤檢測(cè)率。結(jié)合小波去噪的PCA故障檢測(cè)流程如下：

①PCA離線建模。利用歷史正常工況下的數(shù)據(jù)得到負(fù)載矩陣Pa、得分矩陣T、主元個(gè)數(shù)a，以及SPElim和控制限。

②PCA統(tǒng)計(jì)量在線計(jì)算。對(duì)于實(shí)測(cè)數(shù)據(jù)樣本x∈R1×m，計(jì)算其對(duì)應(yīng)的 SPE 和 T2統(tǒng)計(jì)量。

③小波去噪。通過(guò)滑動(dòng)窗對(duì)SPE和T2統(tǒng)計(jì)量進(jìn)行小波去噪。構(gòu)建長(zhǎng)度為wlen的SPE和T2統(tǒng)計(jì)量滑動(dòng)窗，根據(jù)小波函數(shù)、分解層數(shù)lev、閾值、閾值處理方式和閾值重調(diào)規(guī)則進(jìn)行小波去噪，獲得去噪后的統(tǒng)計(jì)量。

2.2.2 小波去噪?yún)?shù)的選取

以上步驟中，小波去噪是關(guān)鍵環(huán)節(jié)，對(duì)故障誤檢測(cè)率影響大。故需要研究小波去噪的參數(shù)選取問(wèn)題。

①小波函數(shù)選?。涸谛〔ǚ治鲋嘘懤m(xù)出現(xiàn)了一系列經(jīng)典小波系，例如Haar小波、Daubechies小波系、Symlets小波系、Coiflets小波系和Spline小波系等[9]。其中，正交小波系可根據(jù)快速算法得到信號(hào)完全重構(gòu)，計(jì)算量較小。故本文采用正交小波系作為基小波。正交小波系有 Haar小波、Daubechies、Symlets和Coiflets小波系。

Haar小波是Daubechies小波系中第一個(gè)小波函數(shù)。Daubechies小波系和Symlets小波系分別有45種小波。小波序號(hào)越大，小波函數(shù)計(jì)算量越大，本文選取了Daubechies小波系和Symlets小波系的前15 個(gè)小波，即 db1，…，db15 和 sym1，…，sym15。

Coiflet小波系包含5種小波函數(shù)，它們所需的計(jì)算量不大。在后續(xù)敏感性分析中選取Coiflet小波系全部5個(gè)小波函數(shù)，即coif1，…，coif5進(jìn)行對(duì)比分析。

②閾值處理方式。小波去噪有2種通用的閾值處理方式：軟閾值法和硬閾值法。設(shè)WT為小波系數(shù)，δ為閾值，軟、硬閾值法可表示為

i）硬閾值法（hard）

ii）軟閾值法（soft）

硬閾值法處理方式簡(jiǎn)單，但會(huì)造成處理后小波系數(shù)的不連續(xù)，引起重構(gòu)信號(hào)振蕩；而軟閾值法處理后小波系數(shù)的連續(xù)性好，在工業(yè)上運(yùn)用更廣泛。

③閾值選取規(guī)則。Matlab工具箱中，有4種閾值選取規(guī)則[15]：固定式閾值（sqtwolog）、Stein 無(wú)偏估計(jì)閾值（rigrsure）、啟發(fā)式閾值（heursure）和極大極小值閾值（minimax）。不同的閾值對(duì)去噪效果造成影響。

④閾值重調(diào)方法：Matlab小波工具箱里，有3種閾值重調(diào)方法：one、sln、mln。它們規(guī)定了各層小波系數(shù)采用的閾值是否重調(diào)的規(guī)則。

⑤分解層數(shù)。分解層數(shù)lev受滑動(dòng)窗長(zhǎng)度wlen影響。例如wlen=8，由于下采樣過(guò)程使小波系數(shù)樣本數(shù)減半，則lev最大為3?；瑒?dòng)窗長(zhǎng)度wlen限制了分解層數(shù)lev的大小。

⑥滑動(dòng)窗參數(shù)?；瑒?dòng)窗參數(shù)包括滑動(dòng)窗長(zhǎng)度wlen和滑動(dòng)窗步長(zhǎng)wstep。本文選取滑動(dòng)窗長(zhǎng)度wlen為2的正整數(shù)倍，以避免邊界失真，保證信號(hào)重建的準(zhǔn)確性[16]?；瑒?dòng)窗步長(zhǎng)wstep應(yīng)小于或等于wlen，還需考慮計(jì)算量因素。

3 應(yīng)用研究

3.1 熱力發(fā)電廠水處理流程

本文以1000 MW熱力發(fā)電廠補(bǔ)給水處理流程為對(duì)象驗(yàn)證本文提出的方法。水處理工藝如圖2所示。原水分別經(jīng)預(yù)處理——陽(yáng)床/陰床——混床，最終出水。陽(yáng)床、陰床、混床中的混合離子樹(shù)脂置換水中的鹽，達(dá)到凈化水質(zhì)的目的。水處理流程配備了大量的傳感器，對(duì)水質(zhì)指標(biāo)進(jìn)行實(shí)時(shí)監(jiān)測(cè)[17]。本文從水處理流程選取了12個(gè)主要傳感器作為故障檢測(cè)對(duì)象，如表1所示。

圖2 陽(yáng)-陰-混床水處理流程Fig.2 Flow chart of cation-anion-mixed bed

表1 傳感器列表Tab.1 List of sensors

3.2 小波PCA故障誤檢測(cè)方法驗(yàn)證

首先從電廠SIS系統(tǒng)中選取正常、穩(wěn)定的運(yùn)行數(shù)據(jù)為學(xué)習(xí)樣本，采樣時(shí)間為5 s，樣本容量100。以該段數(shù)據(jù)為訓(xùn)練樣本建模主元模型。再選擇同一工況下另一組穩(wěn)態(tài)數(shù)據(jù)為檢測(cè)樣本集，采樣時(shí)間5 s，樣本容量1000。分別采用傳統(tǒng)PCA和結(jié)合小波去噪的PCA進(jìn)行水處理過(guò)程故障檢測(cè)，2種方法的驗(yàn)證結(jié)果分別如圖3和圖4所示。

圖4中小波去噪?yún)?shù)設(shè)置如下：db10小波函數(shù)、軟閾值、固定式閾值、mln閾值重調(diào)、wlen=64、wstep=6、lev=3）。

圖3 傳統(tǒng)PCA SPE、T2控制Fig.3 SPE，T2statistics of classical PCA method

圖4 結(jié)合小波去噪PCA SPE、T2控制Fig.4 SPE，T2statistics of PCA combined with wavelet denoising

為了便于書寫，本文將小波去噪?yún)?shù)寫為以下形式：（db10，soft，sqtwolog，mln，wlen=64，wstep=6，lev=3）。

圖3表明傳統(tǒng)PCA方法檢測(cè)樣本SPE、T2控制量存在“毛刺”，且SPE、T2的峰值超出了對(duì)應(yīng)的控制限，導(dǎo)致故障誤報(bào)。雖然故障誤判的持續(xù)時(shí)間短，僅為1～5個(gè)采樣時(shí)刻，但這些短時(shí)間的故障誤判影響了PCA故障檢測(cè)方法的實(shí)際使用效果。圖4表明用小波去噪后，PCA SPE、T2控制量在保留信號(hào)趨勢(shì)的同時(shí)，濾除了毛刺，有效解決了故障誤報(bào)的問(wèn)題。

傳統(tǒng)PCA方法與結(jié)合小波去噪PCA方法的故障誤報(bào)性能如表2所示。

表2 傳統(tǒng)PCA與結(jié)合小波去噪PCA性能統(tǒng)計(jì)Tab.2 Comparison between classical PCA and PCA combined with wavelet denoising

表2表明，通過(guò)對(duì)PCA SPE、T2統(tǒng)計(jì)量進(jìn)行小波去噪，能有效消除故障誤報(bào)，使PCA方法在工業(yè)生產(chǎn)過(guò)程的運(yùn)用更加實(shí)用。圖4為一種去噪?yún)?shù)組合下的特例分析，為了更加徹底地理解小波參數(shù)的影響，以下將對(duì)小波去噪?yún)?shù)的敏感性進(jìn)行分析，以指導(dǎo)PCA小波去噪?yún)?shù)的選取。

3.3 小波去噪?yún)?shù)的敏感性研究

小波去噪是影響PCA故障誤報(bào)的關(guān)鍵環(huán)節(jié)。小波去噪?yún)?shù)多，對(duì)去噪性能都有影響，故需要對(duì)小波去噪?yún)?shù)的敏感性進(jìn)行分析。本文小波去噪最終服務(wù)于降低故障誤報(bào)率，故將系統(tǒng)的故障誤報(bào)次數(shù)作為小波去噪效果的評(píng)判標(biāo)準(zhǔn)。若某一小波參數(shù)降低了故障誤報(bào)次數(shù)，判斷該參數(shù)組合優(yōu)化；反之，則判定參數(shù)組合劣化。

3.3.1 小波系

首先對(duì)小波系進(jìn)行分析，在其他參數(shù)固定時(shí)分別選用 Symlets、Daubechies和Coiflets小波系的小波函數(shù)進(jìn)行對(duì)比研究。各小波系小波函數(shù)濾波效果如圖 5 所示， 其他參數(shù)設(shè)置為（soft，sqtwolog，mln，wlen=64，lev=3，）。圖5中將誤報(bào)次數(shù)作為優(yōu)劣評(píng)判指標(biāo)。

圖5 各小波系濾波效果對(duì)比Fig.5 Comparison among wavelet species

圖5表明，隨著小波函數(shù)序號(hào)的增大，各小波系都能減少誤報(bào)次數(shù)。其中sym2、sym10和db2對(duì)系統(tǒng)性能的提升作用較差，Symlets小波系濾波效果不穩(wěn)定，Coiflet小波系對(duì)性能敏感性較低。Daubechies小波系中小波函數(shù)序號(hào)大于8的小波都能消除誤報(bào)，比Symlets和Coiflet小波系穩(wěn)定；且Daubechies小波系比Symlets小波系計(jì)算量小。故Daubechies小波系，特別是小波函數(shù)序號(hào)大于8的小波函數(shù)具有優(yōu)勢(shì)。以下參數(shù)分析就將重點(diǎn)針對(duì)Daubechies小波系開(kāi)展。

3.3.2 閾值處理方式

當(dāng)其他參數(shù)設(shè)置為（Daubechies小波系，sqtwolog，mln，wlen=64，wstep=4，lev=4），軟、硬閾值濾波效果對(duì)比如圖6所示。

圖6 軟硬閾值濾波效果對(duì)比Fig.6 Comparison between soft and hard threshold

圖6表明，對(duì)Daubechies小波系的各個(gè)小波，軟閾值對(duì)性能的提升都遠(yuǎn)好于硬閾值，且軟閾值濾波效果更加穩(wěn)定。

3.3.3 閾值選取規(guī)則

取其他小波參數(shù)為（Daubechies小波系，soft，mln，wlen=64，wstep=4，lev=4），閾值選取規(guī)則濾波效果對(duì)比如圖7所示。

圖7 各閾值濾波效果對(duì)比Fig.7 Comparison among threshold selection rule

圖7表明，對(duì)Daubechies小波系的各個(gè)小波，sqtwolog閾值濾波效果良好、穩(wěn)定，具有較大優(yōu)勢(shì)。

3.3.4 閾值重調(diào)

取其他小波參數(shù)為（Daubechies小波系，soft，sqtwolog，wlen=64，wstep=4，lev=4），閾值重調(diào)濾波效果對(duì)比如圖8所示。

圖8 閾值重調(diào)濾波效果對(duì)比Fig.8 Comparison among threshold rescaling method

圖8表明，對(duì)Daubechies小波系的各個(gè)小波，mln參數(shù)濾波效果最好。這是由于水處理過(guò)程現(xiàn)場(chǎng)存在非高斯白噪聲，在每一層都需要進(jìn)行閾值重調(diào)。

3.3.5 分解層數(shù)

當(dāng)其他小波參數(shù)選取為（Daubechies小波系，soft，sqtwolog，wlen=64，wstep=4，lev=4），分解層數(shù)濾波效果對(duì)比如圖9所示。

圖9 分解層數(shù)濾波效果對(duì)比Fig.9 Comparison among decomposition level

圖9表明，隨著lev的上升，系統(tǒng)故障誤報(bào)數(shù)明顯下降，濾波效果顯著提升。當(dāng)lev=4時(shí)，已滿足系統(tǒng)正常運(yùn)行要求，繼續(xù)增加lev對(duì)性能的提升不敏感，且增加了計(jì)算量。綜合考慮計(jì)算量與濾波效果，建議lev為4或者5。

3.3.6 滑動(dòng)窗

滑動(dòng)窗長(zhǎng)度wlen與滑動(dòng)窗步長(zhǎng)wstep對(duì)濾波效果無(wú)明顯的影響。分解層數(shù)lev受wlen的限制，當(dāng)lev=4時(shí)，系統(tǒng)故障誤報(bào)率滿足系統(tǒng)要求，綜合考慮建議wlen為64或者128，滑動(dòng)窗步長(zhǎng)wstep為4～10。

4 結(jié)語(yǔ)

傳統(tǒng)PCA在工業(yè)現(xiàn)場(chǎng)故障診斷中易出現(xiàn)故障誤檢測(cè)問(wèn)題，本文將小波去噪應(yīng)用于PCA，旨在提升故障檢測(cè)的性能。以某1000 MW熱力發(fā)電廠補(bǔ)給水處理流程為對(duì)象，以該流程典型工況下運(yùn)行數(shù)據(jù)為基礎(chǔ)驗(yàn)證了小波PCA故障檢測(cè)算法，結(jié)果表明結(jié)合小波去噪的PCA方法能有效減少系統(tǒng)誤報(bào)次數(shù)。

本文還對(duì)小波去噪?yún)?shù)進(jìn)行了敏感性分析，對(duì)小波濾波參數(shù)的設(shè)置提出了建議。建議使用Daubechies小波系小波序號(hào)大于8的小波函數(shù)，采用軟閾值、sqtwolog、mln、4或 5層分解，滑動(dòng)窗長(zhǎng)度為64或128，滑動(dòng)窗步長(zhǎng)為4～10。本文提出了利用小波去噪降低故障誤報(bào)率的方法和步驟，提升了PCA故障檢測(cè)的性能，使其更加實(shí)用化。

[1]張大海，劉宇穗，張世榮，等.發(fā)電廠水處理流程傳感器故障檢測(cè)系統(tǒng)研究[J].自動(dòng)化與儀表，2014，29（5）：9-13.

[2]Chen G，Qian S E.Denoising of hyperspectral imagery using principal component analys is and waveletshrinkage[J].Geoscience and Remote Sensing，IEEE Transactionson，2011，49（3）：973-980.

[3]王健，馮健，韓志艷.基于流形學(xué)習(xí)的局部保持PCA算法在故障檢測(cè)中的應(yīng)用[J].控制與決策，2013（5）：683-687.

[4]周東華，李鋼，李元.數(shù)據(jù)驅(qū)動(dòng)的工業(yè)過(guò)程故障診斷技術(shù)：基于主元分析與偏最小二乘的方法[M].北京：科學(xué)出版社，2011.

[5]Shams B，Budman H M，Duever T A.Fault detection，identification and diagnosis using CUSUM based PCA[J].Chemical Engineering Science，2011，66（20）：4488-4498.

[6]Donoho D L.Adapting to unknown smoothnessvia wavelet shrinkage[J].J Amer Statist Assoc，1995，90：1200-1224.

[7]Donoho D L，Johnstone I.Wavelet shrinkage asymptopia[J].Journal of Royal Statistical Society，1995，57（2）：301-369.

[8]Donoho D L.Denoising by soft-thresholding[J].IEEE Transaction on Information，1995（3）：613-627.

[9]Grossmann A，Morlet J.Decomposition of hardy functions into square integrable wavelets of constant shape[J].SIAM journal on mathematical analysis，1984，15（4）：723-736.

[10]Mallat S，Hwang W L.Singularity detection and processing with wavelets[J].Information Theory，IEEE Transactions on，1992，38（2）：617-643.

[11]Li S，Wen J.A model-based fault detection and diagnostic methodology based on PCA method and wavelet transform[J].Energy and Buildings，2014，68：63-71.

[12]牛征，劉吉臻，牛玉廣.動(dòng)態(tài)多主元模型故障檢測(cè)方法在變工況過(guò)程中的應(yīng)用[J].動(dòng)力工程，2005，25（3）：426-426.

[13]Ding S，Zhang P，Ding E，et al.On the application of PCA technique to fault diagnosis[J].Tsinghua Science& Technology，2010，15（2）：138-144.

[14]Jeng J C.Adaptive process monitoring using efficient recursive PCA and moving window PCA algorithms[J].Journal of the Taiwan Institute of Chemical Engineers，2010，41（4）：475-481.

[15]Michel M，Yves M，Georges O，et al.Wavelet toolbox for use with MATLAB[J].Wavelet Toolbox User’s Guilde，1997.

[16]Strang G，Nguyen T.Wavelets and Filter Banks[M].SIAM，1996.

[17]Faust S D，Aly O M.Adsorption Processes for Water Treatment[M].Elsevier，2013.