含時滯測量值下間歇過程的雙維狀態估計

祁鵬程，趙忠蓋，劉飛

?

含時滯測量值下間歇過程的雙維狀態估計

祁鵬程，趙忠蓋，劉飛

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

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

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

引 言

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

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

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

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

1 問題描述

1.1 間歇過程的雙維動態特性

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

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

圖1 間歇過程雙維動態示意圖

1.2 間歇過程的時滯特性

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

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

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

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

2 間歇過程狀態估計的實現

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

2.1 模型描述

2.1.1 雙維狀態轉移模型

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

對于其他狀態，即＞1，＞1，與當前批次的先前狀態和先前批次相同采樣時間的狀態相關，模型 如下

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

2.2 狀態估計：僅有在線檢測值的時刻

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

2.2.1 單維估計

根據第1批次的狀態空間模型，結合貝葉斯濾波，估計步驟如下[14]：

預測

更新

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

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

考慮粒子退化的問題，需要對粒子進行重采樣。重采樣后，每個粒子的權重變為1/[18]。

2.2.2 雙維估計

由于其他批次的初始狀態與其他狀態的轉移方程不同，本節又分初始狀態和其他狀態進行考慮。根據狀態空間模型和PF算法可以得到初始狀態的預測步驟如下

其中，假設初始狀態與末態是相互獨立的，即

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

（2）當＞2,平滑方式如下

初始狀態的更新步驟如下

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

預測

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

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

那么可得

更新

2.3 狀態估計：兩種檢測值并存的時刻

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

則

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

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

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

（1）第1批次的狀態估計實現

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

（2）其他批次的狀態估計實現

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

3 算法驗證

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

3.1 數字仿真

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

（1）雙維狀態轉移模型

第1批次

初始狀態

除初始狀態之外的其他批次

（2）時滯測量模型

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

圖4 數字仿真無時滯值與含時滯值狀態的初始狀態誤差和均方根誤差

3.2 啤酒發酵過程

發酵是啤酒生產中耗時最長的工序，是一個典型的間歇生產過程，生產過程中兩個關鍵參數為殘糖濃度和菌體濃度。在發酵階段，殘糖濃度不可過高也不可過低，過高則發酵遲緩、啤酒發酵度低、生物穩定性差，過低則微生物耗氧能力下降，溶氧偏高。而菌體濃度高低直接影響啤酒的風味，濃度過高易致酵母衰退、自溶、酒質下降，濃度過低則發酵緩慢、發酵副產物增加。啤酒發酵過程關于殘糖濃度和菌體濃度的狀態轉移模型如下[22]

在實際生產過程中，啤酒發酵是按批量處理的，且在發酵過程中通常會離線采樣分析來得到殘糖濃度和菌體濃度更準確的測量值，對生產過程進行監控、控制和優化。但離線分析有時滯，因此如何融合先前批次和時滯測量值的信息是提高殘糖濃度和菌體濃度估計效果的關鍵。根據圖1和圖2建立啤酒發酵過程含時滯測量值的雙維動態空間模型。

（1）雙維狀態轉移模型

對于啤酒發酵的第1批次，沒有先前批次的信息，當作連續過程處理

對于初始狀態，在上一批次的發酵結束后，會根據啤酒成品的質量來調整本批次原料的糖度和酵母濃度

雖然當前批次的發酵會根據原料變化、過程飄移等問題調整發酵的相應操作，但總體上還是重復上一批次的整個過程。因此當前批次的狀態與當前批次的先前狀態和先前批次相同采樣時間的狀態相關

（2）時滯測量模型

生產過程通常配有在線分析儀，如液相色譜分析儀、YSI2700等，從而得到殘糖濃度和菌體濃度的測量值對生產過程進行監控，但誤差較大。因此為了提高啤酒成品的質量，大概每隔0.5 h對生產過程離線檢測一次，通過實驗室分析得到更準確的測量值。模型如下

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

圖5 啤酒發酵過程無時滯值與含時滯值狀態的初始狀態誤差和均方根誤差

4 結 論

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

References

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

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

[3] 姜蕾蕾. 啤酒發酵過程狀態監測及控制策略的研究 [D]. 哈爾濱: 哈爾濱理工大學, 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

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

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

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