應(yīng)用統(tǒng)計線性回歸的系統(tǒng)誤差最大似然配準

李佳煒，江 晶，吳衛(wèi)華，鄭玉軍

(1.空軍預(yù)警學院 預(yù)警情報系，湖北 武漢 430019；2.空軍預(yù)警學院 空天預(yù)警系，湖北 武漢 430019；3.中國人民解放軍94710部隊，江蘇 無錫 214000)

協(xié)同多傳感器系統(tǒng)通過將空間分離的傳感器組網(wǎng)協(xié)同，可以提高對彈道導(dǎo)彈和衛(wèi)星等快速目標的探測精度、覆蓋范圍和截獲概率，相比單傳感器系統(tǒng)有著更大的潛力和優(yōu)勢[1-2]。然而，系統(tǒng)的融合性能除了受限于傳感器各自的探測能力外，還受到固有的系統(tǒng)誤差影響。在一定數(shù)量級的系統(tǒng)誤差下，傳統(tǒng)的集中式卡爾曼濾波器[3](Centralized Kalman Filter，CKF)的跟蹤性能將嚴重惡化。因此，為改善系統(tǒng)整體的融合性能，需要首先對系統(tǒng)誤差下的多傳感器進行配準和補償[4-8]。

圍繞多傳感器的系統(tǒng)誤差配準問題，國內(nèi)外學者展開了深入的研究。文獻[9-11]針對天基光電傳感器，利用恒星期望位置與實際量測的殘差，構(gòu)建關(guān)于緩變誤差的偏差濾波器和目標狀態(tài)的非線性估計子，從而實現(xiàn)校正傳感器量測的同時估計目標狀態(tài)。在文獻[11]的基礎(chǔ)上，文獻[12-13]分別將導(dǎo)致視線量測偏差的系統(tǒng)誤差描述為姿態(tài)失配角和指向偏移角，建立相應(yīng)偏差矢量的轉(zhuǎn)移模型和已知參考點的量測模型，并通過非線性濾波方法實現(xiàn)對系統(tǒng)誤差的在線配準與補償。上述研究屬于在線估計類算法，有良好的實時性，但是需要覆蓋范圍內(nèi)參考點的位置信息，而且需要增加用于估計目標狀態(tài)的濾波器。在離線批處理類算法方面，文獻[14]提出適用于兩部傳感器相對配準的線性最小二乘(Linear Least Squares，LLS)估計子，并給出了其統(tǒng)計性能邊界。針對傳感器總數(shù)大于2的情況，文獻[15-16]提出最大似然配準(Maximum Likehood Registration，MLR)算法，其通過似然函數(shù)最大化迭代實現(xiàn)系統(tǒng)誤差和目標狀態(tài)的聯(lián)合估計，且適用于多種傳感器觀測體制。文獻[17-18]針對機動平臺的姿態(tài)時變情況，構(gòu)造了能夠同時估計傳感器量測、平臺姿態(tài)角系統(tǒng)誤差的機動最大似然配準算法。文獻[19]推導(dǎo)了WGS-84坐標系下無源傳感器量測對目標狀態(tài)的雅可比矩陣，提出了距離缺失情況下多空基無源傳感器最大似然配準算法。文獻[20]引入衛(wèi)星軌道定向、姿態(tài)角測量與傳感器觀測過程中包含的測量誤差的先驗信息來修正目標狀態(tài)的誤差協(xié)方差，進一步提高了最大似然配準算法對多星載光學傳感器的配準性能。

然而，上述最大似然配準算法均采用泰勒展開的一階項來線性化非線性量測轉(zhuǎn)換，其對非線性函數(shù)的近似僅達到一階多項式精度，限制了算法的配準性能。統(tǒng)計線性回歸(Statistical Linear Regression，SLR)是一種基于點的方法，其利用多個回歸點來近似函數(shù)的非線性，能精確到高階多項式且避免了計算復(fù)雜的雅可比或海塞矩陣[21-22]。因此，筆者提出一種基于統(tǒng)計線性回歸的最大似然配準算法。該算法利用一組不完全相同的回歸點，通過統(tǒng)計線性回歸構(gòu)建目標狀態(tài)關(guān)于去偏量測的回歸方程，并得到投影后目標狀態(tài)的前二階統(tǒng)計特性；同時利用似然最大化迭代準確求解配準方程，實現(xiàn)對系統(tǒng)誤差和目標狀態(tài)的聯(lián)合估計。

1 非隨機系統(tǒng)誤差的最大似然估計

考慮一個協(xié)同多傳感器系統(tǒng)對共同區(qū)域內(nèi)的同一目標進行連續(xù)同步觀測，各傳感器在觀測維度上受到不同程度的系統(tǒng)誤差影響，因此系統(tǒng)的量測方程可建模為

