張 成,韓宏宇,李 元

(1.沈陽化工大學(xué) 理學(xué)院,遼寧 沈陽 110142;2.沈陽化工大學(xué) 計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院,遼寧 沈陽 110142;3.沈陽化工大學(xué) 信息工程學(xué)院,遼寧 沈陽 110142)

0 引 言

隨著現(xiàn)代工業(yè)對產(chǎn)品質(zhì)量和生產(chǎn)安全要求的逐步提高,基于數(shù)據(jù)驅(qū)動的多元統(tǒng)計(jì)過程監(jiān)控(Multivariate Statistical Process Monitoring,MSPM)得到了廣泛的應(yīng)用[1-2]。

針對多變量監(jiān)控問題,傳統(tǒng)的主元分析(Principal Component Analysis,PCA)監(jiān)控方法可以很好地用于穩(wěn)態(tài)線性過程,或者用于監(jiān)控非線性過程的局部操作區(qū)域(近似線性)[3-4]。但是將PCA方法應(yīng)用于非線性或動態(tài)過程可能導(dǎo)致低效和不可靠的過程監(jiān)控,主要原因是線性靜態(tài)的PCA方法不適合描述過程變量中的非線性或動態(tài)特征。為了有效地監(jiān)控非線性過程,基于核主元分析(Kernel PCA,KPCA)的多元統(tǒng)計(jì)過程監(jiān)控方法被應(yīng)用于工業(yè)過程[5]。然而,KPCA并沒有考慮變量的時(shí)間相關(guān)性,即數(shù)據(jù)具有動態(tài)性時(shí),KPCA的監(jiān)控性能就受到制約。為了描述數(shù)據(jù)的動態(tài)特征,Ku等提出了一種動態(tài)主元分析(Dynamic Principal Component Analysis,DPCA)方法[6],即將PCA應(yīng)用于一個(gè)時(shí)間滯后的數(shù)據(jù)矩陣[7-8],從而將樣本之間的相關(guān)性轉(zhuǎn)化為變量之間的相關(guān)性。DPCA監(jiān)控方法也只適用于線性過程,如果將DPCA方法應(yīng)用于動態(tài)非線性過程監(jiān)控,那么該方法的監(jiān)控性能就受到制約。為了解決動態(tài)非線性過程的監(jiān)控問題,Choi等提出了一種基于動態(tài)核主元分析(Dynamic KPCA,DKPCA)[9]的過程監(jiān)控方法。DKPCA的主要思想是在執(zhí)行KPCA之前對數(shù)據(jù)矩陣進(jìn)行了時(shí)滯擴(kuò)展。雖然DKPCA能提取數(shù)據(jù)的非線性特征和動態(tài)特征,檢測結(jié)果也優(yōu)于KPCA和DPCA的結(jié)果,但是通過大量的實(shí)驗(yàn)證明,統(tǒng)計(jì)量T2上仍含有顯著的自相關(guān)性[10],這說明通過時(shí)滯擴(kuò)展并沒有徹底解決數(shù)據(jù)含有動態(tài)性的問題。

針對上述問題,該文提出一種基于去主元相關(guān)性的DKPCA動態(tài)非線性過程故障檢測與診斷方法(Dynamic Kernel Principal Component Analysis basedon Removing Principal Component Correlation,DKPCA-RPCC),其基本思想如下:首先,將低維數(shù)據(jù)引入時(shí)滯參數(shù)構(gòu)成增廣矩陣,再通過KPCA方法計(jì)算增廣矩陣的主成分,稱為真實(shí)得分;其次,假設(shè)當(dāng)前時(shí)刻樣本缺失,利用前移時(shí)刻樣本估計(jì)當(dāng)前值,構(gòu)成新的增廣矩陣,再按照相同的方法計(jì)算新矩陣的主成分,稱為預(yù)估得分;然后,使用真實(shí)得分與預(yù)估得分之間的差異來構(gòu)建統(tǒng)計(jì)量進(jìn)行故障檢測;最后,利用基于變量貢獻(xiàn)圖的方法進(jìn)行故障診斷。該方法在DKPCA的基礎(chǔ)上去除主元之間的相關(guān)性,能夠最大限度地降低樣本之間動態(tài)性的影響。

1 動態(tài)核主元分析

