含時滯測量值下間歇過程的雙維狀態(tài)估計

祁鵬程，趙忠蓋，劉飛

?

含時滯測量值下間歇過程的雙維狀態(tài)估計

祁鵬程，趙忠蓋，劉飛

（江南大學(xué)輕工過程先進控制教育部重點實驗室，江蘇 無錫 214122）

基于粒子濾波研究了間歇過程的狀態(tài)估計問題。根據(jù)間歇過程雙維動態(tài)特性,針對關(guān)鍵參數(shù)在線檢測精度低、離線分析時滯大等問題，分別建立一種雙維狀態(tài)轉(zhuǎn)移模型和時滯測量模型，并利用貝葉斯方法及前/后向平滑，提出一種含時滯測量值下的雙維狀態(tài)估計算法。該算法通過融合先前批次和時滯測量值的信息提高估計精度，并且克服了離線采樣周期和時滯時間不確定的問題。在數(shù)字仿真和啤酒發(fā)酵過程中的仿真應(yīng)用驗證了該算法的有效性。

間歇過程；雙維狀態(tài)空間模型；時滯；粒子濾波；貝葉斯方法；前/后向平滑

引 言

間歇過程也叫批次過程，廣泛應(yīng)用于化工、食品和醫(yī)藥等行業(yè)[1]。生產(chǎn)過程中，關(guān)鍵參數(shù)的變化，如發(fā)酵過程中殘?zhí)菨舛?、菌體濃度等，直接影響產(chǎn)品的質(zhì)量[2]。但是受限于傳感檢測技術(shù)，通常很難獲得穩(wěn)定準確的關(guān)鍵參數(shù)在線檢測值，而是通過估計的方法獲取關(guān)鍵參數(shù)的估計值。在貝葉斯框架下，給定檢測值，通過近似計算狀態(tài)空間的后驗概率密度函數(shù)可以解決間歇過程中關(guān)鍵參數(shù)的估計問題[3]。

目前，許多學(xué)者對間歇過程中的狀態(tài)估計問題進行了大量研究[4-10]?？紤]間歇精餾塔的非線性特性，以及在線成分分析成本高等問題，文獻[5]基于擴展卡爾曼濾波（EKF）設(shè)計最優(yōu)估計器從而得到系統(tǒng)的組成成分。然而，卡爾曼及其擴展形式只能處理高斯噪聲的濾波問題，基于此，文獻[7]中引入粒子濾波(PF)算法并應(yīng)用到間歇過程的狀態(tài)估計中，取得了良好效果。上述研究只考慮了間歇過程的非線性、非高斯特性，狀態(tài)僅隨時間轉(zhuǎn)化，忽略了批次之間的動態(tài)特性。實際上，典型的間歇過程包含許多批次，關(guān)鍵參數(shù)的動態(tài)變化不僅存在于時間維度，還存在于批次維度[8]。與連續(xù)過程不同，間歇過程中所有原料分批次處理，當(dāng)前批次會根據(jù)先前批次的信息來調(diào)整操作步驟，如產(chǎn)品的質(zhì)量。文獻[9]首次考慮了當(dāng)前批次的初始狀態(tài)與先前批次相關(guān)的情況，利用平滑的方式將先前批次的信息融合到當(dāng)前狀態(tài)的估計中，提高了估計精度。文獻[10]在上述文獻的基礎(chǔ)上考慮了所有狀態(tài)，當(dāng)前狀態(tài)既與當(dāng)前批次的先前狀態(tài)相關(guān)也與先前批次相同采樣時間的狀態(tài)相關(guān)，并建立雙維狀態(tài)空間模型對關(guān)鍵參數(shù)進行估計，取得了良好效果。

另一方面，離線檢測在間歇過程中普遍存在，離線檢測一般由人工在實驗室完成且檢測的樣本也是由人工采集。上述過程所用的時間相比于在線采樣時間間隔不可忽略，即離線檢測得到的測量值有很大的時滯。但在實驗室分析得到的數(shù)據(jù)比在線分析儀得到的更準確，因此融合時滯測量值的信息可以進一步提高關(guān)鍵參數(shù)的估計效果[11]。文獻[11]針對線性系統(tǒng)通過卡爾曼濾波和信息融合技術(shù)，實現(xiàn)兩種檢測值的信息融合。鑒于卡爾曼濾波簡單、計算快的優(yōu)點，文獻[12]基于EKF分別分析了全增廣和樣本增廣在處理時滯測量值的性能。文獻[13]根據(jù)時滯測量值與狀態(tài)為積分關(guān)系的特性，建立新模型并提出變維無跡卡爾曼濾波解決此模型的狀態(tài)估計問題。但是，卡爾曼濾波及其擴展形式在使用時噪聲需滿足高斯分布，有很大的局限性。

