基于VBEM的ARFA模型參數(shù)推導(dǎo)和故障檢測(cè)

李吉俊 章智杰 李其操 董自健

摘? 要：文章使用自回歸因子分析模型（ARFA）對(duì)數(shù)據(jù)樣本進(jìn)行動(dòng)態(tài)過(guò)程建模，分析了卡爾曼濾波和EM算法在估計(jì)ARFA模型中回歸矩陣參數(shù)A和載荷矩陣參數(shù)C的方法。在此基礎(chǔ)上，提出了一種使用變分貝葉斯EM（VBEM）故障檢測(cè)方法，對(duì)ARFA模型參數(shù)A和C進(jìn)行推斷和動(dòng)態(tài)過(guò)程故障檢測(cè)。仿真實(shí)驗(yàn)結(jié)果表明，在ARFA模型下，VBEM方法對(duì)下文所述的階躍信號(hào)、斜坡信號(hào)等四類故障的檢測(cè)效果要優(yōu)于EM方法對(duì)該類故障的檢測(cè)效果，并且降低了平均迭代次數(shù)。

關(guān)鍵詞：自回歸因子分析模型;EM算法;變分貝葉斯EM;故障檢測(cè)

中圖分類號(hào)：TP277;TP391.9? ? ? 文獻(xiàn)標(biāo)識(shí)碼：A 文章編號(hào)：2096-4706（2020）24-0001-06

Parameter Derivation and Fault Detection of ARFA Model Based on VBEM

LI Jijun，ZHANG Zhijie，LI Qicao，DONG Zijian

（School of Electronic Engineering，Jiangsu Ocean University，Lianyungang? 222005，China）

Abstract：In this paper，the autoregressive factor analysis（ARFA）model is used to model the dynamic process of data samples，and the methods of Kalman filter and EM algorithm to estimate the regression matrix parameter A and load matrix parameter C in ARFA model are analyzed. On this basis，a fault detection method based on variational Bayes EM（VBEM）is proposed to infer ARFA model parameters A and C and detect dynamic process faults. The simulation results show that under ARFA model，the detection effect of VBEM method is better than that of EM method for four kinds of faults，such as step signal and slope signal，and the average number of iterations is reduced.

Keywords：autoregressive factor analysis model;EM algorithm;variational Bayes EM;fault detection

0? 引? 言

隨著科學(xué)和技術(shù)的發(fā)展，在信息化、數(shù)據(jù)化的今天，工業(yè)系統(tǒng)已逐漸呈現(xiàn)出大規(guī)模、強(qiáng)復(fù)雜性、高集成化等特點(diǎn)，對(duì)生產(chǎn)過(guò)程的安全、平穩(wěn)、產(chǎn)品質(zhì)量等提出了更高的要求。一旦復(fù)雜工程系統(tǒng)發(fā)生故障而又未能得到有效處理，不僅可能導(dǎo)致財(cái)產(chǎn)損失或人員傷亡，還可能導(dǎo)致環(huán)境污染造成生態(tài)災(zāi)難[1]。因此，常常采用有效的工業(yè)過(guò)程故障檢測(cè)與診斷技術(shù)，來(lái)提高復(fù)雜工業(yè)過(guò)程控制系統(tǒng)的可靠性和安全性，及時(shí)有效地發(fā)現(xiàn)和檢測(cè)出系統(tǒng)中存在的故障，以此來(lái)避免由故障造成的人員傷亡和經(jīng)濟(jì)損失[2]。傳統(tǒng)的主元分析方法（Principal Com-ponent Analysis，PCA）是一種廣泛應(yīng)用于過(guò)程監(jiān)控的線性方法，通過(guò)將原始高維空間的數(shù)據(jù)投影到低維空間中來(lái)消除數(shù)據(jù)間的相關(guān)性。然而在實(shí)際工業(yè)過(guò)程短采樣中，還存在著數(shù)據(jù)的自相關(guān)性[3]。為解決這一問(wèn)題，很多學(xué)者提出了多變量動(dòng)態(tài)過(guò)程建模和故障檢測(cè)的方法。KU[4]等人在1995年提出動(dòng)態(tài)PCA（DPCA）算法，該方法通過(guò)把當(dāng)前樣本的前L個(gè)樣本加入到數(shù)據(jù)矩陣中，構(gòu)成時(shí)間滯后的數(shù)據(jù)增廣矩陣，來(lái)提取時(shí)序相關(guān)的關(guān)系，從而消除數(shù)據(jù)的自相關(guān)性。LI[5]等提出了一種間接DPCA模型用于動(dòng)態(tài)過(guò)程監(jiān)測(cè)，該方法可以有效提取數(shù)據(jù)的自相關(guān)性，卻忽略了數(shù)據(jù)的靜態(tài)特性。