動態(tài)核主元分析[9]的核心思想是在應(yīng)用KPCA之前對數(shù)據(jù)矩陣進(jìn)行時(shí)滯擴(kuò)展。該方法不僅可以描述變量之間的非線性特征,而且由于加入了額外的時(shí)移變量,所以也能描述變量之間的自相關(guān)性和滯后交叉相關(guān)性。此外,通過選擇適當(dāng)時(shí)移數(shù)l,使得變量之間的非線性和動態(tài)關(guān)系都出現(xiàn)在方差較小的主元所對應(yīng)的噪聲子空間中[11]。

(1)

(2)

Cp=λp

(3)

通過φ映射之后,原始數(shù)據(jù)在高維空間的維數(shù)較大,無法直接計(jì)算特征值和特征向量,因此引入徑向基核函數(shù),通過核技巧間接計(jì)算高維空間數(shù)據(jù)集的特征值和特征向量,從而簡化計(jì)算。

(4)

(5)

(6)

在高維空間中,樣本協(xié)方差矩陣的特征向量和特征值與核矩陣的特征向量和特征值之間存在以下關(guān)系,如式(7)所示。為了簡化計(jì)算,將協(xié)方差矩陣的特征向量單位化,得到單位向量p,如式(8)所示。高維空間中數(shù)據(jù)集的得分ti可通過式(9)計(jì)算得到,其中P為協(xié)方差的特征向量矩陣,U為核矩陣的特征向量矩陣。

(7)

(8)

(9)

DKPCA的監(jiān)控指標(biāo)[9]如式(10)所示:

T2=tΛ-1tT

(10)

其中,Λ是由特征值λ組成的對角矩陣。

2 基于去主元相關(guān)性的DKPCA故障檢測和診斷方法

DKPCA方法通過時(shí)滯擴(kuò)展并沒有完全解決數(shù)據(jù)具有動態(tài)性的問題,這導(dǎo)致統(tǒng)計(jì)量仍含有較強(qiáng)的自相關(guān)性[10]。為了解決這一問題,將DKPCA與缺失值處理技術(shù)[11]相結(jié)合并利用主成分的差異構(gòu)建統(tǒng)計(jì)量。

具體方法如下:在DKPCA中,建立一個(gè)包含當(dāng)前與過去測量值的矩陣Y。對于每個(gè)新的觀測向量yi可以計(jì)算得分。假設(shè)當(dāng)前的觀測向量缺失,可以使用缺失值處理技術(shù),從過去的數(shù)據(jù)中估計(jì)出當(dāng)前時(shí)刻的測量值并計(jì)算得分,這本質(zhì)上是對得分和觀察值進(jìn)行了一步預(yù)測。因此,得分之間的差異幾乎是連續(xù)非相關(guān)的[8],這意味著可以使用得分差異更好地對樣本進(jìn)行監(jiān)控。所選擇的缺失值處理技術(shù)為條件均值替換法[12],該方法對缺失數(shù)據(jù)的測量矢量進(jìn)行重新排列,不失一般性。

(11)

式中,Y#表示缺失的測量值,Y*表示已知的測量值。

(12)

假設(shè)增廣矩陣已知。為了計(jì)算缺失向量的觀測值,只需要EM算法的期望步驟,把P代入到S的表達(dá)式,可以得到式(13)。使用式(13)可以計(jì)算缺失觀測值的條件期望,如式(14)所示。

(13)

(14)

最終可以得到一個(gè)新的矩陣,如式(15)所示。

(15)

M=K1U

(16)

(17)

(18)

(1)離線建模。

第一步:對訓(xùn)練數(shù)據(jù)集X進(jìn)行標(biāo)準(zhǔn)化。

第二步:根據(jù)時(shí)移l生成增廣矩陣Yx,然后通過KPCA計(jì)算真實(shí)得分Mx。

第四步:根據(jù)式(18)計(jì)算統(tǒng)計(jì)量W,使用核密度估計(jì)方法[15]確定控制限WUCL。

(2)在線檢測。

第一步:對測試樣本應(yīng)用訓(xùn)練數(shù)據(jù)的均值和方差進(jìn)行標(biāo)準(zhǔn)化。

第二步:根據(jù)時(shí)移l生成增廣向量Yfi,使用訓(xùn)練數(shù)據(jù)的核矩陣進(jìn)行中心化,求出測試數(shù)據(jù)的真實(shí)得分Mfi。

第四步:根據(jù)式(18)計(jì)算測試樣本的Wi,若Wi>WUCL,則判定Yfi為故障樣本,否則Yfi為正常樣本。

(3)故障診斷。

(19)

3 仿真實(shí)驗(yàn)

3.1 數(shù)值例子

使用一個(gè)非線性動態(tài)數(shù)值例子[9]測試文中方法的有效性,具體模型如式(20)所示。

(20)

其中,u(t)、y(t)、z(t)分別是輸入變量、輸出變量以及狀態(tài)變量,o(t)是服從均勻分布U～(0,0.5)的隨機(jī)噪聲,w(t)是服從均勻分布U～(0,5)的白噪聲,系統(tǒng)A、B、C、D、w(t)、o(t)分別為:

其中,非線性映射函數(shù)f(u(t))如式(21)所示:

(21)

由以上模型隨機(jī)生成1 000個(gè)正常數(shù)據(jù)作為訓(xùn)練樣本,然后,在變量w1加上一個(gè)幅值為-4.5的階躍變化作為故障數(shù)據(jù)。再由該模型生成100組測試數(shù)據(jù)集,每組測試數(shù)據(jù)前500個(gè)為正常樣本,后500個(gè)為故障樣本。輸出變量u(t)、y(t)被用于過程監(jiān)控。采用DPCA、KPCA、DKPCA以及文中方法對此動態(tài)非線性過程例子進(jìn)行仿真實(shí)驗(yàn),并對上述的四種方法的結(jié)果進(jìn)行分析。

采用DPCA對該數(shù)值例子進(jìn)行故障檢測時(shí),按照95%的累計(jì)方差貢獻(xiàn)率選取主元,因此主元數(shù)設(shè)置為5。按照平行分析法[17],時(shí)滯參數(shù)的選取與特征值有關(guān),特征值為0或接近于0時(shí)表示變量之間線性無關(guān),因而時(shí)滯參數(shù)l的值等于特征值為0或者接近于0的個(gè)數(shù)[6]。DPCA提取的主元之間含有較強(qiáng)的相關(guān)性,特別是前兩個(gè)主元之間的相關(guān)性最大,這將導(dǎo)致計(jì)算的統(tǒng)計(jì)量T2也含有較強(qiáng)的自相關(guān)性,嚴(yán)重影響故障的檢測結(jié)果。

采用KPCA對此數(shù)值例子進(jìn)行故障檢測時(shí),通過經(jīng)驗(yàn)法[18]確定核寬參數(shù)的范圍以及選取最優(yōu)核寬參數(shù)。低維空間的數(shù)據(jù)集通過非線性函數(shù)映射,使得數(shù)據(jù)在高維空間線性可分,此時(shí)故障樣本與正常樣本之間的差異最大,檢測結(jié)果達(dá)到最優(yōu)。采用DKPCA對此數(shù)值例子進(jìn)行故障檢測時(shí),核寬參數(shù)σ和時(shí)滯參數(shù)l采用經(jīng)典的網(wǎng)格搜索和交叉驗(yàn)證[19]方式進(jìn)行選取,該方法選取不同的參數(shù)時(shí),訓(xùn)練數(shù)據(jù)的準(zhǔn)確率如圖1所示。

圖1 校驗(yàn)數(shù)據(jù)的準(zhǔn)確率

采用文中方法進(jìn)行故障檢測時(shí),核寬參數(shù)與時(shí)滯參數(shù)的選取均與DKPCA保持一致。圖2給出了上述四種方法對于100組測試數(shù)據(jù)的檢測結(jié)果,由圖2可以看出,文中方法的檢測結(jié)果明顯高于其他三種方法且平均檢測效率最高。圖3表示四種方法統(tǒng)計(jì)值自相關(guān)性的對比。由圖3可以看出,DPCA、KPCA和DKPCA方法的統(tǒng)計(jì)值都含有較強(qiáng)的一步自相關(guān)性,而文中方法的統(tǒng)計(jì)量的自相關(guān)性明顯降低,因此能夠提高故障檢測率。圖4給出了文中方法對第一組數(shù)據(jù)的故障檢測結(jié)果,由圖4可知,故障樣本與正常樣本幾乎完全分離,證明了文中方法的有效性。使用公式(19)對測試數(shù)據(jù)進(jìn)行故障診斷,測試數(shù)據(jù)集10個(gè)變量對樣本的貢獻(xiàn)程度,如圖5所示,從圖5可以看出第4、9變量的貢獻(xiàn)度較大,所以故障發(fā)生與這兩個(gè)變量的相關(guān)性最大。

圖2 四種方法故障檢測結(jié)果

圖3 統(tǒng)計(jì)量的自相關(guān)性

圖4 DKPCA-RPCC故障檢測結(jié)果

圖5 樣本變量的累計(jì)貢獻(xiàn)

3.2 Tennessee Eastman

本節(jié)的仿真數(shù)據(jù)由最新的TE仿真器生成[20],其基本結(jié)構(gòu)如圖6所示。數(shù)據(jù)集在TE仿真器中通過0.01小時(shí)采樣間隔并持續(xù)運(yùn)行10小時(shí)獲得。正常狀態(tài)下采集1 000個(gè)樣本作為訓(xùn)練集;同時(shí),過程故障在3小時(shí)后引入并持續(xù)到過程結(jié)束。TE過程共包含41個(gè)測量變量和12個(gè)控制變量,仿真器運(yùn)行過程中有三個(gè)變量恒定不變,所以選取剩余的50個(gè)變量進(jìn)行過程監(jiān)控。

圖6 TE過程

仿真使用的計(jì)算機(jī)環(huán)境為Intel(R)Core(TM) i7-7700HQCPU@ 2.80 GHz處理器的64位操作系統(tǒng),且所有實(shí)驗(yàn)結(jié)果均在MATLAB中實(shí)現(xiàn)。

利用文中方法對TE數(shù)據(jù)進(jìn)行仿真驗(yàn)證,并與KPCA、DPCA、DKPCA進(jìn)行對比分析。四種方法均按照95%的累計(jì)貢獻(xiàn)率選取主元,此外DPCA的時(shí)滯參數(shù)選取為1,DKPCA和文中方法的核寬參數(shù)和時(shí)滯參數(shù)均通過網(wǎng)格搜索和交叉驗(yàn)證的方式選取。表1給出了上述各種方法的檢測結(jié)果,由于第六組數(shù)據(jù)缺失,因此不做測試。以第27組測試數(shù)據(jù)為例進(jìn)行故障檢測,圖7分別給出DKPCA-RPCC、KPCA、DPCA和DKPCA對于故障27的檢測結(jié)果。對比四種方法的檢測結(jié)果可以看出,文中方法的檢測率最高。為了進(jìn)一步確定故障產(chǎn)生的原因,使用式(19)對第27組數(shù)據(jù)進(jìn)行故障診斷,圖8給出了故障27的診斷結(jié)果。從圖中可以看出,第9個(gè)和第49個(gè)變量的貢獻(xiàn)度最大,即故障27是由變量9和49引起的。變量9和49分別代表反應(yīng)堆溫度和反應(yīng)堆冷卻水流量,圖9給出了反應(yīng)堆溫度和反應(yīng)堆冷卻水流量的變化趨勢,由圖可以看出與計(jì)算的結(jié)果恰好相互印證,與實(shí)際相符。

圖7 故障27的檢測結(jié)果

圖8 故障27的貢獻(xiàn)圖

圖10表示四種方法統(tǒng)計(jì)量的自相關(guān)性,明顯可以看出,DKPCA-RPCC方法構(gòu)建的統(tǒng)計(jì)量具有較低水平的自相關(guān)性。根據(jù)故障檢測率以及構(gòu)建統(tǒng)計(jì)量的自相關(guān)性檢驗(yàn)所提出方法的性能。統(tǒng)計(jì)量的自相關(guān)性越低,說明主成分之間含有的動態(tài)性越低,故障檢測性能更好。

圖10 四種方法統(tǒng)計(jì)量的自相關(guān)

4 結(jié)束語

提出了一種基于去主元相關(guān)性的DKPCA動態(tài)非線性過程故障檢測與診斷方法。該方法相對傳統(tǒng)的DKPCA方法有兩點(diǎn)改進(jìn),一是通過缺失值處理技術(shù)預(yù)測與原始矩陣保留相同特征的矩陣,二是通過真實(shí)得分與預(yù)估得分作差消除相同特征并在差異空間構(gòu)造新的統(tǒng)計(jì)量。該方法進(jìn)一步降低了動態(tài)性的影響,使得統(tǒng)計(jì)量的自相關(guān)性顯著降低。另外,給出了基于變量累計(jì)貢獻(xiàn)率的故障診斷方法。通過數(shù)值例子和TE過程驗(yàn)證該方法在故障檢測與診斷中具有更好效果,對于動態(tài)非線性監(jiān)控過程具有一定的指導(dǎo)意義。

由于該方法在計(jì)算過程中多次使用核函數(shù),因而出現(xiàn)檢測過程運(yùn)行時(shí)間長,占用內(nèi)存空間大的問題,接下來將考慮如何優(yōu)化算法的問題,提高算法的運(yùn)行效率。