基于模糊有序聚類算法的間歇過程故障檢測

張 成，趙海濤，孫韶媛

(1. 華東理工大學 信息科學與工程學院，上海 200030；2. 東華大學 信息科學與技術學院，上海 201620)

由于間歇生產過程具有高附加值、小批量和品種多等優點，現已成為現代工業中的一種重要的生產方式[1]。因此，對間歇過程進行有效的故障檢測，及時發現故障變得十分重要。常用的方法有多向主成分分析(multi-way principal component analysis， MPCA)和多向偏最小二乘 (multi-way partial least squares， MPLS)法[2-4]。間歇過程的一個重要特性為多階段性，每一個生產過程都由幾個子階段組成。但是，傳統MPCA和MPLS方法沒有考慮到間歇過程的多階段特性，只是對所有的批次建立單一的檢測模型，這樣會造成較高的故障誤報率和漏報率。

針對以上問題，研究學者進行了大量的研究。Lu等[5]提出了一種基于階段劃分的故障檢測策略，利用K均值聚類(K-means)算法對過程數據進行階段劃分，對不同的子階段分別建立檢測模型，取得了較好的檢測效果。該聚類算法利用樣本與初始類別之間的相似性，確定樣本所屬的類別，并根據所有類別的平方誤差和來確定最終的聚類結果。Liu等[6]將近鄰傳播聚類算法(affinity propagation clustering algorithm， APCA)應用于間歇過程的故障檢測。APCA不需要事先指定聚類中心以及聚類個數，而且對數據的初始值不敏感。Yew等[7]將模糊C均值聚類算法(fuzzyC-means，FCM)應用于間歇過程的子階段劃分。FCM由K-means算法改進而來，引入了一個隸屬度函數，通過各個樣本與每個類別的隸屬度來決定該樣本所屬的類別[8]。除此之外，研究學者也對階段劃分做了大量的工作[9-13]。在實際生產中，過程數據是按時間順序來進行采集的，每個子階段內的數據在時間上必須是連續的。然而，以上幾種算法都沒有考慮到過程數據的連續性，導致子階段劃分出現斷續問題，即子階段內的數據不是某一時間段內監測到的連續數據，出現了不連續的問題。

圖1 基于時間片的聚類過程Fig.1 Clustering process based time slice

此外，很多方法利用三維數據矩陣相鄰時間片數據矩陣之間的相似性或者其負載矩陣的相似性來進行聚類[14-16]，其過程如圖1所示。其中:I為批次總數;K為每個批次的采樣總數;J為過程監測變量總數;Xi為第i個時間片的數據矩陣，即第i個采樣點處所有批次的監測數據。使用該方法就必須保證每個矩陣維度相同，這就要求每一個批次過程的生產時間要一致，由于對于所有的批次來說，J都是相同的，要使每個矩陣維度相同，必須要求所有批次的K相同。但在實際生產中，受于初始值和環境等因素的影響，每一個批次的生產時長都會有所不同，即出現數據不等長問題。

為了更準確地對間歇過程進行階段劃分，本文基于有序聚類算法提出了一種模糊有序聚類算法(fuzzy order clustering algorithm， FOCA)。在初始化聚類時，FOCA規定每個類別中的樣本按照采集時的先后順序保持連續，直到聚類完成。為了降低噪聲、離群點以及數據波動等對聚類結果的影響，FOCA借鑒了FCM中的隸屬度函數[8]，利用模糊的方法來降低數據波動的影響。此外，FOCA應用于階段劃分時，并未采用圖1中的方式，而是對每一個批次過程單獨進行子階段劃分，然后再將每個子階段的數據按變量方向展開。因此，在對批次過程進行子階段劃分時不會受數據不等長問題的影響。階段劃分完成之后，建立每個子階段的主成分分析(PCA)檢測模型，利用T2和平方預測誤差(squared prediction error, SPE)統計量對新采集樣本進行故障檢測。

1 基于模糊有序聚類的多階段建模

1.1 模糊有序聚類算法

為了保證對批次過程進行子階段劃分時不出現斷續的問題，在初始化類別時不能隨機選取數據作為類別樣本，可以對原批次過程數據進行截取，截取之后的類別數據屬于批次過程某一時間段內的數據。該初始化之后的聚類一定不存在斷續問題，每個聚類中的數據都是連續的。

設有一組樣本數據X=(x1，x2， …，xK)，包含了K個連續的樣本數據，xk∈RJ，k=1， 2， …，K。FOCA在初始化類別時只在樣本X中截取連續的一段樣本數據作為初始化類別數據。

(1)

式中：xk為類別Gc中的樣本數據；xic和xic+1-1分別為類別Gc的起始樣本和末尾樣本；ic和ic+1為相鄰類別之間的分割點序號；Nc為類別Gc中的樣本個數。

