噪聲相關系統中估計誤差方差奇異的濾波融合算法

楊鈺晗，李 智，牛頓標，宋恩彬*

（1.四川大學數學學院，成都 610065；2.電子信息控制重點實驗室，成都 610063）

R.E.Kalman于1960年提出Kalman濾波［1］，自此Kalman濾波便在目標跟蹤、遙感、導航定位、戰場監視、制造過程監控和復雜機械狀態維護等實際領域［2］有著廣泛的應用。然而單個傳感器的Kalman濾波難以滿足高精度和高容錯性等要求［3］，而多傳感器Kalman濾波融合接收信息更全面，能適應不同的環境變化，獲得更準確的估計值，因此在眾多領域得到廣泛的應用，并逐漸成為研究的熱點。

多傳感器系統的濾波融合方式主要有中心式融合與分布式融合兩種。分布式融合相比中心式融合，雖然信息量有一定的損失，但計算量小、通信負擔低、網絡結構靈活、可靠性高，是目前應用最為廣泛的信息處理方法。所以，使信息損失達到最小，甚至等價于中心式融合的分布式Kalman濾波融合算法受到了越來越多的關注。

文獻［4］-［6］提出了經典分布式Kalman濾波最優融合算法，并證明了該分布式算法融合結果與中心式Kalman濾波等價，使得分布式Kalman濾波融合既具有分布式系統的眾多優點，又具有與中心式Kalman濾波一樣優良的濾波性能。然而在實際應用中，由于連續系統離散化，系統或信號處于同一噪聲源污染的環境以及模型轉化等原因，各傳感器之間的觀測噪聲往往是相關的，由此而形成的噪聲相關系統得到廣泛研究。例如，文獻［7］針對這類系統，提出了與中心式等價的一類融合算法。文獻［8］通過對測量值和相關參數進行線性變換，使傳感器噪聲相關性解耦，實現與中心式等價的估計性能。此外，在有些隨機動態系統中，狀態向量可能會受到各種約束［9］，例如，車輛在固定道路中移動，鐘擺沿圓周進行擺動等。這些系統固有的約束信息限制了目標狀態的運動軌跡，同時導致融合過程中實際的過程噪聲方差是不可逆的，估計誤差方差陣也是奇異的。文獻［10］針對帶有不等式約束線性系統提出了一類無味Kalman濾波算法；文獻［2］和文獻［11］針對估計誤差方差是奇異的動態系統提出了一類融合算法，其估計性能與中心式的融合算法一致，解決了這類問題。

但在實際應用中，多傳感器線性估計系統復雜多樣，可能出現估計誤差方差是奇異的噪聲相關系統，而上文提及的融合算法都只針對解決其中一類問題，無法適用于兩種問題同時存在的情況。本文針對這類特殊系統，提出了一個基于最優加權矩陣的融合算法，該算法放寬了已有算法對系統的觀測噪聲和過程噪聲的限制，即允許：（a）觀測噪聲方差陣和過程噪聲是奇異的；（b）不同傳感器之間的觀測噪聲是相關的。最后通過模擬仿真，將本文所提出的融合算法與中心式融合算法和文獻［2］中對應算法進行對比分析，比較三種算法的估計誤差，數值結果驗證了所提融合算法的有效性。

1 問題描述

1.1 系統描述

本文考慮由l個傳感器構成的離散分布式動態系統，即

其中，F k∈R n×n為狀態轉移矩陣，H ik∈R li×n為第i個傳感器的狀態觀測矩陣，vk，xk∈R n分別為k時刻的狀態與過程噪聲，yik，wik∈R li分別為第i個傳感器k時刻的觀測向量和觀測噪聲。為方便后文描述，可將l個傳感器的觀測組合起來，得到聯合的觀測方程，即

其中，yk=(y1k'，…，ylk')'，H k=(H1k'，…，Hlk')'和wk=(w1k'，…，wlk')'。

1.2 過程噪聲和觀測噪聲統計特性描述

