EM算法在廣播星歷計算衛星位置中的應用

肖琴琴，宋迎春，杜 琨

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

在衛星導航定位中，為了確定用戶的位置，必須首先通過衛星星歷計算衛星的位置，因此，獲取GPS衛星的位置是求得接收機坐標的必要環節。目前，獲得衛星位置的途徑主要有兩種：一種是通過廣播星歷計算得到，一種是通過事后的精密星歷直接獲取。廣播星歷是在接收機觀測時接收衛星發播的導航電文經譯碼取得的實時星歷，而精密星歷是由有關單位提供的事后處理星歷。其中廣播星歷的精度偏低［1］，約為2m，而由IGS提供的精密星歷精度較高，優于5cm，但是精密星歷需觀測工作結束11d后才可用［2］。GPS廣播星歷雖然精度較精密星歷低，但它具有實時、易獲取的特點，常在動態實時導航定位中得到應用。因此，對實時性要求較高的定位工作而言，只能通過廣播星歷獲得衛星位置。

使用廣播星歷直接計算衛星位置可以通過文獻［3］中的公式進行求解，然而當采樣頻率較高需多次計算衛星的位置時，直接計算的方法就會占用很多的內存并耗費很長時間?？紤]到衛星位置是一個連續變化的過程，在一定時間段內可將衛星坐標表示為時間多項式，在內存中僅保存多項式的系數，那么在計算衛星坐標時，只要調出多項式系數即可。這與按固定公式計算衛星坐標的方法相比，大大提高了數據處理的效率。

由于衛星遮擋、衛星信號傳播限制等原因，可能造成部分有用信息丟失，導致數據不完全，或者計算出來的星歷參數值具有一定的誤差，從而使廣播星歷計算出的衛星坐標與其真值之間有較大的偏差，另外每計算一個觀測歷元的衛星坐標也需要較大的計算量，因此如何提高衛星坐標計算的效率及精度是值得探討的問題。EM（expectation maximization）算法是近年發展很快且應用很廣的一種算法，它能夠在觀測數據的基礎上增加一些“潛在數據”，從而簡化計算并完成一系列簡單的極大化或模擬，因此可以考慮將EM算法應用到衛星的軌道標準化。由于不同衛星的精度等級不一樣，因此在選擇擬合點的個數以及擬合階數時也會有所差別。當發現某些擬合點的誤差太大時，如果直接刪除進行降階處理必然導致誤差增大，這時如利用EM算法在觀測數據的基礎上加一些“潛在數據”，可使得擬合的精度和可靠性仍能滿足要求。

1 常用的插值法和擬合法

利用多項式進行插值擬合通常有很多種方法，如拉格朗日多項式插值、埃爾米特插值、普通的多項式擬合、勒讓德多項式擬合和切比雪夫多項式擬合。但目前最常用的多項式是切比雪夫多項式，關于切比雪夫擬合的原理已在文獻［4］詳細闡述。

用n階切比雪夫多項式來逼近時間段［t0，t0＋Δt］中的衛星星歷時，先要將變量t∈［t0，t0＋Δt］變換為變量τ∈［－1，1］：

于是，衛星坐標可表示為

式中：n為多項式的階數，Cxi，Cyi，Czi為切比雪夫多項式的系數。

切比雪夫多項式Ti的遞推公式為

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

2 EM算法及其在衛星軌道計算中的應用

廣播星歷文件經常會由于衛星信號限制等原因，導致部分數據失真或缺失，在數據缺失或數據量較少的情況下，會導致坐標的計算精度降低。采用EM算法，添加與衛星軌道信息相關的“潛在數據”，可以有效解決這一問題，大大提高衛星坐標擬合的精度。

2.1 EM算法的原理

EM算法［5］是一種迭代算法，它的每一次迭代由兩步組成：E步（求期望）和M步（極大化）［6］。一般地，以表示θ的基于觀測數據的后驗分布密度，稱為觀測數據后驗分布；以P（θ/Y，Z）表示添加數據（缺失數據）Z后得到的關于θ的后驗分布密度函數，稱為完全數據后驗分布。P（Z/θ，Y）表示在給定θ和觀測數據Y下潛在數據（缺失數據）Z的條件分布密度函數，其目的是計算觀測后驗分布的參數。

EM算法的進行如式（4）、式（5）所示。記θi為第i＋1次迭代開始時后驗參數的估計值，則第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）是在整個參數空間Ω內使得Q（θ/θi，Y）最大的θ的值所組成的集合。將上述E步和M步進行 迭 代 直 至 ‖θi＋1－θi‖ 或 者充分小時為止［7－9］。

2.2 EM算法在廣播星歷計算衛星軌道中的應用

這里以8階切比雪夫多項式擬合為例，衛星的坐標可以表示為8階切比雪夫多項式：

則誤差方程為

廣播星歷文件是每隔2h進行一次更新，在用多項式擬合其函數模型時，可先對某個時刻的廣播星歷進行外推1h，把某些加密時刻的數據當成是缺失數據，應用EM算法進行處理，可提高擬合的精度。假設在式（8）中lk為缺失數據（加密時刻的衛星坐標數據），已知誤差向量V服從高斯正態分布，利用EM算法建立似然函數方程：