(2)

式中：μks為樣本xk與類別s的隸屬度；b為控制聚類結果模糊程度的常數，且b>1。

本文隸屬度的計算借鑒FCM中的隸屬函數，要求一個樣本與所有類別的隸屬度之和為1，即

(3)

隸屬度的計算公式如式(4)所示。

(4)

式中：S為需要計算隸屬度的類別總數。若類別c不是第一個或者最后一個類別，那么S為3，即將類別c之前的所有類別視為一類，類別c之后的所有類別視為一類，類別c單獨視為一類，總共3類；若類別c是第一個或者最后一個類別，則S為2。

類別初始化完成之后，需要根據準則函數來確定最終的聚類結果。設一聚類結果如式(5)所示。

P(K，C):{xi1，xi1+1， …，xi2-1}， {xi2，xi2+1， …，xi3-1}， …， {xiC，xiC+1， …，xK}

(5)

式中：P(K，C)表示將K個樣本聚成C類，其分割點為xi1，xi2， …，xiC，其中i1=1。每個類別中的樣本保持連續和有序。由式(2)可知，類別的誤差和越小，該類別的數據劃分越合理。當所有類別的誤差和達到最小值時，得到的便是最優聚類結果。因此，定義該聚類結果的準則函數如式(6)所示。

(6)

當K和C固定時，e值越小，表示所有類別的誤差和越小，聚類結果也就越合理。e達到最小值的P(K，C)將會是最優的聚類結果。則

e[P(K， 2)]=min{D(G1)+D(G2)}

(7)

D(Gc)}

(8)

式(7)表示將K個樣本聚成兩類。式(8)表示將K個樣本聚成C類。后者需先找到一個類別，即GC= (xiC，xiC+1， …，xK)，然后將剩下的iC-1個樣本聚成C-1類，依次類推找到所有的類別，所得到的聚類結果便是最優聚類。FOCA的具體步驟如下:

(1) 從第一個樣本開始，迭代截取連續的一段數據作為初始化聚類，此步驟需遍歷出所有的聚類結果；

(2) 根據式(2)計算每一個聚類的誤差和，用于表示該聚類的合理性，其中模糊策略的加入可優化聚類結果；

(3) 根據式(6)和計算出的聚類誤差和，計算出所有聚類結果的準則函數值；

(4) 根據式(7)和(8)以及準則函數值，計算出樣本聚類分割點，得到最優聚類結果。

1.2 算法分析

FOCA在確定初始化類別時，需遍歷所有可能的聚類結果。以第一個樣本為類別起始點，將第1到第2之間的樣本視為一類，將第1到第3之間的樣本視為一類，依次類推，將第1到第K之間的所有樣本視為一類；然后以第二個樣本為類別的起始點，同上推出所有的類別，直到遍歷出所有可能的聚類結果。由于每一個類別中的樣本前后有序，而且在時間上是連續的，這樣就保證了數據的連續性。

間歇過程的每一批次都是由幾個子階段組成，由一個子階段到相鄰的后一個子階段之間的過渡不可能瞬間完成，這中間有一個過渡過程[17]。過渡階段的數據不穩定，波動比較大，對聚類會產生一定的影響。過渡階段的數據一般位于子階段數據集的兩端，這一區域的樣本數據沒有確定的子階段屬性，因此可加上隸屬度進行軟聚類。由式(2)可知，隸屬度的加入可改變該類別的誤差和的大小，最終會通過式(6)改變聚類結果。但是，對于FOCA來說，每一個樣本點都有可能成為相鄰子階段之間的分割點，隸屬度的加入主要會對過渡區域中樣本的類別歸屬產生影響，并對該區域的聚類結果進行優化。而對于子階段中間部分即穩定區域樣本數據的類別歸屬影響不大。

當使用FOCA對批次過程數據進行子階段劃分時，對每一個批次數據單獨進行階段劃分，聚類完成之后，再將每一個子階段內的數據按變量方向進行展開。整個子階段的劃分過程如圖2所示。由圖2可知，FOCA不會受數據不等長問題的影響。

圖2 FOCA劃分子階段過程Fig.2 Phase partition process of FOCA

1.3 基于多階段的故障檢測模型

每個子階段的數據按變量方向展開成二維矩陣的形式：X(KI×J)，然后建立每個子階段的PCA模型。設有一子階段c的數據矩陣為Xc，其PCA模型如式(9)所示。

(9)

式中：Pc(J×R)為Xc的負載矩陣；R為主元個數；Tc為得分矩陣；Ec為殘差矩陣。

PCA模型建立之后，通常使用T2以及SPE兩個統計量來檢測過程的運行是否正常。對于一個待檢測樣本xnew來說，其T2統計量計算如式(10)所示。

(10)