對于不同的時刻k，過程噪聲與觀測噪聲wik為不相關的零均值高斯隨機變量，但當i≠j時，因本文考慮不同傳感器的觀測噪聲wik與wjk在同一時刻k存在相關性，所以聯合的觀測噪聲方差陣為R k=cov(wk)，該矩陣第(i，j)塊為，可以為非零陣。其中主對角塊為第i個傳感器的觀測噪聲方差陣過程噪聲協方差矩陣為Q k=cov(vk)∈R n×n，由于本文放寬了對噪聲方差的限制，所以只需求Q k是一個半定陣。

1.3 經典Kalman濾波

根據文獻［12］中經典Kalman濾波的結果，利用所有傳感器觀測的中心式濾波迭代式，即

其中，(·)?表示矩陣的偽逆（Moore-Penrose逆），也稱廣義逆。為方便下文融合算法的描述與公式推導，類似中心式濾波的迭代，可得到第i個局部傳感器的濾波迭代式，即

1.4 問題難點

因本文考慮同時存在噪聲相關和估計誤差方差陣奇異的系統，所以其分布式融合存在兩個難點。第一，因為各傳感器的觀測噪聲具有相關性，所以觀測噪聲協方差矩陣R k不是對角塊的，無法利用文獻［2］的方法，通過對R k進行分塊直接得到分布式的融合算法。第二，因為過程噪聲協方差矩陣是奇異陣，估計誤差協方差矩陣也是奇異的，而文獻［9］所提融合算法需要用到上述兩個矩陣的逆矩陣，所以無法利用該算法融合各傳感器的估計。

因為上述兩個難點，現有文獻的方法都不再適用，所以下文考慮采用類似文獻［15］的最優加權算法進行融合。

2 分布式Kalman濾波融合算法

本文提出的分布式融合算法基于文獻［13］和［14］中最優加權融合矩陣的結論。假定x為系統要進行估計的n維狀態，是關于x的“局部”無偏估計，且下列的估計誤差協方差矩陣已知，即

所對應的誤差協方差矩陣為

其中，

根據上述公式（6）和（7）及文獻［11］的融合方法，為形成遞推式的融合公式，選擇2l+1個“局部”估計，即

根據公式（6），對于具有l個傳感器的動態系統（1）和（2），中心式Kalman濾波的最優線性融合估計由下式給出，即

其中2l+1個“局部”估計的估計誤差協方差陣為

類似于公式（7），上述融合估計的誤差協方差矩陣為

為遞推給出的誤差協方差矩陣中的元素，首先確定初值如下：

由公式（4）和式（5），可得到如下“局部”估計誤差的遞推式，即

其中，第二個等號由公式（4）代入得到，第三個等號由合并同類項得到。

類似的，有

其中，第二個等號由公式（5）代入得到，第三個等號由合并同類項得到。

其中，第二個等號由公式（5）代入得到，第三個等號通過合并同類項得到。

更進一步，由公式（13）和（15）可得Ck中各元素的迭代式，即

類似的，由公式（13）和（14）可得

類似的，由公式（14）和（15）可得

類似的，由公式（15）可得

類似的，由公式（15）可得

這樣就得到了Ck中的所有元素，從而獲得了Ck的遞推式，得到最優線性融合估計。

由公式（9）-（20）可以發現，此融合方法不涉及計算矩陣的逆，同時也不要求觀測噪聲協方差矩陣是對角塊的，因此可以很好地解決本文所考慮的分布式估計融合問題。

3 仿真實例

為驗證本文所提出的融合算法的性能，下面以實踐中廣泛應用的兩個目標運動模型作為例子進行模擬仿真，例1為目標在三維空間內做圓周運動，例2為目標在二維空間內做直線運動，且目標狀態都受到等式約束，從而導致融合過程中估計誤差協方差矩陣是奇異矩陣。在本例當中，假設各局部傳感器之間的觀測噪聲是相關的，因此系統的觀測噪聲方差矩陣是非對角塊的。

3.1 例1圓周運動

