EM算法在衛(wèi)星軌道計算中的應用＊

杜 琨，宋迎春，肖琴琴

（中南大學地球科學與信息物理學院，湖南 長沙410083）

0 引 言

利用GPS信號進行導航定位時，需要已知GPS衛(wèi)星在空間的瞬時位置，衛(wèi)星星歷的精度將直接影響定位的精度和結果，因此，精確計算GPS衛(wèi)星在空間的瞬時位置是精密定位的前提和基礎［1－2]。隨著現(xiàn)代測量、導航、氣象等領域對衛(wèi)星星歷的精度和實時性的要求越來越高，對利用GPS廣播星歷和精密星歷來準確計算衛(wèi)星在空間的瞬時位置，變得越來越重要。目前，已有很多文獻研究了利用精密星歷文件計算任意時刻的衛(wèi)星空間位置，且取得了很好的效果，然而這些研究工作都是基于數(shù)據(jù)完全，即星歷文件數(shù)據(jù)準確、信息量較大無缺失數(shù)據(jù)時進行研究的［3－4]。在實際應用中，精密星歷文件并不完全準確，由于衛(wèi)星遮擋、衛(wèi)星信號傳播限制等多方面的原因，可能造成部分有用信息丟失，導致數(shù)據(jù)不完全。在這種情況下，如果仍采用常規(guī)的計算方法，則可能受到數(shù)據(jù)不完全的影響。對于信息量不足，通常只能降階對其擬合或插值，這兩種情況都會造成很大的誤差，影響擬合或插值的效果。已有研究表明，對于15min歷元間隔的精密星歷而言，如果要精確至10－8，用8階朗格朗日內(nèi)插已足夠保證精度，而切比雪夫多項式擬合至少需要12階才能達到厘米級精度［3]。利用切比雪夫多項式擬合要達到厘米級精度，至少需要12組數(shù)據(jù)才能對其進行擬合。但由于精密星歷是每15min一組，數(shù)據(jù)量相對較少，所以對于短弧段用切比雪夫多項式擬合不能滿足要求。EM算法是一種常用的缺失數(shù)據(jù)處理方法，目前，已被廣泛應用于測量領域。2010年，惠沈盈等把EM算法應用于數(shù)據(jù)缺失時的動態(tài)定位解算［5]，2011年，林東方等又把EM算法應用于數(shù)據(jù)缺失時的GPS高程擬合［6]，都取得了很好的效果。通過對EM算法的研究，對不完全數(shù)據(jù)添加與衛(wèi)星軌道信息相關的“潛在數(shù)據(jù)”，能大大提高衛(wèi)星擬合的精度。

1 常用的插值法和擬合法

精密星歷文件是以離散的形式、按一定的時間間隔（通常為15min）給出衛(wèi)星在空間的三維坐標、三維改正速度以及衛(wèi)星鐘改正數(shù)等信息。而在GPS事后數(shù)據(jù)處理中，經(jīng)常要用到任意時刻的衛(wèi)星坐標，因此需要對精密星歷進行插值或者擬合以獲得采樣間隔更高的衛(wèi)星軌道信息。利用精密星歷文件計算衛(wèi)星位置以及運動速度通常采用拉格朗日插值和切比雪夫擬合。

關于拉格朗日插值和切比雪夫擬合的原理已經(jīng)在大量文獻中有過說明，僅以切比雪夫多項式擬合為例。用n階切比雪夫多項式來逼近時間段［t0，t0＋Δt]中的衛(wèi)星星歷時，先要將變量t∈［t0，t0＋Δt]變換為變量τ∈［－1，1]：

式中：n為多項式的階數(shù)；Cxi，Cyi，Czi為切比雪夫多項式Ti的系數(shù)。切比雪夫多項式Ti的遞推公式如下：

根據(jù)已知的衛(wèi)星坐標，用最小二乘法擬合出多項式系數(shù)Cxi，Cyi，Czi后，就可以計算該時段任意時刻的衛(wèi)星位置。

2 EM算法及其在衛(wèi)星軌道計算中的應用

IGS精密星歷文件經(jīng)常因為計算失誤、文件傳播等原因導致部分數(shù)據(jù)失真或缺失。在數(shù)據(jù)缺失或數(shù)據(jù)量較少的情況下，無法采用高階的插值法或擬合法。這時，只能降階對其插值或擬合，從而導致計算精度降低。采用EM算法，添加與衛(wèi)星軌道信息相關的“潛在數(shù)據(jù)”，可以有效解決這一問題，大大提高衛(wèi)星擬合的精度。

2.1 EM算法的原理

EM算法是一種迭代算法。它的每一次迭代由兩步組成：E步（求期望）和 M步（極大化）。一般，以P（θ／Y）表示θ的基于觀測數(shù)據(jù)的后驗分布密度，稱為觀測數(shù)據(jù)后驗分布。以P（θ／Y，Z）表示添加數(shù)據(jù)（缺失數(shù)據(jù)）Z后得到的關于θ的后驗分布密度函數(shù)，稱為完全數(shù)據(jù)后驗分布。P（Z／θ，Y）表示在給定θ和觀測數(shù)據(jù)Y下潛在數(shù)據(jù)（缺失數(shù)據(jù)）Z的條件分布密度函數(shù)。我們的目的是計算觀測后驗分布P（θ／Y）的參數(shù)，EM 算法如下進行，記θi為第i＋1次迭代開始時后驗參數(shù)的估計值，則第i＋1次迭代的兩步為：

E步：將P（θ／Y，Z）或logP（θ／Y，Z）關于Z 的條件分布求期望，從而把Z積掉，即