本文研究間歇過程中關(guān)鍵參數(shù)的估計問題,針對關(guān)鍵參數(shù)在線檢測精度低、離線分析時滯大等問題，結(jié)合間歇過程的雙維動態(tài)特性，建立新的狀態(tài)空間模型，基于PF提出一種含時滯測量值下的雙維狀態(tài)估計算法。該算法不僅考慮了先前批次的信息，還融合了離線檢測的時滯測量值。分別針對僅有在線檢測值和兩種檢測值并存兩種情況，利用前/后向平滑(FBS)對先前批次進行平滑從而融合先前批次的信息。處理時滯測量值時，運用貝葉斯方法更新粒子權(quán)重實現(xiàn)兩種檢測值的信息融合，通過調(diào)整平滑步數(shù)解決時滯時間不確定的情況，并進行了數(shù)字仿真和啤酒發(fā)酵過程的仿真應(yīng)用。

1 問題描述

1.1 間歇過程的雙維動態(tài)特性

一個典型的間歇過程包含許多批次。在單個批次中，與連續(xù)過程一樣，狀態(tài)隨時間演化，具有時間動態(tài)特性。而在批次間，由于原料變化、運行調(diào)整等原因，狀態(tài)在批次維度上同樣具有動態(tài)特性。考慮到批次之間開始時加入的原料有緩慢的變化，且一般會根據(jù)最后產(chǎn)品的質(zhì)量來調(diào)節(jié)下一批次的原料，因此考慮當(dāng)前批次的初始狀態(tài)與先前批次的初態(tài)和末態(tài)相關(guān)的情況。對于其他狀態(tài)，同樣會根據(jù)先前批次相同采樣時間的運行狀態(tài)調(diào)整操作，所以考慮其他狀態(tài)既與當(dāng)前批次的先前狀態(tài)相關(guān)也與先前批次相同采樣時間的狀態(tài)相關(guān)的情況[10]。

假設(shè)間歇過程包含個批次，每個批次有個采樣時刻，則雙維動態(tài)模型示意圖如圖1所示。其中，和代表第批次在采樣時刻的狀態(tài)和檢測值。水平和豎直的箭頭代表狀態(tài)的轉(zhuǎn)換方向，斜箭頭表示每個采樣時刻得到的檢測值。實線表示兩個狀態(tài)之間沒有其他的狀態(tài)，而虛線代表有其他的狀態(tài)。由于第一批次沒有先前批次，與連續(xù)過程一樣，只有時間上的動態(tài)特性。而其他批次不僅具有時間維度的動態(tài)關(guān)系，還有批次維度的動態(tài)關(guān)系。

圖1 間歇過程雙維動態(tài)示意圖

1.2 間歇過程的時滯特性

間歇過程中，在線檢測值直接由傳感器測得，即時但精度低，離線檢測值通過人工離線分析得到，有時滯但精度高，且兩種檢測值的采樣周期不同，給間歇過程中關(guān)鍵參數(shù)的估計帶來了困難，如何融合時滯測量值是提高估計效果的關(guān)鍵。間歇過程含時滯測量值的示意圖如圖2所示，其中，表示離線檢測時刻，表示第批次在時刻的在線檢測值，表示第批次在時刻的離線檢測值，表示離線檢測值的時滯時間。垂直的箭頭表示在線檢測值沒有時滯。虛線箭頭表示離線檢測值在時刻采樣，檢測值經(jīng)過時滯在時刻得到。因為離線采樣周期較大，通常最近一次離線采樣的檢測值在進行下一次離線采樣之前就已經(jīng)得到。

圖2 間歇過程含時滯測量值示意圖

另外，在實際生產(chǎn)過程中，在線采樣頻率固定且周期小，而離線采樣由人工完成，采樣頻率不固定且周期大，且由于離線分析的操作時間不是確定的，因此時滯時間也可能是時變的。

定義1：第批次時刻所有在線檢測值的集合為，離線檢測值的集合為。

2 間歇過程狀態(tài)估計的實現(xiàn)

首先根據(jù)間歇過程的雙維動態(tài)特性和時滯特性，分別建立新的雙維狀態(tài)轉(zhuǎn)移模型和時滯測量模型?？紤]兩種檢測值的采樣周期不同，且離線采樣周期不固定，本文從兩方面進行分析：僅有在線檢測值的時刻和在線檢測值、離線檢測值同時存在的時刻。

2.1 模型描述

2.1.1 雙維狀態(tài)轉(zhuǎn)移模型

第2及以后批次需考慮先前批次的信息。對于初始狀態(tài)，即＞1，，與上一批次的初態(tài)和末態(tài)相關(guān)，其模型如下

對于其他狀態(tài)，即＞1，＞1，與當(dāng)前批次的先前狀態(tài)和先前批次相同采樣時間的狀態(tài)相關(guān)，模型 如下

2.1.2 時滯測量模型 根據(jù)圖2可以得到間歇過程含時滯測量值的測量模型如下

2.2 狀態(tài)估計：僅有在線檢測值的時刻

由模型描述可知，第1批次狀態(tài)僅隨時間變化，屬于單維情況，而其他批次屬于雙維情況，所以分單維和雙維兩種情況分析。

2.2.1 單維估計

根據(jù)第1批次的狀態(tài)空間模型，結(jié)合貝葉斯濾波，估計步驟如下[14]：

預(yù)測

更新

對于非線性非高斯過程，PF利用Monte Carlo（MC）方法估計目標狀態(tài)后驗分布，取得了良好效果。其原理是用樣本均值代替積分運算[15]

其中，為粒子數(shù)。但是真實的狀態(tài)分布情況未知，通過重要性采樣分布可以獲取先驗粒子[16]，則

考慮粒子退化的問題，需要對粒子進行重采樣。重采樣后，每個粒子的權(quán)重變?yōu)?/[18]。

2.2.2 雙維估計

由于其他批次的初始狀態(tài)與其他狀態(tài)的轉(zhuǎn)移方程不同，本節(jié)又分初始狀態(tài)和其他狀態(tài)進行考慮。根據(jù)狀態(tài)空間模型和PF算法可以得到初始狀態(tài)的預(yù)測步驟如下

其中，假設(shè)初始狀態(tài)與末態(tài)是相互獨立的，即

因為隨著得到的測量值信息增多，初態(tài)對末態(tài)的影響較小，且每個批次的長度較長，可以忽略初態(tài)與末態(tài)的關(guān)系[10]。而在前一批次估計中已得到，只要計算可得，很明顯是一個平滑分布。本文采用FBS平滑算法，該方法保存每個采樣時刻的粒子作為平滑粒子，在當(dāng)前時刻估計狀態(tài)的基礎(chǔ)上遞推計算先前時刻的平滑狀態(tài)，且有效地解決了平滑過程中粒子退化的問題[19]。具體如下

（2）當(dāng)＞2,平滑方式如下

初始狀態(tài)的更新步驟如下

對于其他狀態(tài)，根據(jù)狀態(tài)空間模型和PF算法可以得到其他狀態(tài)的估計步驟如下：

預(yù)測

與式（11）的處理方式一樣，同樣需要假設(shè)

表示兩個狀態(tài)在時間維度上相差一個批次的長度就相互獨立。通過計算平滑分布即可得到。而在計算前必須先計算，具體如下

那么可得

更新

2.3 狀態(tài)估計：兩種檢測值并存的時刻

在處理時滯測量值前，首先根據(jù)在線檢測值的信息更新當(dāng)前時刻的粒子權(quán)重，然后再通過貝葉斯方法融合含時滯的離線檢測值，再次更新粒子權(quán)重。具體如下

則

對于單維估計，平滑的處理方式與式（15）相同，對于雙維估計與式（17）相同，只是在時，是和同時作為時刻檢測值時的權(quán)重，通過遞推計算即可得到分布對應(yīng)的粒子和權(quán)重。

實際操作中，由于離線檢測值是需要經(jīng)過人工實驗分析才能得到，操作時間具有不確定性，即時滯時間是不確定的。另外，離線采樣周期不固定，即離線檢測時刻不固定。從本節(jié)可以看出融合時滯測量值方法的關(guān)鍵一步是計算，在本文式是采用FBS算法進行平滑處理。因為當(dāng)前時刻（時滯值到達時刻）等于離線檢測時刻與時滯時間之和，在知道和下就可以確定，而和在實際過程中很容易得到，所以一旦確定和就可以根據(jù)當(dāng)前調(diào)整平滑步數(shù)計算出。因此本文提出的算法可以克服時滯時間和采樣周期不確定的問題。

綜上，含時滯測量值下間歇過程的雙維狀態(tài)估計示意圖如圖3所示，具體步驟如下。

（1）第1批次的狀態(tài)估計實現(xiàn)

⑤重采樣：根據(jù)歸一化后的權(quán)重得到個近似服從分布的粒子，并重置所有粒子權(quán)重為1/。

（2）其他批次的狀態(tài)估計實現(xiàn)

⑦重采樣：根據(jù)歸一化后的權(quán)重得到個近似服從分布的粒子，并重置所有粒子權(quán)重為1/。

3 算法驗證

為了比較間歇過程中考慮時滯測量值與不考慮時滯測量值的估計效果，選擇每個批次的初始狀態(tài)誤差（ISE）和均方根誤差（RMSE）[21]作為性能指標

3.1 數(shù)字仿真

根據(jù)文獻[8]得到間歇過程的數(shù)字仿真模型如下。

（1）雙維狀態(tài)轉(zhuǎn)移模型

第1批次

初始狀態(tài)

除初始狀態(tài)之外的其他批次

（2）時滯測量模型

圖4為數(shù)字仿真的效果圖。其中，黑色虛線是沒有考慮時滯測量值的誤差，紅色點線代表融合了時滯測量值的誤差。由圖可知，隨著批次的增加，狀態(tài)的ISE和RMSE都越來越小，趨于收斂，且在考慮時滯測量值的情況下狀態(tài)的ISE和RMSE比不考慮時滯測量值時的小，表明在考慮時滯測量值的情況下關(guān)鍵參數(shù)的估計效果較好。

圖4 數(shù)字仿真無時滯值與含時滯值狀態(tài)的初始狀態(tài)誤差和均方根誤差

3.2 啤酒發(fā)酵過程

發(fā)酵是啤酒生產(chǎn)中耗時最長的工序，是一個典型的間歇生產(chǎn)過程，生產(chǎn)過程中兩個關(guān)鍵參數(shù)為殘?zhí)菨舛群途w濃度。在發(fā)酵階段，殘?zhí)菨舛炔豢蛇^高也不可過低，過高則發(fā)酵遲緩、啤酒發(fā)酵度低、生物穩(wěn)定性差，過低則微生物耗氧能力下降，溶氧偏高。而菌體濃度高低直接影響啤酒的風(fēng)味，濃度過高易致酵母衰退、自溶、酒質(zhì)下降，濃度過低則發(fā)酵緩慢、發(fā)酵副產(chǎn)物增加。啤酒發(fā)酵過程關(guān)于殘?zhí)菨舛群途w濃度的狀態(tài)轉(zhuǎn)移模型如下[22]

在實際生產(chǎn)過程中，啤酒發(fā)酵是按批量處理的，且在發(fā)酵過程中通常會離線采樣分析來得到殘?zhí)菨舛群途w濃度更準確的測量值，對生產(chǎn)過程進行監(jiān)控、控制和優(yōu)化。但離線分析有時滯，因此如何融合先前批次和時滯測量值的信息是提高殘?zhí)菨舛群途w濃度估計效果的關(guān)鍵。根據(jù)圖1和圖2建立啤酒發(fā)酵過程含時滯測量值的雙維動態(tài)空間模型。

（1）雙維狀態(tài)轉(zhuǎn)移模型

對于啤酒發(fā)酵的第1批次，沒有先前批次的信息，當(dāng)作連續(xù)過程處理

對于初始狀態(tài)，在上一批次的發(fā)酵結(jié)束后，會根據(jù)啤酒成品的質(zhì)量來調(diào)整本批次原料的糖度和酵母濃度

雖然當(dāng)前批次的發(fā)酵會根據(jù)原料變化、過程飄移等問題調(diào)整發(fā)酵的相應(yīng)操作，但總體上還是重復(fù)上一批次的整個過程。因此當(dāng)前批次的狀態(tài)與當(dāng)前批次的先前狀態(tài)和先前批次相同采樣時間的狀態(tài)相關(guān)

（2）時滯測量模型

生產(chǎn)過程通常配有在線分析儀，如液相色譜分析儀、YSI2700等，從而得到殘?zhí)菨舛群途w濃度的測量值對生產(chǎn)過程進行監(jiān)控，但誤差較大。因此為了提高啤酒成品的質(zhì)量，大概每隔0.5 h對生產(chǎn)過程離線檢測一次，通過實驗室分析得到更準確的測量值。模型如下

圖5為啤酒發(fā)酵過程估計效果圖，分別為殘?zhí)菨舛群途w濃度的ISE和RMSE。黑色虛線是沒有考慮時滯測量值的誤差，紅色點線代表融合了時滯測量值的誤差。從圖中可以看出，隨著批次的增加，狀態(tài)的ISE和RMSE都越來越小，趨于收斂。由于每個批次只離線檢測兩次，即只能得到兩個時滯測量值的信息，導(dǎo)致估計效果沒有明顯提高，但在考慮時滯測量值的情況下狀態(tài)估計的總體效果要好于不考慮時滯測量值時的情況。

圖5 啤酒發(fā)酵過程無時滯值與含時滯值狀態(tài)的初始狀態(tài)誤差和均方根誤差

4 結(jié) 論

本文針對間歇過程，建立一種雙維狀態(tài)空間模型，研究了含時滯測量值下關(guān)鍵參數(shù)的估計問題。分別分析了僅有在線檢測值和兩種檢測值并存兩種情況?？紤]到間歇過程非線性、非高斯特性，引入PF算法進行狀態(tài)估計。為充分利用離線檢測值，以及間歇過程沿批次維度上歷史批次信息，通過貝葉斯方法對粒子權(quán)重進行更新，并采用FBS算法進行平滑，實現(xiàn)了多源信息的融合。數(shù)字仿真和啤酒發(fā)酵的應(yīng)用仿真結(jié)果說明：在間歇過程中融合時滯測量值可以有效地提高關(guān)鍵參數(shù)的估計效果。

References

[1] 侯迪波. 間歇生產(chǎn)過程中的知識發(fā)現(xiàn)方法研究 [D]. 杭州: 浙江大學(xué), 2004.
HOU D B. Research on knowledge discovery method in batch process [D]. Hangzhou: Zhejiang University, 2004.

[2] 李慧娟. 青霉素發(fā)酵過程建模與優(yōu)化 [D]. 沈陽: 東北大學(xué), 2010.
LI H J. Modeling and optimization of penicillin fermentation process [D]. Shenyang: Northeastern University, 2010.

[3] 姜蕾蕾. 啤酒發(fā)酵過程狀態(tài)監(jiān)測及控制策略的研究 [D]. 哈爾濱: 哈爾濱理工大學(xué), 2014.
JIANG L L. Study on state monitoring and control strategy of beer fermentation process [D]. Harbin: Harbin University of Science and Technology, 2014.

[4] KOZUB D J, MACGREGOR J F. State estimation for semibatch polymerization reactors [J]. Chemical Engineering Science, 1992, 47 (5): 1047-1062.

[5] OISIOVICI R M, CRUZ S L. State estimation of batch distillation columns using an extended Kalman filter [J]. Chemical Engineering Science, 2000, 55 (20): 4667-4680.

[6] WANG J L, ZHAO L Q, YU T. On-line estimation in fed-batch fermentation process using state space model and unscented Kalman filter [J]. Chinese Journal of Chemical Engineering, 2010, 18 (2): 258-264.

[7] CHEN T, JULIAN M, ELAINE M. Particle filters for state and parameter estimation in batch processes [J].Journal of Process Control, 2005, 15 (6): 665-673.

[8] YAO Y, GAO F. Batch process monitoring in score space of two-dimensional dynamic principal component analysis (PCA) [J]. Industrial & Engineering Chemistry Research, 2007, 46 (24): 8033-8043.

