張曉玲, 曹玉蘋(píng), 鄧曉剛
(1.中國(guó)石油大學(xué) 勝利學(xué)院, 山東 東營(yíng) 257097;2.中國(guó)石油大學(xué)(華東) 控制科學(xué)與工程學(xué)院, 山東 青島 266580)
間歇生產(chǎn)方式在現(xiàn)代工業(yè)界所占比重越來(lái)越高,其過(guò)程監(jiān)控和故障診斷也引起廣泛關(guān)注[1]。長(zhǎng)期以來(lái),過(guò)程監(jiān)控的主流統(tǒng)計(jì)方法[2-5]如主元分析(principal component analysis, PCA)、偏最小二乘(partial least squares, PLS)、獨(dú)立元分析(independent component analysis, ICA)、典型相關(guān)分析(canonical correlation analysis, CCA)等在應(yīng)用于間歇過(guò)程故障診斷時(shí),都是采用將三維多批次歷史數(shù)據(jù)展開(kāi)成二維的多向(multi-way)處理方法,即使針對(duì)間歇過(guò)程的不同特點(diǎn)涌現(xiàn)出一系列改進(jìn)算法[6-14],也大都對(duì)數(shù)據(jù)進(jìn)行展開(kāi)處理,如按批次展開(kāi)、按變量展開(kāi)成二維矩陣的方式等。此類(lèi)多向展開(kāi)的處理方法將三維多批次數(shù)據(jù)中每一批次的二維矩陣展開(kāi)變成向量看待,勢(shì)必會(huì)破壞原有數(shù)據(jù)的相關(guān)性,導(dǎo)致部分信息的丟失和特征不完全提取。
近年來(lái),一些基于張量的數(shù)據(jù)降維方法,如張量PCA (tensor PCA, TPCA)[15-16]、張量鄰域保持嵌入算法(tensor neighborhood preserving embedding, TNPE)[17]等,使用張量模式的特征提取方法對(duì)張量形式的數(shù)據(jù)進(jìn)行特征提取,保證了張量形式的特征空間不被破壞,張量分析逐漸在計(jì)算機(jī)視覺(jué)、圖像處理、模式識(shí)別等領(lǐng)域成為研究熱點(diǎn)[18]。本文借鑒以上方法中張量分析的思想,對(duì)間歇過(guò)程在線監(jiān)控時(shí),將三維數(shù)據(jù)的每一批次看作一個(gè)二階張量,應(yīng)用TPCA算法對(duì)三維數(shù)據(jù)直接提取特征信息,不破壞三維樣本中每一批次的原始矩陣空間結(jié)構(gòu),并定義R2統(tǒng)計(jì)量和SPE統(tǒng)計(jì)量建立故障檢測(cè)模型;另外,為能更及時(shí)檢測(cè)出間歇過(guò)程中常出現(xiàn)的漸變故障(如間歇精餾過(guò)程中催化劑活性的緩慢變化、間歇反應(yīng)釜換熱器由于長(zhǎng)時(shí)間使用結(jié)垢現(xiàn)象越來(lái)越嚴(yán)重而導(dǎo)致?lián)Q熱性能的逐漸下降),本文利用累積和控制圖(cumulative sum,CUSUM)分析思想,將其應(yīng)用于 TPCA算法中,提出累加和的張量主元分析(summed tensor principal component analysis, STPCA),將三維歷史數(shù)據(jù)按一定的滑動(dòng)窗口累加求和,累積歷史信息,再做TPCA算法處理,建立故障檢測(cè)模型,從而提高對(duì)漸變故障在初始出現(xiàn)階段檢測(cè)的靈敏度。
可通過(guò)以下迭代方法實(shí)現(xiàn)式(1)中U、V的求解:
Step2:設(shè)置循環(huán)變量n初始值為1,進(jìn)入循環(huán)迭代運(yùn)算。
Step4:對(duì)CY1進(jìn)行特征值分解得到特征值矩陣Λ和特征向量矩陣Vn+1,只選取前J'個(gè)最大的特征值組成,Vn+1則是由前J'個(gè)特征值對(duì)應(yīng)的特征向量矩陣。
Step5:令Y2m=XmVn,求Y2的協(xié)方差陣
Step8:輸出最終的U∈RK×K'、V∈RJ×J'。
圖1 三維數(shù)據(jù)預(yù)處理Fig.1 Schematic diagram of three-way dataset preprocessing
注意,TPCA變換處理的數(shù)據(jù)X是經(jīng)過(guò)標(biāo)準(zhǔn)化處理的,其處理過(guò)程如圖1所示,第1步將初始數(shù)據(jù)展開(kāi)成大小為M×KJ的二維矩陣,沿著每列進(jìn)行均值為0方差為1的歸一化處理,第2步再恢復(fù)成三維結(jié)構(gòu)數(shù)據(jù)X(M×K×J)。此處的展開(kāi)僅用于數(shù)據(jù)的預(yù)處理,應(yīng)用TPCA變換時(shí)使用的仍是三維矩陣X(M×K×J)。因此,預(yù)處理后的三維數(shù)據(jù)X(M×K×J)通過(guò)以上Step1~8的迭代算法獲得U和V后,經(jīng)過(guò)映射公式Y(jié)m=UTXmV,則可被投影到低維張量特征空間得到Y(jié)(M×K'×J' ) = {Ym∈K'×J',m= 1,2,… ,M}。
各類(lèi)間歇過(guò)程中的漸變故障很常見(jiàn),漸變故障是某變量操作值緩慢偏離正常變化區(qū)間、不斷累積的過(guò)程,故障初期變化量較小,很多故障檢測(cè)方法對(duì)該類(lèi)故障存在預(yù)報(bào)大大延遲的問(wèn)題。為提高檢測(cè)故障靈敏度,本文利用PAGE等[19]提出的累積和(CUSUM)控制圖思想,提出一種基于累加和的張量主元分析算法(STPCA),將經(jīng)過(guò)預(yù)處理后的正常批次數(shù)據(jù)X(M×K×J)按照每批次內(nèi)滑動(dòng)窗口大小進(jìn)行歷史信息的疊加,數(shù)據(jù)處理公式如下:
summed_xm(k)表示第m批次中按照式(2)累加和后的第k采樣時(shí)刻觀測(cè)值,h為累加求和的寬度,構(gòu)成移動(dòng)的滑動(dòng)窗口,當(dāng)k>h,每一時(shí)刻的數(shù)據(jù)值都由前h個(gè)多變量樣本求和得到,從而得到求和后的summed_X(M×K×J) = { summed_Xm∈RK?RJ,m= 1 ,2,…,M},然后用第1節(jié)中TPCA算法對(duì)其進(jìn)行特征提取,求得變換矩陣U和V,為下一步故障診斷打下基礎(chǔ)。
在線監(jiān)控過(guò)程主要包含2個(gè)階段。
第1階段:離線模型的建立
② 對(duì)X(M×K×J)按式(2)進(jìn)行累加和處理,獲得summed_X(M×K×J)。
③ 運(yùn)用 TPCA算法,確定投影矩陣U∈RK×K'和V∈RJ×J',計(jì)算經(jīng)過(guò)投影后的低維張量特征空間數(shù)據(jù)Y={Ym∈RK'?RJ',m= 1,2,… ,M} ,Ym是累加和處理后summed_Xm的低維投影,變換公式為:Ym=UT* s ummed_Xm*V。
④ 建立故障檢測(cè)模型。
對(duì)基于TPCA和STPCA 2種算法的間歇過(guò)程故障檢測(cè),本文都采用 R2和SPE統(tǒng)計(jì)量。
R2是基于支持向量數(shù)據(jù)域描述(support vector domain description,SVDD)[20]方法建立的統(tǒng)計(jì)量,表示每批次Xm的投影Ym到特征空間中心的距離,其計(jì)算式如下[21]:
拉格朗日乘子iα通過(guò)下式的最大化問(wèn)題求解[21]:
C為懲罰因子,可取1。
SPE統(tǒng)計(jì)量用來(lái)測(cè)量殘差空間的變化。先通過(guò)U、V重構(gòu)summed_Xm的估計(jì)值,得到殘差矩陣從而得SPE統(tǒng)計(jì)量的計(jì)算式:
R2和SPE統(tǒng)計(jì)量皆可通過(guò)核密度估計(jì)法[22]確定各自95% 或99% 的置信限(也稱(chēng)控制限、閾值),在每個(gè)采樣時(shí)刻k分別利用核密度估計(jì)法計(jì)算對(duì)應(yīng)采樣點(diǎn)所有批次統(tǒng)計(jì)量的閾值,因此可得到動(dòng)態(tài)變化的控制限,當(dāng)統(tǒng)計(jì)量超出各自控制限所判定的正常變化范圍,則意味著故障的出現(xiàn)。
第2階段:在線故障檢測(cè)
①采集新批次第k采樣時(shí)刻數(shù)據(jù)根據(jù)離線建模標(biāo)準(zhǔn)化處理時(shí)得到的均值和方差進(jìn)行預(yù)處理。
②對(duì)預(yù)處理后的數(shù)據(jù)xnew,k按照式(2)進(jìn)行累加求和處理,得到 s ummed_xnew,k。
③計(jì)算新的投影矩陣和殘差:
其中uk表示矩陣U中的第k行向量。
④參考式(3)和(5)分別計(jì)算新的在線2R和SPE統(tǒng)計(jì)量,得到
并判斷統(tǒng)計(jì)量是否偏離各自閾值,以確定間歇過(guò)程運(yùn)行狀況。
分批補(bǔ)料的盤(pán)尼西林發(fā)酵過(guò)程是典型的間歇式生產(chǎn)過(guò)程,本文基于BIROL等[23]開(kāi)發(fā)的標(biāo)準(zhǔn)仿真軟件Pensim 2.0產(chǎn)生的盤(pán)尼西林發(fā)酵數(shù)據(jù)進(jìn)行在線監(jiān)控算法仿真研究。本仿真由該軟件在正常操作范圍內(nèi)對(duì)初始條件加入一定的微小偏移和適當(dāng)噪聲生成30批次數(shù)據(jù)作為訓(xùn)練樣本,每批次含11個(gè)觀測(cè)變量,反應(yīng)時(shí)間相等,均為400 h,采樣間隔為0.5 h,從而得到用于建立MPCA、TPCA、STPCA監(jiān)控模型。另外,利用Pensim2.0產(chǎn)生6個(gè)故障批次,表1是3種跳變故障,表2是漸變類(lèi)型故障。
表1 跳變故障Table 1 Step faults used for detection
表2 漸變故障Table 2 Gradual faults used for detection
建模時(shí),3種方法都需要對(duì)數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化,不同的是,MPCA方法如圖1中Step 1按批次展開(kāi)成二維,按列作歸一化處理,得到二維矩陣X(30×8800),即X(M×KJ),在此基礎(chǔ)上建立MPCA模型。TPCA和STPCA算法則需將標(biāo)準(zhǔn)化后的二維矩陣X(30×8800)再按圖1中Step 2還原成X(30× 8 00×11),在保持和不破壞原來(lái)每批次數(shù)據(jù)結(jié)構(gòu)基礎(chǔ)上建立各自的監(jiān)控模型。因此,傳統(tǒng)PCA算法會(huì)把每批次的數(shù)據(jù)Xm∈RK?RJ轉(zhuǎn)換為向量xm∈RKJ看待,特征值分解的矩陣大小是KJ×KJ,基于張量的 TPCA方法其特征值分解的矩陣大小僅為K×K或J×J,在特征值分解計(jì)算復(fù)雜度上要更低。另外,STPCA相比TPCA方法在數(shù)據(jù)處理上多一步累加和處理,以增大對(duì)漸變故障檢測(cè)的靈敏度。最后的監(jiān)控模型所用到統(tǒng)計(jì)量也有差異:MPCA使用T2和SPE,TPCA和STPCA則運(yùn)用本文定義的 R2和SPE統(tǒng)計(jì)量。
圖2 基于MPCA方法的故障1檢測(cè)結(jié)果Fig.2 Fault 1 detection results based on MPCA
圖3 基于TPCA方法的故障1檢測(cè)結(jié)果Fig.3 Fault 1 detection results based on TPCA
圖4 基于MPCA方法的故障2檢測(cè)結(jié)果Fig.4 Fault 2 detection results based on MPCA
當(dāng)新批次中出現(xiàn)故障 1、2、3跳變故障時(shí),基于 MPCA的在線故障檢測(cè)結(jié)果如圖2、4、6所示,MPCA的SPE統(tǒng)計(jì)量能準(zhǔn)確預(yù)報(bào)故障的發(fā)生(140~300采樣點(diǎn)),但在非故障時(shí)段有較高誤報(bào),從圖中可明顯看到有多點(diǎn)超出控制限,特別是對(duì)故障3達(dá)到了5.12% 的誤報(bào)率(見(jiàn)表3)。相對(duì)而言,圖3、5、7所示TPCA故障檢測(cè)結(jié)果中,SPE統(tǒng)計(jì)圖既能準(zhǔn)確檢出故障,又具有較低誤報(bào)率,由表3可見(jiàn),對(duì)前兩個(gè)故障的誤報(bào)率為0,僅對(duì)故障3有約1.63% 的誤報(bào),但仍低于MPCA方法的結(jié)果。與SPE統(tǒng)計(jì)量相比,MPCA方法的T2統(tǒng)計(jì)圖檢測(cè)效果較差,只能對(duì)故障1有效,對(duì)故障2和3完全不敏感,無(wú)任何超限警報(bào);反觀TPCA方法的 R2統(tǒng)計(jì)圖,在采樣點(diǎn)140~300僅無(wú)法檢測(cè)出故障3,能完全準(zhǔn)確檢測(cè)到故障2,對(duì)故障1僅有少量點(diǎn)的漏報(bào)。從誤報(bào)率情況分析,TPCA的 R2統(tǒng)計(jì)量也優(yōu)于MPCA的T2統(tǒng)計(jì)量。總體分析,TPCA兩種統(tǒng)計(jì)圖對(duì)3個(gè)故障的檢出率高于 MPCA方法,且具有較低的誤報(bào)率??梢?jiàn),將利用張量分解直接處理三維數(shù)據(jù)的TPCA方法引入到間歇過(guò)程故障檢測(cè)中,構(gòu)造 R2和SPE統(tǒng)計(jì)量進(jìn)行實(shí)時(shí)監(jiān)控是可行并具有一定效果的。
表3 MPCA和TPCA方法的故障檢測(cè)誤報(bào)率對(duì)比Table 3 Comparison of false positive rates (FPR) by MPCA and TPCA
圖5 基于TPCA方法的故障2檢測(cè)結(jié)果Fig.5 Fault 2 detection results based on TPCA
圖6 基于MPCA方法的故障3檢測(cè)結(jié)果Fig.6 Fault 3 detection results based on MPCA
圖7 基于TPCA方法的故障3檢測(cè)結(jié)果Fig.7 Fault 3 detection results based on TPCA
為提高檢測(cè)漸變故障的靈敏度,本文在TPCA基礎(chǔ)上提出基于STPCA的間歇過(guò)程故障診斷方法,同樣使用 R2和SPE統(tǒng)計(jì)圖進(jìn)行故障分析。仿真實(shí)驗(yàn)設(shè)置了表2所述3種故障,分別是通風(fēng)速率、攪拌功率、底物流加速率從第200采樣點(diǎn)以一定速度呈現(xiàn)斜坡變化,模擬間歇過(guò)程中可能出現(xiàn)的漸變類(lèi)故障。圖8到圖13反映了MPCA、STPCA方法對(duì)3個(gè)故障批次的檢測(cè)效果。
從圖8、10和12所示結(jié)果觀察,MPCA的T2和SPE統(tǒng)計(jì)量均在一定時(shí)刻越過(guò)閾值不斷增加,表明漸變型故障的發(fā)生;同樣,圖9、11、13中STPCA方法的統(tǒng)計(jì)結(jié)果圖也給出了故障警報(bào)。但明顯的是,兩種方法的SPE統(tǒng)計(jì)圖比基于得分的T2和基于特征空間的 R2統(tǒng)計(jì)量趨于更靈敏,給出故障預(yù)報(bào)要早。為直觀對(duì)比不同方法的檢測(cè)效果優(yōu)劣,本文通過(guò)表4列出兩種方法各統(tǒng)計(jì)量所確定的故障出現(xiàn)時(shí)刻,并假設(shè)統(tǒng)計(jì)量連續(xù)8個(gè)采樣點(diǎn)超過(guò)控制限則認(rèn)定該采樣時(shí)刻為故障檢測(cè)點(diǎn)。對(duì)比SPE統(tǒng)計(jì)量,針對(duì)漸變故障4、5、6,STPCA分別領(lǐng)先MPCA方法9個(gè)、28個(gè)、30個(gè)采樣點(diǎn)提前給出故障指示;STPCA的2R統(tǒng)計(jì)圖檢測(cè)故障的能力雖遜于SPE統(tǒng)計(jì)量,卻仍比MPCA的2T統(tǒng)計(jì)量分別早117、336、56個(gè)采樣時(shí)刻超出控制限。綜合看來(lái),兩種方法檢測(cè)漸變故障時(shí),相對(duì)故障實(shí)際發(fā)生時(shí)刻都有或多或少的滯后,原因在于故障變量剛開(kāi)始是緩慢變化的,偏離正常范圍需要一定的時(shí)間間隔,只有當(dāng)變量值漸漸漂移累積超出一定量后,監(jiān)測(cè)系統(tǒng)才能判斷過(guò)程出現(xiàn)了非正常的變化。但STPCA方法將歷史批次數(shù)據(jù)進(jìn)行疊加處理,放大故障變量的初始偏移,并在張量空間進(jìn)行特征值分解、不破壞原有數(shù)據(jù)結(jié)構(gòu)基礎(chǔ)上得到監(jiān)控統(tǒng)計(jì)量,結(jié)合以上結(jié)果圖和表4,不難看出該方法能更早檢測(cè)出漸變故障,體現(xiàn)該方法的優(yōu)越性。
圖8 基于MPCA方法的故障4檢測(cè)結(jié)果Fig.8 Fault 4 detection results based on MPCA
圖9 基于STPCA方法的故障4檢測(cè)結(jié)果Fig.9 Fault 4 detection results based on STPCA
圖10 基于MPCA方法的故障5檢測(cè)結(jié)果Fig.10 Fault 5 detection results based on MPCA
圖11 基于STPCA方法的故障5檢測(cè)結(jié)果Fig.11 Fault 5 detection results based on STPCA
圖12 基于MPCA方法的故障6檢測(cè)結(jié)果Fig.12 Fault 6 detection results based on MPCA
圖13 基于STPCA方法的故障6檢測(cè)結(jié)果Fig.13 Fault 6 detection results based on STPCA
表4 MPCA和STPCA檢出的故障發(fā)生采樣點(diǎn)Table 4 Fault sampling points based on MPCA and STPCA
需要說(shuō)明的是,STPCA方法更適用于漸變故障存在的場(chǎng)合,尤其是故障出現(xiàn)早期,通過(guò)CUSUM累積過(guò)程中的微小變化,提高對(duì)漸進(jìn)、緩慢變化類(lèi)故障初期較小故障偏移的檢測(cè)。同時(shí),該方法中累加求和寬度h的取值與檢測(cè)效果關(guān)系密切,若h取值太小,利用 CUSUM 累積的歷史信息不夠豐富,對(duì)漸變故障出現(xiàn)初期檢測(cè)準(zhǔn)確率會(huì)較低,反之一味追求快速檢測(cè)出故障而將h取得過(guò)大,會(huì)因帶入部分噪聲信息而增加誤報(bào)率,導(dǎo)致故障檢測(cè)效果下降。目前總體上,h大小的選擇沒(méi)有統(tǒng)一標(biāo)準(zhǔn),多是根據(jù)經(jīng)驗(yàn)或交叉驗(yàn)證方法確定。本文通過(guò)反復(fù)試驗(yàn),選擇適當(dāng)?shù)膆值,既能提高對(duì)漸變故障的檢出率,較早發(fā)出故障報(bào)警,又防止出現(xiàn)較高誤報(bào)的情況。另外,STPCA算法在故障檢測(cè)時(shí),對(duì)新批次數(shù)據(jù)引入CUSUM的處理,不可避免會(huì)產(chǎn)生故障檢測(cè)延時(shí),因此,本方法并不適用于故障檢測(cè)實(shí)時(shí)性要求較高的場(chǎng)合。
傳統(tǒng)將三維數(shù)據(jù)變成二維的間歇過(guò)程故障診斷方法,破壞數(shù)據(jù)原有的空間結(jié)構(gòu)和相關(guān)特性,必然會(huì)對(duì)過(guò)程監(jiān)控效果產(chǎn)生不利影響。本文將基于張量分解的TPCA算法引入到間歇過(guò)程監(jiān)控中,直接將三維數(shù)據(jù)作為張量目標(biāo)進(jìn)行運(yùn)算,保留了數(shù)據(jù)完整性,定義2R 和SPE統(tǒng)計(jì)量建立故障檢測(cè)模型。經(jīng)仿真驗(yàn)證,基于TPCA的間歇過(guò)程故障診斷方法適合檢測(cè)突變故障,檢測(cè)效果較好。針對(duì)間歇過(guò)程中常發(fā)生的漸變故障,本文提出STPCA方法,通過(guò)對(duì)三維數(shù)據(jù)多個(gè)連續(xù)時(shí)刻的采樣值求和,不斷放大歷史信息,同時(shí)又保證不破壞三維原始數(shù)據(jù)結(jié)構(gòu),因此在檢測(cè)漸變類(lèi)型故障的仿真實(shí)驗(yàn)中,STPCA方法有更好效果,提高故障檢測(cè)靈敏度,相應(yīng)降低了漏檢率。
符號(hào)說(shuō)明: