林東方 宋迎春 金 昊
不完全測量數(shù)據(jù)的 E M處理算法*
林東方 宋迎春 金 昊
推導(dǎo)了 EM算法用于不完全測量數(shù)據(jù)的實(shí)現(xiàn)方法。結(jié)果表明,應(yīng)用 EM算法處理不完全測量數(shù)據(jù),改善了測量精度,與完全數(shù)據(jù)下的平差結(jié)果極其接近,明顯優(yōu)于不進(jìn)行平差處理的結(jié)果。
不完全測量數(shù)據(jù);平差;G-M模型;最小二乘原理;EM算法
在測量中,由于儀器設(shè)備精度、人為因素和外界環(huán)境的限制,測量中往往伴隨著誤差。為了保證測量的精度,將測量誤差控制在一定的范圍內(nèi),在必要觀測的基礎(chǔ)上常常進(jìn)行多余觀測,應(yīng)用最小二乘原理平差計(jì)算,得到最優(yōu)的估計(jì)參數(shù)。但是由于工程施工、農(nóng)業(yè)生產(chǎn)或意外事故常常使一些控制點(diǎn)被破壞或找不到,無法進(jìn)行觀測?;蛘哂^測后,由于儀器故障、人為因素造成某些測站的觀測數(shù)據(jù)被破壞或丟失。這些情況的發(fā)生,都會造成觀測數(shù)據(jù)的不完整,使完全觀測數(shù)據(jù)變?yōu)椴煌耆^測數(shù)據(jù)。由于數(shù)據(jù)的不完全,在進(jìn)行平差計(jì)算時就必須將破壞或丟失的數(shù)據(jù)刪除之后再進(jìn)行平差,這嚴(yán)重影響了觀測質(zhì)量和測量精度。當(dāng)數(shù)據(jù)缺失過于嚴(yán)重時,就無法進(jìn)行平差計(jì)算,測量數(shù)據(jù)便無法使用,觀測資料也就廢置不用,這就浪費(fèi)了大量人力、物力和財力。
鑒于此,一些學(xué)者提出,對不完全數(shù)據(jù)進(jìn)行插值處理,將缺失數(shù)據(jù)通過插值方法補(bǔ)充進(jìn)去,而后進(jìn)行平差計(jì)算,如時間序列處理的插值方法[1],這種方法雖然有一定的效果,但插值數(shù)據(jù)是對剩余數(shù)據(jù)的進(jìn)一步處理所得,并沒有利用更多的潛在信息。而在這方面,EM(Expectation-Maximization)算法就有著很大的優(yōu)勢。EM算法是用于非完全數(shù)據(jù)參數(shù)估計(jì)的一種有效方法,它是一種數(shù)據(jù)添加算法,即在觀測數(shù)據(jù)的基礎(chǔ)上添加一些“潛在數(shù)據(jù)”,利用對數(shù)據(jù)處理有益的潛在信息,得到數(shù)據(jù)缺失下參數(shù)的最優(yōu)估計(jì)。當(dāng)測量數(shù)據(jù)不完全時,可以利用觀測數(shù)據(jù)所服從的一些規(guī)律,對缺失數(shù)據(jù)的取值范圍加以限制,即采用 EM算法,利用缺失數(shù)據(jù)在現(xiàn)有條件下的期望值,對不完全數(shù)據(jù)進(jìn)行處理計(jì)算,最終得到一個較好的參數(shù)估計(jì)結(jié)果。
2.1 EM算法
EM算法的最大特點(diǎn)是通過對完全數(shù)據(jù)的處理來解決不完全數(shù)據(jù)問題[2,3]。它主要應(yīng)用于下面兩種不完全數(shù)據(jù)參數(shù)估計(jì):第一,觀測數(shù)據(jù)不完全,這是由于觀測過程的局限性或觀測后的意外事故所導(dǎo)致;第二,似然函數(shù)不是解析的,或者似然函數(shù)的表達(dá)式過于復(fù)雜從而導(dǎo)致極大似然函數(shù)難以解算。第一種情況在測量數(shù)據(jù)處理中經(jīng)常遇到。
EM算法是一種迭代算法,它的每一步迭代由兩步組成:E步(求期望)和M步(極大化)[4]。E步在給定己觀測到的數(shù)據(jù)和現(xiàn)有參數(shù)下,求“缺失數(shù)據(jù)”的條件期望;M步計(jì)算參數(shù)的MLE估計(jì),這與己知似然求參數(shù)的MLE估計(jì)的計(jì)算方法一致。
具體地講,以 P(θ/Y)表示θ的基于觀測數(shù)據(jù)的后驗(yàn)分布密度,稱為觀測數(shù)據(jù)后驗(yàn)分布。以 P(θ/Y, Z)表示添加數(shù)據(jù) (缺失數(shù)據(jù))Z后得到的關(guān)于θ的后驗(yàn)分布密度函數(shù),稱為完全數(shù)據(jù)后驗(yàn)分布。P(Z/ θ,Y)表示在給定θ和觀測數(shù)據(jù) Y下潛在數(shù)據(jù) (缺失數(shù)據(jù))Z的條件分布密度函數(shù)。我們的目的是計(jì)算觀測后驗(yàn)分布 P(θ/Y)的參數(shù)。EM的算法為,記θi為第 i+1次迭代開始時后驗(yàn)參數(shù)的估計(jì)值,則第 i +1次迭代的兩步為:
E步:將 P(θ/Y,Z)或 logP(θ/Y,Z)關(guān)于 Z的條件分布求期望,從而把 Z積掉,即
M步:將 Q(θ/θi,Y)極大化,即找一個點(diǎn)θi+1,使
EM算法在每一次迭代后均提高觀測極大似然密度函數(shù)值,具有良好的全局收斂性[4,5]。
2.2 (G-M)模型下不完全測量觀測數(shù)據(jù)的EM處理方法
測量中常采用 Gauss-Markov(G-M)模型[6]
在測量平差中,間接平差的誤差方程為[7]
這里假設(shè)觀測數(shù)據(jù) ln由于意外事故丟失,應(yīng)用 EM算法建立似然函數(shù)方程 P(θ/Y,Z)
式中θ為未知參數(shù) X,Y為不完全觀測數(shù)據(jù) (l1,l2,…,ln-1),Z為缺失數(shù)據(jù) ln。
由于測量誤差改正數(shù) V服從期望為零的正態(tài)分布,因此可以得到缺失數(shù)據(jù) ln的條件分布概率密度函數(shù)為
由式(6)、(7)便得到 EM算法的期望步
以水準(zhǔn)測量為例。在圖 1中,A、B為已知水準(zhǔn)點(diǎn),其高程 HA=12.013 m,HB=10.013 m,可視為無誤差,為了確定 C及D點(diǎn)的高程,共觀測了 4個高差,h1=-1.004 m,h2=1.516 m,h3=2.512 m,h4=1.520 m。
當(dāng)測量數(shù)據(jù)完全時,即高差觀測值h1、h2、h3、h4沒有丟失,我們采用最小二乘原理,應(yīng)用間接平差法,計(jì)算得到 HC、HD的高程為:HC=11.008 3 m,HD=12.525 7 m。
圖1 水準(zhǔn)測量路線Fig.1 Leveling line
假設(shè)由于意外事故,h2、h4數(shù)據(jù)丟失,測量數(shù)據(jù)變?yōu)椴煌耆珨?shù)據(jù),由于 h2、h4的丟失,沒有了多余觀測,傳統(tǒng)的平差方法失效,直接計(jì)算得到 HC、HD的高程為:HC=11.009 m,HD=12.525 m。
在非完全數(shù)據(jù)的情況下,兩種計(jì)算方法所得結(jié)果與完全數(shù)據(jù)情況下平差計(jì)算所得結(jié)果進(jìn)行比較(表 2)。不平差所得參數(shù) HC、HD的值與完全數(shù)據(jù)下間接平差所得參數(shù) HC、HD的值均相差 0.7 mm,而 EM算法所得參數(shù) HC、HD的值與完全數(shù)據(jù)下間接平差所得參數(shù)值均相差 0.2 mm。由此可見,在測量數(shù)據(jù)不完全時,采用 EM算法進(jìn)行數(shù)據(jù)處理能夠得到更好的結(jié)果。
表1 高程點(diǎn)計(jì)算結(jié)果(單位:m)Tab.1 Elevation results(un it:m)
表2 高程點(diǎn)計(jì)算結(jié)果比較(單位:m)Tab.2 Comparison between elevation results(un it:m)
推導(dǎo)了在測量數(shù)據(jù)不完全時,采用 EM算法進(jìn)行數(shù)據(jù)處理的實(shí)現(xiàn)方法,通過實(shí)例,驗(yàn)證了 EM算法在處理不完全測量數(shù)據(jù)時的有效性。結(jié)果表明,在測量數(shù)據(jù)缺失,沒有多余觀測的情況下,使用 EM算法的處理結(jié)果與完全數(shù)據(jù)時平差結(jié)果極其接近,效果十分明顯。在 EM算法下,許多以往廢置的數(shù)據(jù)均可以重新啟用,例如,往年廢置變形監(jiān)測數(shù)據(jù),往年廢置水準(zhǔn)網(wǎng)觀測數(shù)據(jù),有了這些歷史數(shù)據(jù)的加入,充分利用過去的觀測數(shù)據(jù)[8],對分析高層建筑物、水壩形變,地殼形變具有不可估量的意義。
EM算法是一種迭代收斂算法,它的收斂速度較慢[9],參數(shù)初值的選取對迭代的收斂速度影響較大,因而在初值選取時,我們應(yīng)利用參數(shù)的先驗(yàn)信息,選取一個較接近參數(shù)估計(jì)值的初值,以提高 EM算法的收斂速度。EM算法在處理單個缺失數(shù)據(jù)或連續(xù)缺失數(shù)據(jù)時的效果顯著,但是 EM算法在處理重度缺失數(shù)據(jù)時,處理效果不明顯[10],在處理離散缺失數(shù)據(jù)時,須分批處理,因此如何提高重度缺失數(shù)據(jù)時 EM算法的處理效果和系統(tǒng)整體處理離散缺失數(shù)據(jù)尚需進(jìn)一步的研究。應(yīng)用 EM算法要求觀測值與參數(shù)之間具有一定的分布關(guān)系,當(dāng)觀測值與參數(shù)之間無分布關(guān)系時則不適用 EM算法。
1 武艷強(qiáng),黃立人.時間序列處理的新插值方法[J].大地測量與地球動力學(xué),2004,(4):43-45.(Wu Yanqiang and Huang Liren.A new interpolationmethod in time series analyzing[J].Journal of Geodesy and Geodynamics,2004, (4):43-45)
2 DempsterA P,Laird N M and Rubin D B.Maximum likelihood from incomplete data via the EM algorithm[J].Journal of the Royal Statistifcal SocietyB,1977,39:1-38.
3 錢俊,舒寧.基于 EM算法和單幅雷達(dá)圖像陰影的控制點(diǎn)坡度校正[J].武漢大學(xué)學(xué)報 (信息科學(xué)版),2004,(12):1 089-1 092.(Qian Jun and Shu Ning.Correction of control point slope based on EM algorithm and shading of single SAR image[J].Geomatics and Information Science ofWuhan University,2004,(12):1 089-1 092)
4 Graham C G and Juan C A.Approximate EM algorithms for parameter and state estimation in nonlinear stochastic models [A].Proceedings of the 44th IEEE Conference on Decision and Control,and the European Control Conference[C]. 2005,12-15.
5 王兆軍.EM算法收斂的必要條件[J].南開大學(xué)學(xué)報(自然科學(xué)版),1994,(2):85-88.(Wang Zhaojun.The necessary condition on the convergence of the EM algorithm [J].Acta Scientiarum Naturalium Universitatis Nankaiensis,1994,(2):85-88)
6 郭金運(yùn),靳奉祥,劉國林.不完全測量數(shù)據(jù)的應(yīng)用研究[J].測繪通報,2002,(2):7-8.(Cuo Jinyun,Jin Fengxiang and Liu Guolin.Study ofApplication of incomplete data in surveying[J].Bulletin of Surveying andMapping,2002, (2):7-8)
7 武漢大學(xué)測繪學(xué)院測量平差學(xué)科組.誤差理論與測量平差基礎(chǔ) [M].武漢:武漢大學(xué)出版社,2003.(Research Group of Surveying Adjustment,School of Surveying and Mapping,Wuhan University.The base of errors theory and surveying adjustment[M].Wuhan;WuhanUniversity Press, 2003)
8 熊俊楠,等.高層建筑物變形監(jiān)測數(shù)據(jù)處理方法研究及工程應(yīng)用[J].測繪科學(xué),2010,(7):69-71.(Xiong Junnan, et al.Research and application for high-rise deformation observation data processing method[J].Science of Surveying andMapping,2010,(7):69-71)
9 高旅端,陳志,王家潤.一種加速 EM算法收斂的方法[J].數(shù)理統(tǒng)計(jì)與應(yīng)用概率,1998,(4):342-348.(Gao Luduan,Chen Zhi andWang Jiarui.A method for accelerating convergence of fhe E M algorithm[J].Mathematical Statistics and Applied Probability,1998,(4):342-348)
10 張香云,張秀偉.不同缺失率下 EM算法的參數(shù)估計(jì)[J].數(shù)理統(tǒng)計(jì)與管理,2008,(3):428-431.(Zhang Xiangyun and Zhang Xiuwei. Parametrical estimation for EM algorithm in the case of different losing-rate[J].Application of Statistics and Management,2008,(3):428-431)
EM PROCESSING ALGORITHM OF INCOM PLETE SURVEY ING DATA
Lin Dongfang,Song Yingchun and Jin Hao
(School of Geosciences and Info-Physics,Central South University,Changsha 410083)
For the incomplete surveying data caused by various factors,traditional adjustment can not be achieved.It is har mful to the accuracy of the surveying,even leads to abandoned observation data.Through EM algorithm,we can take the potential infor mation which is helpful to i mprove the accuracy of data processing effectively.EM processing algorithm for incomplete surveying data is proved.The results show that the surveying accuracy is improved by using the EM algorithm to process the incomplete data.The results are extremely si milar to the adjustment resultswith complete data and significantly better than the resultswithout adjustment.
incomplete surveying data;adjustment;G-M model;principle of least squares;EM algorithm
(中南大學(xué)地球科學(xué)與信息物理學(xué)院測繪與國土信息工程系,長沙 410083)
1671-5942(2011)04-0112-04
2011-01-16
國家自然科學(xué)基金 (40874005);教育部博士點(diǎn)基金(200805331086)
林東方,男,1986年生,碩士研究生,主要研究方向?yàn)楝F(xiàn)代測量數(shù)據(jù)處理.E-mail:lindongfang223@163.com
P207;P203
A