KIM[6]等提出了一種基于概率主成分分析（PPCA）的多變量過(guò)程監(jiān)控方法。并指出當(dāng)概率生成模型應(yīng)用于過(guò)程監(jiān)控時(shí)，大多數(shù)的統(tǒng)計(jì)量如單變量Shewhart圖、多變量Shewhart圖、T2和SPE統(tǒng)計(jì)量都可以被統(tǒng)一到概率模型的框架中。楊沛武[7]等提出了一種基于動(dòng)態(tài)PPCA（Dynamic PPCA）的故障檢測(cè)技術(shù)，通過(guò)增廣時(shí)間序列變量來(lái)擴(kuò)大數(shù)據(jù)矩陣，進(jìn)而將數(shù)據(jù)的自相關(guān)性納入相應(yīng)的增廣矩陣中來(lái)進(jìn)行DPPCA建模。然而該方法無(wú)法全部獲取數(shù)據(jù)的互相關(guān)特征。LI[8]等提出了一種新的動(dòng)態(tài)隱變量模型（DLV）。在DLV中，采用AR模型創(chuàng)建了一個(gè)新的殘差，并通過(guò)PCA分解來(lái)提取剩余部分的靜態(tài)變化。而? 和? 統(tǒng)計(jì)量則分別用來(lái)監(jiān)測(cè)動(dòng)態(tài)和靜態(tài)變化。因此，DLV能夠同時(shí)提取數(shù)據(jù)的自相關(guān)性和互相關(guān)性。WEN[9]等提出了一種基于線性高斯?fàn)顟B(tài)空間模型（LGSSM）的動(dòng)態(tài)過(guò)程監(jiān)測(cè)方法。LGSSM通過(guò)一階馬爾科夫特性來(lái)描述過(guò)程的動(dòng)態(tài)關(guān)系，并對(duì)過(guò)程數(shù)據(jù)做了降維處理。雖然該方法也能夠同時(shí)提取數(shù)據(jù)的自相關(guān)性和互相關(guān)性，但對(duì)高階動(dòng)態(tài)過(guò)程的特征提取能力較差。周樂(lè)[10]提出了一種自回歸因子分析模型（ARFA），它具有普適性的動(dòng)態(tài)模型，能夠同時(shí)提取數(shù)據(jù)間的動(dòng)態(tài)和靜態(tài)關(guān)系，有效提髙動(dòng)態(tài)過(guò)程建模和故障檢測(cè)效果。

本文在ARFA的等價(jià)模型和卡爾曼濾波的基礎(chǔ)上，提出了一種新的變分貝葉斯EM的動(dòng)態(tài)過(guò)程故障檢測(cè)方法。該方法推導(dǎo)模型回歸矩陣參數(shù)A和載荷矩陣參數(shù)C的期望和故障檢測(cè)的總體思路是：由參數(shù)和數(shù)據(jù)的聯(lián)合分布，根據(jù)變分理論求出A和C的邊緣分布，再求出參數(shù)A和C的后驗(yàn)期望，并把這兩個(gè)參數(shù)期望代入到ARFA的標(biāo)準(zhǔn)化模型和卡爾曼濾波方程中進(jìn)行迭代，直到模型參數(shù)A和C收斂為止。把收斂后的模型參數(shù)代入T2和SPE統(tǒng)計(jì)量的相關(guān)公式中進(jìn)行計(jì)算。當(dāng)樣本的兩個(gè)統(tǒng)計(jì)量均在控制限以下時(shí)，則認(rèn)為沒(méi)有故障發(fā)生;均在控制限以上時(shí)，則認(rèn)為有故障發(fā)生[3]。

本文安排為，第一節(jié)介紹ARFA模型及其EM算法;第二節(jié)為利用變分貝葉斯EM方法對(duì)ARFA模型進(jìn)行推斷;第三節(jié)為仿真實(shí)驗(yàn);第四節(jié)為結(jié)論。

1? ARFA模型及其EM算法

1.1? ARFA模型及其等價(jià)系統(tǒng)

ARFA模型形式為：

fk=Azk-1+wk

jk=Cfk+vk? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?（1）

其中jk∈RM為測(cè)量變量，fk∈RD為動(dòng)態(tài)因子，zk-1= [fk-1T，fk-2T，…，fk-LT]T∈RDL包含了動(dòng)態(tài)因子的前L個(gè)值。A∈RD×DL為回歸矩陣，C∈RM×D為載荷矩陣。系統(tǒng)的動(dòng)態(tài)噪聲和測(cè)量噪聲分別為wk和vk，它們分別服從均值為零，方差為Q和R的高斯分布[10]。綜上所述，ARFA的模型有效的參數(shù)為{A，C，Q，R}。

為了利用卡爾曼濾波來(lái)估計(jì)系統(tǒng)狀態(tài)，ARFA模型可以等價(jià)為下述線性系統(tǒng)[10]：

（2）

系統(tǒng)（2）的定義與ARFA相似，因此，我們可以得到兩者之間的關(guān)系為和? 。相應(yīng)的動(dòng)態(tài)噪聲，其中，

1.2? 利用EM算法來(lái)估計(jì)ARFA模型參數(shù)A和C

K個(gè)輸出序列J=（j1，j2，…，jk）T和它們相應(yīng)的動(dòng)態(tài)因子F=（f1，f2，…，fk）T的聯(lián)合對(duì)數(shù)似然函數(shù)為：

（3）

在M步中，利用對(duì)數(shù)極大似然函數(shù)分別對(duì)A和C求一階偏導(dǎo)數(shù)，得到的模型參數(shù)更新值為：

（4）

（5）

其中，矩陣、和? 分別由式（6）、式（7）、式（8）給出：

（6）

（7）

（8）

1.3? 基于ARFA模型的動(dòng)態(tài)過(guò)程故障檢測(cè)

我們構(gòu)建T2和SPE兩種統(tǒng)計(jì)量來(lái)檢測(cè)動(dòng)態(tài)過(guò)程是否發(fā)生故障。與傳統(tǒng)的PCA構(gòu)建T2統(tǒng)計(jì)量的方法不同，ARFA中動(dòng)態(tài)因子fk由于包含了數(shù)據(jù)的動(dòng)態(tài)關(guān)系而不再相互獨(dú)立，因此，fk已經(jīng)不再適用于直接構(gòu)建統(tǒng)計(jì)量以檢測(cè)隱變量空間的變化[10]。為此，我們利用動(dòng)態(tài)因子的估計(jì)殘差來(lái)構(gòu)建T2統(tǒng)計(jì)量為：

（9）

（10）

其中， 由式（10）給出：

（11）

此外，根據(jù)模型的預(yù)測(cè)誤差，我們可以構(gòu)造SPE統(tǒng)計(jì)量為[10]：

（12）

SPE=ejkTejk? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? （13）

根據(jù)正常樣本的統(tǒng)計(jì)量，我們可以估計(jì)故障檢測(cè)的統(tǒng)計(jì)限。T2統(tǒng)計(jì)量的控制限可以由χ2分布估計(jì)為[11]：

T2=T2lim=χ2（D）? ? ? ? ? ? ? ? ? ? ? ? （14）

其中，α為顯著性水平，D為動(dòng)態(tài)因子的維度。

SPE統(tǒng)計(jì)量的控制限同樣由近似的χ2分布估計(jì)為[12]：

SPE=SPElim～g·χ2（h）? ? ? ? ? ? ? ? ?（15）

g和h滿足以下方程，其中，α為顯著性水平，g為常數(shù)，h為自由度。

gh=mean（SPEnormal）? ? ? ? ? ? ? ? ? ?（16）

2g2h=var（SPEnormal）? ? ? ? ? ? ? ? ? ?（17）

2? 利用變分貝葉斯EM方法對(duì)ARFA模型進(jìn)行推斷

參數(shù)和數(shù)據(jù)的聯(lián)合分布為：

P（A，C，fk，jk）=P（jk|A，C，fk）·P（A，C，fk）（18）

變分近似推導(dǎo)的目的就是為了尋找一個(gè)可以逼近待求后驗(yàn)分布、可以處理的函數(shù)。根據(jù)變分理論[13，14]有：

P（A，C，fk）≈q（A，C，fk）=q（A）q（C）q（fk）（19）

A的邊緣分布為[15]：

（20）

其中，a（j）為A矩陣的第j列， 為協(xié)方差矩陣，SA，（j）為數(shù)量矩陣SA的第j列。

C的邊緣分布為[15]：

（21）

其中，c（j）為C矩陣的第j列， 為協(xié)方差矩陣，SC，（j）為數(shù)量矩陣SC的第j列。

假定A和C矩陣中的超參數(shù)α和γ均如服從以下形式的伽馬分布：

（22）

（23）

相關(guān)數(shù)據(jù)統(tǒng)計(jì)矩陣的表達(dá)式由式（24）、式（25）、式（26）、式（27）給出：

（24）

（25）

（26）

（27）

定義和的逆矩陣由式（28）和式（29）給出[15]：

（28）

（29）

qθ（θ）關(guān)于參數(shù)向量θ=（A，C）的變分后驗(yàn)的期望由式（30）和式（31）給出：

（30）

（31）

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

本文在處理器為CPU Intel i5-5200U 2.20 GHz的筆記本中，采用Matlab（2016b版本）編程語(yǔ)言來(lái)對(duì)一個(gè)二階動(dòng)態(tài)系統(tǒng)的數(shù)值例子進(jìn)行仿真實(shí)驗(yàn)。數(shù)值例子為：

mk=Azk-1+wk

nk=Cmk+vk

其中，nk∈R4為k時(shí)刻的觀測(cè)值，mk∈R2是其相應(yīng)的動(dòng)態(tài)因子，zk-1=[mk-1T，mk-2T]T∈R4包含了動(dòng)態(tài)因子的前兩個(gè)值。系統(tǒng)的動(dòng)態(tài)噪聲和測(cè)量噪聲分別為wk和vk，它們均服從均值為0，方差為0.01的高斯分布?；貧w矩陣和載荷矩陣分別為：

我們從標(biāo)準(zhǔn)正態(tài)分布中隨機(jī)取出100個(gè)2維數(shù)據(jù)作為正常樣本。為了進(jìn)行故障檢測(cè)實(shí)驗(yàn)，我們產(chǎn)生了四種測(cè)試集，每個(gè)測(cè)試集當(dāng)中仍然包括100個(gè)樣本點(diǎn)。這些測(cè)試集分別引入了四種不同的常見(jiàn)故障，它們分別是：故障1：在第51～100個(gè)樣本點(diǎn)上引入幅度為8的階躍信號(hào)故障;故障2：在第51～100個(gè)樣本點(diǎn)上引入幅度為1.6×（50-t）的斜坡信號(hào)故障;故障3：在第51～100個(gè)樣本點(diǎn)上引入幅度為0.48×（50-t）2的拋物信號(hào)故障;故障4：在第51～100個(gè)樣本點(diǎn)上引入均值為12，方差為4的高斯噪聲。

為了驗(yàn)證本文所提方法的有效性，我們?cè)诿糠N方法下對(duì)各種故障分別重復(fù)3 000次實(shí)驗(yàn)。并記錄每次實(shí)驗(yàn)中T2統(tǒng)計(jì)量和SPE統(tǒng)計(jì)量的檢測(cè)正確率以及所用方法的迭代次數(shù)，最后再把少部分的異常數(shù)據(jù)剔除掉并求取剩余數(shù)據(jù)在兩種方法下的均值和方差，統(tǒng)計(jì)結(jié)果如表1所示。

由表1可以看出，EM方法中的T2統(tǒng)計(jì)量對(duì)故障1檢測(cè)正確率的均值為0.277 3，而在VBEM方法中的均值為0.284 6。兩種方法下的SPE統(tǒng)計(jì)量對(duì)故障1檢測(cè)正確率的均值幾乎相等。VBEM中的T2統(tǒng)計(jì)量和SPE統(tǒng)計(jì)量對(duì)故障2、3、4檢測(cè)正確率的均值均高于EM方法。另外，從最后一列的平均值中也可以看到VBEM中的T2統(tǒng)計(jì)量和SPE統(tǒng)計(jì)量對(duì)總的故障檢測(cè)正確率的均值也均高于EM方法。所用VBEM方法的平均迭代次數(shù)為4.061 8低于EM方法的平均迭代次數(shù)4.666 2?？偟膩?lái)說(shuō)，VBEM方法對(duì)本文所述四類故障的檢測(cè)效果要優(yōu)于EM方法。

如圖1～8所示為在某次實(shí)驗(yàn)中所得到的兩種方法對(duì)各故障的檢測(cè)結(jié)果。圖中水平虛線為正常數(shù)據(jù)下的統(tǒng)計(jì)量的控制限，實(shí)線曲線分別為數(shù)據(jù)的T2統(tǒng)計(jì)量和SPE統(tǒng)計(jì)量，當(dāng)統(tǒng)計(jì)量超過(guò)它們的控制限時(shí)，則表明檢測(cè)到系統(tǒng)有故障發(fā)生。

如表2所示，在EM方法中T2統(tǒng)計(jì)量的漏報(bào)率較為嚴(yán)重。在VBEM方法中圖2和圖8的T2統(tǒng)計(jì)量的漏報(bào)率也較為嚴(yán)重。在EM方法中SPE統(tǒng)計(jì)量的漏報(bào)率較小。在圖3中SPE統(tǒng)計(jì)量的漏報(bào)率僅為2%。在圖5中SPE統(tǒng)計(jì)量的漏報(bào)率僅為14%。在VBEM方法中SPE統(tǒng)計(jì)量的檢測(cè)正確率為100%，漏報(bào)率為0%。

4? 結(jié)? 論

本文針對(duì)ARFA模型下的動(dòng)態(tài)過(guò)程故障檢測(cè)問(wèn)題，構(gòu)建了相應(yīng)的T2和SPE兩種統(tǒng)計(jì)量，并提出了一種新的變分貝葉斯EM故障檢測(cè)方法。在對(duì)上文所述四類故障的正確檢測(cè)率上，本方法要優(yōu)于傳統(tǒng)的EM故障檢測(cè)方法。但本文未對(duì)高階動(dòng)態(tài)系統(tǒng)和實(shí)際的工業(yè)數(shù)據(jù)進(jìn)行仿真實(shí)驗(yàn)，因此這也是今后的一個(gè)研究方向。

參考文獻(xiàn)：

[1] 何瀟，郭亞琦，張召，等.動(dòng)態(tài)系統(tǒng)的主動(dòng)故障診斷技術(shù) [J].自動(dòng)化學(xué)報(bào)，2020，46（8）：1557-1570.

[2] 王世林.基于多元統(tǒng)計(jì)分析的工業(yè)過(guò)程故障檢測(cè)研究 [D].北京：華北電力大學(xué)（北京），2018.

[3] 劉春燕.基于PCA的過(guò)程工業(yè)故障檢測(cè) [D].綿陽(yáng)：西南科技大學(xué)，2017.

[4] KU W，STORER R H，GEORGAKIS C. Disturbance detection and isolation by dynamic principal component analysis [J].Chemometrics and Intelligent Laboratory Systems，1995，30（1）：179-196.

[5] LI W H，QIN S J. Consistent dynamic PCA based on errors-in-variables subspace identification [J].Journal of Process Control，2001，11（6）：661-678.

[6] KIM D，LEE I B. Process monitoring based on probabilistic PCA [J].Chemometrics and Intelligent Laboratory Systems，2003，67（2）：109-123.

[7] 楊沛武.概率PCA多元統(tǒng)計(jì)方法在過(guò)程監(jiān)控中的應(yīng)用研究 [D].無(wú)錫：江南大學(xué)，2008.

[8] LI G，QIN S J，ZHOU D H. A New Method of Dynamic Latent-Variable Modeling for Process Monitoring [J].IEEE Transactions on Industrial Electronics，2014，61（11）：6438-6445.

[9] WEN Q J，GE Z Q，SONG Z H. Data-based linear Gaussian state-space model for dynamic process monitoring [J].AlChE Journal，2012，58（12）：3763-3776.

[10] 周樂(lè).基于概率的工業(yè)過(guò)程數(shù)據(jù)建模與故障檢測(cè) [D].杭州：浙江大學(xué)，2015.

[11] QIN J S. Statistical Process Monitoring：Basics and Beyond [J].Journal of Chemometrics，2003，17（8-9）：480-502.

[12] BOX G E P. Some Theorems on Quadratic Forms Applied in the Study of Analysis of Variance Problems，I. Effect of Inequality of Variance in the One-Way Clas-sification [J].The Annals of Mathematical Statistics，1954，25（2）：290-302.

[13] JAAKKOLA T S. Tutorial on variational approximation methods [M]//OPPER M，SAAD D. Advanced Mean Field Methods：Theory and Practice.Massachusetts：MIT Press，2001：129-159.

[14] JORDAN M I，GHAHRAMANI Z，JAAKKOLA T S，et al. An Introduction to Variational Methods for Graphical Models [J].Machine Learning，1999，37（2）：183-233.

[15] BEAL M J. Variational algorithms for approximate Bayesian inference [D].London：University College London，2003.

作者簡(jiǎn)介：李吉?。?995—），男，漢族，河南鄧州人，碩士研究生在讀，研究方向：故障檢測(cè)、圖像處理;章智杰（1996—），男，漢族，江蘇連云港人，碩士研究生在讀，研究方向：圖像處

理;李其操（1997—），男，漢族，浙江諸暨人，碩士研究生在讀，研究方向：圖像處理;通訊作者：董自?。?973—），男，漢族，江蘇連云港人，教授，博士研究生，研究方向：檢測(cè)與控制、通信技術(shù)。