M 步：將Q（θ／θi，Y）極大化，即找一個點θi＋1，使

如此形成了一次迭代θi→θi＋1，θi＋1∈M（θi），這里M（θi）是在整個參數(shù)空間Ω 內(nèi)使得Q（θ／θi，Y）最大的θ的值所組成的集合。將上述E步和M步進行迭代直至‖θi＋1－θi‖或者‖（Q（θi＋1／θi，Y）－Ω（θ／Qi，Y）‖充分小時為止［5－7]。

2.2 EM算法在衛(wèi)星軌道計算中的應用

以12階切比雪夫多項式擬合為例。衛(wèi)星的坐標可以表示為k階切比雪夫多項式為

則誤差方程為

式中：V為誤差向量；C為擬合系數(shù)；T表示時間參數(shù)切比雪夫遞推公式；l為觀測向量（衛(wèi)星坐標）。可以表示為

精密星歷文件是以離散的形式、按一定的時間間隔（通常為15min）給出衛(wèi)星在空間的三維坐標、三維改正速度以及衛(wèi)星鐘改正數(shù)等信息。在用多項式擬合其函數(shù)模型時，可以先對精密星歷進行加密，而把加密時刻的數(shù)據(jù)當成是缺失數(shù)據(jù)，應用EM算法進行處理，從而提高擬合的精度。假設在上式中l(wèi)n為缺失數(shù)據(jù)（加密時刻的衛(wèi)星坐標數(shù)據(jù)），已知誤差向量V服從高斯正態(tài)分布，利用EM算法建立似然函數(shù)方程P（θ／Y，Z），

缺失數(shù)據(jù)ln的條件分布概率密度函數(shù)為

式中，θi為經(jīng)過i次迭代后的參數(shù)值。

則得到EM算法的期望步：

EM 算法的極大化步，將Q（θ／θi，Y）極大化，積分后對各參數(shù)求偏導數(shù)，計算得到參數(shù)θi＋1，將θi＋1代入E步進行迭代，循環(huán)直至‖θi＋1θi‖或者‖Q（θi＋1／θi，Y）－Q（θi／θi，Y）‖充分小時為止。

3 算例分析

本例選取的數(shù)據(jù)是從IGS數(shù)據(jù)中心下載的2012年1月29日的精密星歷。為了使擬合結果達到厘米級精度，這里選取1號衛(wèi)星0時0分0秒每隔30min共12個歷元的X方向的坐標值作為擬合節(jié)點的數(shù)據(jù)，而將15min間隔的數(shù)據(jù)作為已知數(shù)據(jù)進行對比。

分別采用12階切比雪夫多項式擬合、同時假設擬合節(jié)點中第五個歷元的數(shù)據(jù)為粗差或者丟失，在這種情況下降階采用11階切比雪夫多項式擬合以及采用EM算法處理后，再采用12階切比雪夫多項式擬合，計算擬合弧段內(nèi)每15min歷元節(jié)點的結果如表1所示。

表1 擬合結果精度比較（單位：m）

在缺失數(shù)據(jù)下，降階應用11階擬合與采用EM算法處理后，再采用12階切比雪夫擬合，其絕對誤差曲線如圖1所示。

圖1 絕對誤差曲線

從圖1可以看出，由于4號和5號節(jié)點間的歷元數(shù)據(jù)缺失，從而導致4號節(jié)點的擬合結果產(chǎn)生明顯的振蕩。而5號節(jié)點處于整個擬合弧段的中間，因此擬合結果較好。

4 結 論

由于精密星歷是每15min一組，數(shù)據(jù)量相對較少，所以對于短弧段用切比雪夫多項式擬合不能滿足要求。當不能采用12階切比雪夫多項式擬合時，只能采用降階擬合的方法，即采用11階或者更低階數(shù)的切比雪夫多項式來求取任意時刻的衛(wèi)星空間坐標，此時，可能導致擬合精度無法滿足要求，采用EM算法，通過對其添加與計算衛(wèi)星位置有關的“潛在數(shù)據(jù)”，能有效解決這一問題，大大提高擬合的精度。由實驗結果可以看出，在用切比雪夫多項式擬合內(nèi)插擬合弧段內(nèi)任意時刻的衛(wèi)星位置時，擬合弧段兩端節(jié)點的擬合效果較差，使用EM算法也出現(xiàn)了同樣的問題。在數(shù)據(jù)缺失的情況下，與缺失數(shù)據(jù)相鄰的節(jié)點擬合精度要相對較差。

［1]徐紹銓，張華海，楊志強，等.GPS測量原理及應用［M].武漢：武漢大學出版社，2007：11－86.

［2]李征航，張小紅.衛(wèi)星導航定位新技術及高精度數(shù)據(jù)處理方法［M].武漢：武漢大學出版社，2009：1－24.

［3]茍長龍.廣播星歷插值和精密星歷外推方法研究［D].長沙：中南大學，2009：18－22.

［4]萬亞豪，張書畢，侯東陽.GPS精密星歷插值法與擬合法的精度比較［J].全球定位系統(tǒng)，2011，36（2）：1－4.

［5]惠沈盈.觀測數(shù)據(jù)不完全的動態(tài)定位算法研究［D].長沙：中南大學，2010：12－14.

［6]林東方，宋迎春，肖琴琴.缺失數(shù)據(jù)下基于EM算法的GPS高程擬合［J].大地測量與地球動力學，2012，32（2）：1－4.

［7]楊基棟.EM算法理論及其應用［J].安慶師范學院學報·自然科學版，2009：1－5.