zk=h(xk)+b+wk，

(1)

(2)

由于目標狀態(tài)xk未知，在b和xk均未知的情況下，僅能得到似然函數(shù)p(Z|X，b)，其中目標狀態(tài)集合X={xk∶k=1，2，…，K}。通過聯(lián)合最大化關(guān)于X和b的似然函數(shù)，并根據(jù)噪聲項wk在時域內(nèi)相互獨立，可得X和b的最大似然估計值為

(3)

式(3)也稱為系統(tǒng)的配準方程，其中p(zk|xk，b)為量測矢量zk的條件概率密度。

2 統(tǒng)計線性誤差傳播理論框架

(4)

由式(4)可知，目標狀態(tài)投影誤差是由高斯量測噪聲經(jīng)非線性量測轉(zhuǎn)換函數(shù)傳播后產(chǎn)生的。在系統(tǒng)誤差bn已知的條件下，目標狀態(tài)xk，n的概率密度可近似為一個均值為xk的高斯分布。利用統(tǒng)計線性回歸對非線性量測轉(zhuǎn)換進行線性化處理，可近似求得目標狀態(tài)的一階、二階統(tǒng)計特性。

(5)

(6)

其中，

(7)

(8)

(9)

非線性函數(shù)g的線性化形式可表示為

xk=gL(yk)=Ayk+c+e。

(10)

式(10)也稱為xk關(guān)于yk的回歸方程，其中殘差e的均值和協(xié)方差矩陣為

(11)

(12)

其中，

(13)

因此，經(jīng)線性化函數(shù)gL傳播后的前二階矩分別為

(14)

var(gL(yk))=AΣykAT+(Σxk-AΣykAT)=Σxk。

(15)

3 基于統(tǒng)計線性回歸的最大似然配準算法

(16)

(17)

(18)

3.1 目標狀態(tài)估計

使得條件概率密度p(zk|xk，b)在目標狀態(tài)空間中最大化的xk就是目標狀態(tài)的最大似然估計?？紤]到噪聲序列wk，n在傳感器間具有相互獨立性，則可將p(zk|xk，b)表示為緊湊形式：

(19)

聯(lián)合式(10)，并忽略線性化誤差ek，n的影響，則式(19)在目標狀態(tài)空間中的近似形式可表示為

p(zk，1，zk，2，…，zk，N|xk，b)≈

(20)

根據(jù)式(13)，目標狀態(tài)誤差協(xié)方差Σxk，n可近似為

(21)

其中，

(22)

(23)

3.2 系統(tǒng)誤差估計

(24)

(25)

其中，

(26)

(27)

其中，

(28)

(29)

(30)

其中，

(31)

(32)

代入式(27)，式(32)可化簡為式(20)的形式：

(33)

(34)

3.3 算法執(zhí)行流程

SLR-MLR算法主要由系統(tǒng)誤差估計和目標狀態(tài)估計兩個部分組成。算法的具體執(zhí)行流程概括如下。

輸入：歷史量測集合Z={zk∶k=1，…，K}。

誤差估計：

fork=1，…，Kdo

forn=1，…，Ndo

end for

end for

end while

狀態(tài)估計：

根據(jù)式(27)，計算所有傳感器配準后量測在目標狀態(tài)空間中的投影Xk；

4 仿真實驗與性能分析

4.1 場景設(shè)置與評價指標

(35)

(36)

假設(shè)系統(tǒng)在目標發(fā)射60 s后開始進行觀測，量測采樣間隔為1 s，各地面站的地理坐標分別為?1=(120°E，38°N)、?2=(112.5°E，50°N)和?3=(105°E，38°N)。整個系統(tǒng)的觀測性能參數(shù)如下：各雷達的系統(tǒng)誤差分別設(shè)置為b1=[2.5 km，3.5°，-1.5°]T、b2=[-1.8 km，-3°，1°]T和b3=[2 km，-2.5°，-1°]T；量測誤差標準差均設(shè)置為σρn=0.1 km，σθn=0.2°和σφn=0.25°。仿真場景示意如圖1所示。

(a) 雷達位置與目標軌跡

(b) 量測在狀態(tài)空間中的投影

將系統(tǒng)誤差估計的平均誤差(Mean Error，ME)和誤差標準差(Error Standard Deviation，ESD)作為性能評價指標，來衡量算法對系統(tǒng)誤差的配準性能。此外，采用文獻[15]定義的克拉美羅下界(Cramér-Rao Lower Bound，CRLB)為性能參考邊界。該邊界是一個保守值，對于數(shù)量不大的蒙特卡羅仿真，能夠較準確地評估算法理論上的最優(yōu)估計性能。注意：系統(tǒng)誤差估計的克拉美羅下界的計算需要利用目標狀態(tài)的真實值xk，其表達式為