缺失數據lk的條件分布概率密度函數為

式中：θi為經過i次迭代后的參數值。

則得到EM算法的期望步：

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

3 算例分析

由于GPS廣播星歷的外推時間一般為1h［10］，因此，只對廣播星歷中的某一歷元前、后各1h內的衛星位置進行擬合。本文采用從Internet網上下載的GPS衛星2010－09－05的廣播星歷abpo2480.10n數據（RINEX格式），對19號衛星在1：00－3：00時段內X方向上的坐標值進行擬合，并將擬合出的結果分別與直接由廣播星歷和IGS精密星歷計算出的軌道進行對比。由于內插點較多，本文僅選取靠近缺失值附近的8個點值來進行對比。計算結果見表1、表2和圖1、圖2，其中表1與圖1是擬合值與廣播星歷計算值之間的比較結果，表2與圖2是擬合值與IGS精密星歷的比較結果。

圖2 與精密星歷計算值的絕對誤差曲線

表2 與IGS擬合結果精度比較 m

由計算結果可以看出：

1）當用于切比雪夫多項式擬合的數據精度較好且數據完整的情況下，使用切比雪夫多項式內插獲得的衛星位置和直接使用廣播星歷用戶算法計算出的結果非常接近，且其誤差在毫米級以內，這充分說明使用切比雪夫多項式來擬合衛星軌道是合理的，如表1和圖1所示。

圖1 與廣播星歷計算值的絕對誤差曲線

表1 與廣播星歷直接擬合結果精度比較 m

2）在缺失數據的情況下，使用降階的切比雪夫多項式計算出的衛星位置將會產生明顯的偏差，且在缺失點附近尤為顯著，偏差大小達到數據完整時的1000多倍，如表1所示。這也說明，在數據缺失時有必要采取一定的措施來提高切比雪夫多項式擬合的精度和可靠性。

3）通過引入EM算法來生成缺失數據的潛在值，解決了因數據不足需要對切比雪夫多項式進行降階處理的問題，且使用生成的潛在數據與其它數據一起擬合出的衛星位置的精度和可靠性與數據完整時基本相當，如表1和圖1所示。同時也表明，在數據缺失的情況下，使用EM算法來生成缺失數據的潛在值是可行的，且其精度和可靠性能夠滿足需求。

4）將表2、圖2與表1、圖1進行比較可知，無論是將直接使用GPS廣播星歷用戶算法計算出的衛星位置作為參考值，還是將IGS精密星歷內插出的衛星位置作為參考值，其結論都是相同的，只不過前者的效果尤為顯著，其原因是廣播星歷本身的精度要明顯地低于IGS精密星歷。

4 結束語

在動態導航定位中，數據的采樣頻率越來越高，如果直接使用廣播星歷用戶算法來計算衛星的位置，必然會大大增加導航定位計算的內存或負荷，從而導致實時性難以滿足。此時，選擇較少的采樣點并使用廣播星歷用戶算法計算相應時刻的衛星位置，再用切比雪夫多項式擬合，將會節約導航定位計算占用的內存并減少計算負荷。但是，由于各方面因素的影響，直接使用廣播星歷計算出的衛星位置個數及精度會受到一定的影響，導致部分所需節點的觀測數據缺失或不可用，此時必須使用降階的切比雪夫多項式進行擬合，其精度和可靠性將會大大降低。本文引入的EM算法能夠生成缺失或不可用節點位置的潛在值，彌補了因數據缺少切比雪夫多項式需要降階處理的不足，且使用EM算法后的切比雪夫多項式不需降階，其擬合結果的精度和可靠性與數據完整時基本相同。因此，EM算法確實是一種解決部分衛星軌道數據缺失或不可用時獲取任意時刻的衛星位置的較好方法。

［1］邱蕾，廖遠琴，花向紅.基于IGS精密星歷的衛星坐標插值［J］.測繪工程，2008，17（4）：15－18.

［2］樓益棟，劉萬科，張小紅.GPS衛星星歷精度分析［J］.測繪信息與工程，2003，28（6）：4－6.

［3］徐紹銓，張華海，楊志強.GPS測量原理及應用［M］.武漢：武漢大學出版社，2002.

［4］余鵬，孫學金，趙世軍.GPS定位中衛星坐標計算的切比雪夫多項式擬合法［J］.氣象科技，2004，32（3）：198－201.

［5］DEMPSTER，A.P，LAIRD，N.M.，and Rubin，D.B.Maximum Likelihood From Incomplete Data Via the EM Algorithm［J］.Journal of the Royal Statistifcal Society B，1977，39（1）：1－38.

［6］GRAHAM C.G，JUAN C.A，Approximate EM Algorithms for Parameter and State Estimation in Nonlinear Stochastic Models［J］.Proceedings of the 44th IEEE Conference on Decision and Control，and the European Control Conference，2005：368－373

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

［8］楊基棟.EM算法理論及其應用［J］.安慶師范學院學報：自然科學版，2009，15（4）：30－35

［9］曾傳璜，張鑫，張晶晶，等.EM算法的研究［J］.軟件導刊，2008，7（9）：97－98

［10］羅力，黃聲享.單組廣播星歷精度分析及其衛星軌道擬合研究［J］.測繪信息與工程，2010，35（1）：21－22.