式中：Pc，Sc分別表示xnew所屬子階段數據矩陣對應的負載矩陣以及特征值矩陣。

對于每個子階段來說，T2統計量近似服從F分布。所以每個子階段的T2統計量控制限可以通過式(11)來進行估算。

(11)

式中：Rc為主元個數；Nc為子階段c的樣本個數；FRc，Nc-Rc，α為自由度為Rc，Nc-Rc的F分布；α為置信度。

SPE統計量計算如式(12)所示。

(12)

式中：Qnew為SPE統計量；I為單位矩陣，Pc為xnew所屬子階段數據矩陣對應的負載矩陣。

每個子階段的SPE統計量近似服從χ2分布，每個子階段的SPE統計量控制限可以通過式(13)進行計算。

(13)

式中：ec和fc分別為子階段c中所有樣本數據SPE統計量值的均值和方差。

由以上分析可知，要計算一個待檢測樣本的T2和SPE統計量的值需要某一階段c的Pc、Rc、Nc等參數，該階段c是當前待檢測樣本所屬的子階段，即要對待檢測樣本進行檢測就需要知道該樣本屬于哪一個子階段。Qin等[18]設計了一個指標φ用于表示當前待檢測樣本與所有子階段的相似性，其計算式如式(14)所示。

(14)

1.4 基于FOCA的間歇過程故障檢測流程

基于FOCA的間歇過程故障檢測可以分為離線建模以及在線檢測兩個部分。

離線建模的步驟如下：

步驟1 采集正常生產過程數據作為訓練樣本，組成數據集X(I×J×K)。其中：I為總批次個數；J為過程監控的變量個數；K為批次運行過程中采集到的數據個數，每個批次會所有不同。

步驟2 對每個批次數據Xi(J×K) (i=1，2， …，I)進行標準化處理，然后通過FOCA對每個批次的數據分別進行子階段劃分，得到每個批次的子階段。

步驟3 將所有批次中相同的子階段進行整合，并按變量方向進行展開，整合之后的新數據代表了每個子階段的訓練樣本數據。

步驟4 對于每個子階段的訓練樣本數據，通過PCA算法確定每個子階段的檢測模型，并確定每個子階段的SPE和T2控制限，該控制限會被用于后續新樣本的檢測。

在線檢測的步驟如下：

步驟1 采集新的樣本數據xnew(J×1)，對其進行標準化處理，作為待檢測樣本數據。

步驟2 通過式(14)確定當前待檢測樣本屬于哪一個子階段。

步驟3 通過待檢測樣本所屬的子階段檢測模型，來計算待檢測樣本的SPE和T2統計量值。

步驟4 判斷SPE和T2統計量值是否超過控制限。如果超過，則表示過程出現故障，發出警報；如果未超過，則表示過程運行正常。

2 仿真試驗

青霉素的生產過程是一個典型的多階段、非線性且動態的間歇過程。2002年，由Birol等[19]設計出了青霉素的仿真模型Pensim 2.0，這一模型被廣泛應用于間歇過程的故障檢測研究。

本文試驗仿真采用Pensim 2.0生成仿真數據，包括50組正常批次和3組故障批次。每個批次的采樣周期為0.5 h，共產生800個左右的采樣值。共選取10個生產變量進行監控，包括通風率、攪拌功率、底料流速、pH值、物料溫度、二氧化碳濃度、氧氣濃度、產熱量、冷水流速和青霉素濃度。在進行仿真試驗時，將K-means方法作為比較對象。

2.1 子階段劃分

K-means和FOCA對正常批次29的子階段劃分結果如圖3所示。由從圖3(a)可知，階段劃分結果出現了很嚴重的斷續問題，除第1個階段外，第2和3階段的數據完全斷續。由圖3(b)可知，3個子階段內的數據連續，沒有出現斷續的問題。由此可以表明，FOCA對批次階段劃分的效果要優于K-means。

在FOCA對批次29聚類過程中，準則函數值的變化曲線如圖4所示。由1.2節分析可知，要將批次聚成3類需要兩個分割點，首先需要求第2個分割點然后再求第1個分割點。圖4便為求第2個分割點時準則函數值的變化曲線，求第1個分割點時的變化曲線與此類似。由圖4可知，準則函數值先急劇減小后緩慢增大，其最小值對應的樣本點便為第2個分割點。

(a) K-means

(b) FOCA圖3 K-means和FOCA對于批次29的階段劃分結果Fig.3 Phase partition result of K-means andFOCA for batch 29

圖4 FOCA迭代過程中準則函數值的變化Fig.4 Variation of criterion function value during the iterative process of FOCA

本次試驗采用的計算機配置：CPU為Intel Core i3-4500，內存為8G，系統為Win7 64位操作系統。經計算可知，FOCA對批次進行階段劃分的平均耗時為7 min 23 s。由于這是對訓練樣本數據的處理耗時，所以不會對在線故障檢測造成影響。