[9] ZHAO Z G, HUANG B, LIU F. Bayesian method for state estimation of batch process with missing data [J]. Computers and Chemical Engineering, 2013, 53 (11): 14-24.

[10] ZHAO Z G, HUANG B, LIU F. State estimation in batch process based on two-dimensional state-space model [J]. Industrial & Engineering Chemistry Research, 2014, 53 (50): 19573-19582.

[11] GUO Y F, ZHAO Y, HUANG B. Development of soft sensor by incorporating the delayed infrequent and irregular measurements [J]. Journal of Process Control, 2014, 24 (11):1733-1739.

[12] GOPALAKRISHNAN A, KAISARE N S, NARASIMHAN S. Incorporating delayed and infrequent measurements in extended Kalman filter based nonlinear state estimation [J]. Journal of Process Control, 2011, 21 (1): 119-129.

[13] GUO Y F, HUANG B. State estimation incorporating infrequent, delayed and integral measurements [J]. Automatica, 2015, 58: 32-38.

[14] FERDINAND E U. State-space estimation with a Bayesian filter in a coupled PDE system for transient gas flows [J]. Original Research Article, 2015, 39 (2): 682-692.

[15] ANDRIEU C, DOUCET A, SINGH S,. Particle methods for change detection, system identification, and control [J]. Proceedings of the IEEE, 2004, 92 (3): 423-443.

[16] DOUCET A, GODSILL S, ANDRIEU C. On sequential Monte Carlo sampling methods for Bayesian filtering [J]. Statistics and Computing, 2000, 10 (3): 197-208.

[17] STORVIK G. Particle filters for state-space models with the presence of unknown static parameters [J]. IEEE Trans. Signal Processing, 2002, 50 (2): 281-289.

[18] RIGATOS G G. particle filtering for state estimation in nonlinear industrial systems [J]. Instrumentation and Measurement, 2009, 58 (11): 3885-3900.

[19] NURMINEN H, RISTIMAKI A, ALI-LOYTTY S,. Particle filter and smoother for indoor localization[C]//4th International Conference on Indoor Positioning and Indoor Navigation. Montbéliard, France: IEEE, 2013: 1-10.

[20] ORGUNER U, GUSTAFSSON F. Storage efficient particle filters for the out of sequence measurement problem[C]//11th International Conference on Information Fusion. Piscataway, NJ, USA: IEEE, 2008: 1-8.

[21] KOTTAKKI K K, BHARTIYA S, BHUSHAN M. State estimation of nonlinear dynamical systems using nonlinear update based unscented Gaussian sum filter [J]. Journal of Process Control, 2014, 24 (9): 1425-1443.

[22] ZHAO Z G, HUANG B, LIU F. Parameter estimation in batch process using EM algorithm with particle filter [J]. Computers and Chemical Engineering, 2013, 57 (15): 156-172.

Two dimensional state estimation in batch process with delayed measurements

QI Pengcheng, ZHAO Zhonggai, LIU Fei

(Key Laboratory of Advanced Control for Light Industry Processes, Ministry of Education, Jiangnan University,
Wuxi 214122, Jiangsu, China)

This work investigates the state estimation in a batch process based on the particle filter method. Considering that the two dimensional dynamics and the key parameters are obtained online with low accuracy or offline with large time delay, a two dimensional state transition model and delayed measurement model are developed. In addition, a two dimensional state estimation algorithm is proposed by using Bayesian method and the forward-backward smoothing algorithm. The proposed algorithm improves the estimation accuracy by fusing information of previous batches and delayed measurements, and overcomes the influence of the uncertainty of offline sampling period and time delay. The applications in a numerical example and a beer fermentation show the effectiveness of the proposed method.

batch process; two dimensional state space model; time delay; particle filter; Bayesian method; forward-backward smoothing

supported by the National Natural Science Foundation of China (61134007, 61573169).

date: 2015-12-14.

Prof. ZHAO Zhonggai, gaizihao@ jiangnan.edu.cn

TQ 920.6

A

0438—1157（2016）09—3784—09

10.11949/j.issn.0438-1157.20151894

國家自然科學(xué)基金項目（61134007，61573169）。

2015-12-14收到初稿，2016-04-28收到修改稿。

聯(lián)系人：趙忠蓋。第一作者：祁鵬程（1991—），男，碩士研究生。