(37)

4.2 實驗結(jié)果與性能分析

仿真實驗首先通過單次蒙特卡羅運行，給出了量測采樣點數(shù)為400時的SLR-MLR配準效果，系統(tǒng)誤差估計和目標狀態(tài)融合結(jié)果如圖2與圖3所示。從配準結(jié)果可以看出，算法在迭代次數(shù)遞增至3～5次后，能較好地收斂到系統(tǒng)誤差的預(yù)設(shè)值附近，并且在任意一個雷達不同觀測維度上的估計精度能夠達到95%以上。此外，在量測集合經(jīng)配準補償后，算法的狀態(tài)融合結(jié)果穩(wěn)定且接近目標的真實狀態(tài)。

圖2 系統(tǒng)誤差估計結(jié)果

圖3 目標狀態(tài)融合結(jié)果

然后，通過100次蒙特卡羅運行，圖4和圖5給出了迭代總數(shù)為5、采樣總數(shù)為200～400情況下SLR-MLR的配準性能，分別為系統(tǒng)誤差估計的平均誤差和誤差標準差。由圖可知，隨著采樣點數(shù)的遞增，算法在各個觀測維度上的平均誤差趨于平穩(wěn)且收斂于0，誤差標準差接近于克拉美羅下界的開平方；在采樣總數(shù)增加至300后，其穩(wěn)態(tài)平均誤差絕對值的算術(shù)平均分別為0.153 8 km、0.266 3 km、0.230 7 km，0.007 7°、0.002 9°、0.003 3°，以及0.015 2°、0.010 1°、0.007 2°。結(jié)果表明：SLR-MLR能有效地解決協(xié)同多傳感器系統(tǒng)的配準問題，并且對系統(tǒng)誤差的估計是漸進無偏有效估計。

基于上述仿真場景和相同的性能參數(shù)，將筆者所提SLR-MLR算法與文獻[15-20]采用的基于一階泰勒多項式近似的最大似然配準算法(Taylor-MLR)進行比較。圖6給出了迭代總數(shù)為5、采樣總數(shù)為250～300情況下兩種算法對雷達1的配準性能，待測試的SLR-MLR中參數(shù)k取值為0、10、20、30、40。圖中，SLR-MLR的ESD明顯低于Taylor-MLR，并且增大參數(shù)k可以明顯地改善其配準性能。這是因為Taylor-MLR只能精確到一階多項式，而利用不敏變換生成加權(quán)回歸點的SLR-MLR可以達到三階多項式精度，因此可以更精確地捕捉觀測過程的非線性特征以及不確定性。

圖4 系統(tǒng)誤差估計的平均誤差

圖5 系統(tǒng)誤差估計的誤差標準差

圖6 算法性能比較

為進一步分析參數(shù)k對算法性能的影響，圖7給出了采樣總數(shù)為300時不同k下SLR-MLR的配準性能。由圖可知，當60≤k≤100時，SLR-MLR的距離ESD和俯仰ESD達到最優(yōu)，而其方位ESD在k增大到80后開始收斂，因此建議k的取值區(qū)間為[80，100]。此外，通過仿真實驗可知，兩種算法單個配準回路的平均運行時間分別為4.52 s、3.12 s，因此SLR-MLR的計算復(fù)雜度約為Taylor-MLR的1.45倍，計算復(fù)雜度沒有明顯增加。

圖7 不同參數(shù)k下SLR-MLR的配準性能

5 總 結(jié)

筆者對協(xié)同多傳感器系統(tǒng)探測彈道目標情況下的系統(tǒng)誤差配準問題進行了研究，提出了一種基于統(tǒng)計線性回歸的最大似然配準算法。通過理論分析和仿真實驗，可以得到如下結(jié)論：(1) SLR-MLR算法能較好地實現(xiàn)對系統(tǒng)誤差和目標狀態(tài)的聯(lián)合估計，其中獲得的系統(tǒng)誤差估計是漸進無偏有效估計；(2) 利用不敏變換生成加權(quán)回歸點，SLR-MLR能達到三階多項式精度，并且適當增大縮放比例參數(shù)可以提升算法的配準性能；(3) SLR-MLR的配準精度高于Taylor-MLR，且計算復(fù)雜度沒有明顯增加，具有一定的實際應(yīng)用價值。

需要注意的是，所提算法僅適用于靜態(tài)系統(tǒng)誤差的配準場景。對于時變系統(tǒng)誤差下的傳感器配準問題，采用偽量測方法得到系統(tǒng)誤差的閉式表達，再應(yīng)用遞歸最小二乘[24]或者遞歸期望最大化[25]算法，可實現(xiàn)對系統(tǒng)誤差估計的動態(tài)更新。這也是下一步協(xié)同多傳感器配準的研究難點和重點。