2.2 正常批次的檢測

為了驗證本文方法的有效性，選取一個正常樣本數據進行檢測，并與基于K-means的故障檢測方法進行對比，檢測結果如圖5所示。由圖5(a)可知，基于K-means檢測方法的T2和SPE統計量都誤檢測出了故障，且后者的誤報率較高，在大約90個樣本點處出現誤報。由圖5(b)可知，在過程正常運行期間，基于FOCA的故障檢測方法的T2和SPE統計量都處于控制限以下，沒有出現誤報。由此可知，基于FOCA的故障檢測方法對正常樣本的檢測效果要優于基于K-means的故障檢測方法。

(a) K-means

(b) FOCA圖5 基于K-means和FOCA的檢測方法對正常批次的檢測結果Fig.5 Monitoring result of normal batch processbased on K-means and FOCA

2.3 故障批次的檢測

通過Pensim 2.0仿真平臺可以生成不同變量(通風率、攪拌功率和培養基進料流速)、不同類型(階躍或斜坡)、不同時長、不同大小的故障類型。現對3組故障的檢測結果進行分析。

故障1為通風率的階躍故障，通風率以8%的階躍增加，故障時間為150～350 h。基于K-means和FOCA的故障檢測方法的檢測結果如圖6所示。由圖6可知，兩種方法都及時的檢測出了故障，而且沒有出現誤報。圖中橢圓包圍起來的部分為漏報部分。通過兩圖對比可知，基于FOCA的故障檢測方法的漏報率要比基于K-means的故障檢測方法的漏報率低很多。對于故障1的檢測，基于FOCA的故障檢測方法要優于基于K-means的方法。

(a) K-means

(b) FOCA圖6 基于K-means和FOCA的故障檢測方法對故障1的檢測結果Fig.6 Monitoring result of fault 1 based onK-means and FOCA

故障2為通風率階躍故障，通風率以10%的階躍增加，故障時間為85～335 h。基于K-means和FOCA的故障檢測方法的檢測結果如圖7所示。由7(a)可知，基于K-means故障檢測方法的T2和SPE統計量都及時檢測到了故障，但在故障發生期間，在大量采樣點處T2位于控制限以下，造成了很高的漏報率。而由7(b)可知，基于FOCA的故障檢測方法不僅能夠及時地檢測出故障，而且T2統計量沒有出現漏報，SPE統計量也僅在少量采樣點處位于控制限以下。綜合兩種統計量的檢測結果，基于FOCA的方法效果較好。

(a) K-means

(b) FOCA圖7 基于K-means和FOCA的故障檢測方法對故障2的檢測結果Fig.7 Monitoring result of fault 2 based onK-means and FOCA

故障3為攪拌功率的斜坡故障，攪拌功率以0.05 W/h增加，故障時間為100～300 h。基于兩種方法的檢測結果如圖8所示。由圖8可知，兩種方法對故障3進行檢測時都沒有出現誤報。但是，對比圖8(a)和8(b)可知，基于FOCA方法的T2統計量在第141.5 h處檢測到了故障，SPE統計量在第166 h處檢測到了故障；基于K-means方法的T2統計量在第147 h處檢測到了故障，SPE統計量在第232 h處檢測到了故障。由此可見，基于FOCA故障檢測法的第一次檢測出故障的時間要早于基于K-means的故障檢測方法。

經統計，基于FOCA的故障檢測方法對每個樣本檢測所耗費的時間為0.095 ms，在故障發生時，可以準確迅速地檢測出故障并發出警報，有助于現場操作人員及時做出相應的處理。

由以上3組故障檢測結果對比可知，在對不同故障類型、故障大小和故障時長的檢測上，基于FOCA的故障檢測方法要優于基于K-means的方法。這主要是由于FOCA對階段的劃分更加準確，基于階段劃分所建立的故障檢測模型更加符合生產過程的實際情況。因此，基于FOCA的故障檢測方法能夠得到更好的故障檢測結果。

(a) K-means

(b) FOCA圖8 基于K-means和FOCA的故障檢測方法對故障3的檢測結果Fig.8 Monitoring result of fault 3 based onK-means and FOCA

3 結 語

針對傳統聚類算法在對過程批次數據進行聚類時出現的斷續問題，本文提出的模糊有序聚類算法(FOCA)對過程數據的連續性進行了充分考慮，并加入了隸屬度，有效避免了斷續問題的發生，同時降低了噪聲對聚類的影響，優化了聚類結果。而且，FOCA只針對單一的批次進行子階段劃分，這樣既可減少對其他批次數據的影響，又能解決數據不等長問題。對青霉素生產過程的仿真試驗表明，基于FOCA的故障檢測方法比基于K-means的故障檢測方法有著更優的檢測效果。