估計目標在三維空間內做勻速圓周運動，位置向量xk的隨機動態系統為

同時假設狀態xk受到等式約束，即被約束在矩陣D的零空間中，即

帶有約束公式（23）的動態系統公式（21）（22）就轉化為無約束系統，即

其中，P=I-D?D，(·)?表示矩陣的廣義逆，由于P是奇異陣，所以此時轉化后的“轉移矩陣”PFk也變為奇異陣。在此例中，假設

假設系統由兩個傳感器進行觀測，其觀測矩陣和帶有相關性的觀測噪聲協方差矩陣以下列形式給出，即

動態系統的初始值設定如下：

利用蒙特卡洛法進行500次模擬，每次模擬Kalman濾波迭代300步，利用平均目標位置向量的均方誤差損失Ek評估算法性能，即

其中，表示第i次模擬第k步迭代的融合后的估計值，表示第i次模擬第k步迭代的狀態值。

中心式濾波、本文所提出的最優加權分布式融合濾波和文獻［2］的融合算法的估計誤差數值結果如圖1所示。

圖1 圓周運動時不同Kalman濾波融合方法的估計誤差Fig.1 Estimation errors for different Kalman filter fusion methods of circle moving model

由圖1可見，文獻［2］算法估計誤差曲線在本文所提的最優加權分布式融合算法和中心式融合算法的誤差曲線上方，估計誤差較大；而本文所提的算法對應誤差曲線非常接近中心式融合誤差曲線，表明本文所提的算法具有較好的性能。

3.2 例2直線運動

假設目標在二維空間內做勻加速直線運動，狀態向量xk∈R4的前兩維表示估計目標的位置向量，后兩維表示估計目標的速度向量，目標依然服從等式約束，假設觀測、轉移矩陣如下：

觀測噪聲方差陣保持不變，等式約束矩陣和過程噪聲方差陣設為

初始值設定為

此時，中心式濾波、文獻［2］中融合濾波與本文所提的分布式融合濾波的數值結果如圖2所示。

由圖2可見，文獻［2］所提算法的估計誤差漸漸發散；而本文所提算法對應誤差曲線依然非常接近中心式融合誤差曲線，表現出穩健的性能，顯示了該算法的優越性。

圖2 直線運動時不同Kalman濾波融合方法的估計誤差Fig.2 Estimation errors for different Kalman filter fusion methods of straight line moving model

3.3 簡要分析

由上述數值結果可以看出，本文所提算法的融合結果誤差損失小，遠遠低于文獻［2］中算法的估計誤差，具有較高的精確度。仔細分析可以發現，文獻［2］所提算法之所以在例1和例2中性能較差，是因為文獻［2］中算法要求矩陣R k是對角塊的，再由如下等式算出，即［2］

其中，Kk(i)是公式（4）中矩陣K k的第i個列塊。利用R k可以分塊為，中心式濾波中帶觀測的項就可以由各傳感器的“局部”估計和一步預報線性給出，從而得到與中心式Kalman濾波等價的分布式融合算法。在本文所考慮的系統中，由于不滿足對角塊的條件，故文獻［2］中算法誤差較大。

而本文所提出的分布式融合算法性能雖然稍遜于中心式融合，但上述仿真結果表明兩者相差不大。此外，在實際應用中，多傳感器估計系統可能受到通信量的限制，或者觀測出現異常數據，此時分布式融合算法相比中心式算法具有更強的適應性和可靠性。

4 結 語

本文針對噪聲相關系統中狀態估計受狀態等式約束的情況，提出了一個基于線性最小均方誤差估計（MMSE）的分布式Kalman濾波融合算法，且融合算法的加權矩陣可以遞推得到。最后通過數值實驗與中心式濾波進行估計誤差對比，驗證了該算法具有較高的精確度。另外，實際應用中可能出現觀測值異常（比如丟包、延遲等），以及出現觀測噪聲相關且估計誤差協方差奇異的情況。針對該場景，考慮先剔除異常值再融合估計，以得到較好的估計性能，這也是接下來的一個研究方向。