抗差估計(jì)在自由落體式絕對(duì)重力儀中的應(yīng)用

胡　明　張為民　楊　萌　田　蔚　鐘　敏

1　中國(guó)科學(xué)院測(cè)量與地球物理研究所大地測(cè)量與地球動(dòng)力學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，武漢市徐東大街340號(hào)，430077

?

抗差估計(jì)在自由落體式絕對(duì)重力儀中的應(yīng)用

胡明1張為民1楊萌1田蔚1鐘敏1

1中國(guó)科學(xué)院測(cè)量與地球物理研究所大地測(cè)量與地球動(dòng)力學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，武漢市徐東大街340號(hào)，430077

針對(duì)自由落體式絕對(duì)重力儀中的異常數(shù)據(jù)，提出抗差參數(shù)估計(jì)的數(shù)據(jù)處理方法。對(duì)比分析傳統(tǒng)最小二乘擬合與抗差最小二乘擬合計(jì)算的重力值，結(jié)果表明,抗差參數(shù)估計(jì)通過合理調(diào)整權(quán)重因子抑制了異常數(shù)據(jù)的影響，提高了重力測(cè)量精度。

傳統(tǒng)最小二乘估計(jì)；抗差參數(shù)估計(jì)；絕對(duì)重力儀；自由落體

FG5是目前廣泛使用的絕對(duì)重力儀產(chǎn)品，結(jié)合g9軟件，其獲得的重力測(cè)量精度達(dá)μGal級(jí)。文獻(xiàn)[1]給出了FG5重力儀g9軟件數(shù)據(jù)處理與參數(shù)擬合的一般流程；文獻(xiàn)[2]從時(shí)域上對(duì)自由落體各數(shù)據(jù)權(quán)重進(jìn)行了詳細(xì)分析與討論；文獻(xiàn)[3-6]在參數(shù)擬合公式中增加正弦項(xiàng)與余弦項(xiàng)，抑制激光器調(diào)制頻率對(duì)輸出的影響；文獻(xiàn)[7]將Lp范數(shù)理論應(yīng)用于絕對(duì)重力儀數(shù)據(jù)處理中，并詳細(xì)討論了各種噪聲模型的影響。自由落體式絕對(duì)重力儀中的時(shí)間-距離測(cè)量數(shù)據(jù)難免存在觀測(cè)粗差，一般采取的舍棄部分初始和結(jié)束段數(shù)據(jù)的預(yù)處理方法很難剔除落體數(shù)據(jù)中的全部粗差。傳統(tǒng)最小二乘參數(shù)擬合方法適用于近似服從正態(tài)分布的測(cè)量數(shù)據(jù)，參數(shù)擬合結(jié)果對(duì)非正態(tài)分布的觀測(cè)粗差極為敏感[7-8]。因此，如何抑制絕對(duì)重力儀中異常數(shù)據(jù)的影響，提高重力儀測(cè)量精度，是一個(gè)值得深入研究的問題。

本文針對(duì)自由落體式絕對(duì)重力儀，將抗差參數(shù)估計(jì)應(yīng)用到等距離間隔的數(shù)據(jù)處理中，并與傳統(tǒng)最小二乘擬合結(jié)果和g9軟件處理結(jié)果進(jìn)行對(duì)比。結(jié)果表明，抗差參數(shù)估計(jì)采用變權(quán)重的方法抑制了異常數(shù)據(jù)對(duì)重力加速度測(cè)量的影響；同等條件下，抗差參數(shù)估計(jì)得到的重力值相對(duì)傳統(tǒng)最小二乘擬合更接近g9軟件的處理結(jié)果，提高了重力值測(cè)量精度。

1　傳統(tǒng)最小二乘擬合

對(duì)于自由落體式的高精度絕對(duì)重力儀，若忽略重力梯度效應(yīng)并考慮激光器的波長(zhǎng)調(diào)制效應(yīng)，位移測(cè)量值xi與時(shí)間ti之間的關(guān)系可表示為[4-5]：

(1)

式中，x0、v0、g分別為初始位移、初始速度、當(dāng)?shù)刂亓铀俣?，M為干涉條紋分頻數(shù)，λ為激光波長(zhǎng)，ωd為激光器波長(zhǎng)調(diào)制頻率，B和C為正弦項(xiàng)和余弦項(xiàng)相關(guān)系數(shù)。

根據(jù)傳統(tǒng)最小二乘擬合標(biāo)準(zhǔn)進(jìn)行參數(shù)估計(jì)：

(2)

式中，

(3)

自由落體式FG5重力儀以等距離間隔模式工作，落體棱鏡經(jīng)過相同位移觸發(fā)相應(yīng)時(shí)刻信號(hào)，且隨落體棱鏡下落速度的增加，時(shí)刻數(shù)據(jù)越來越密集，時(shí)間間隔越來越短。若設(shè)置重力儀分頻數(shù)為4 000(即4 000個(gè)條紋產(chǎn)生一個(gè)時(shí)間測(cè)量脈沖信號(hào))，單次落體約采集180個(gè)時(shí)刻點(diǎn)，典型的時(shí)刻與時(shí)間間隔曲線如圖1所示。圖1(a)為時(shí)刻測(cè)量曲線，時(shí)刻隨點(diǎn)數(shù)增加而逐漸增大；圖1(b)為4 000分頻距離對(duì)應(yīng)的時(shí)間間隔，由于落體棱鏡下落速度越來越快，時(shí)間間隔會(huì)隨之逐漸變小。在120數(shù)據(jù)段之后的時(shí)間間隔異常跳躍，可能是自由落體系統(tǒng)震動(dòng)或時(shí)間計(jì)數(shù)器跳躍等所致，在數(shù)據(jù)處理時(shí)需要預(yù)先剔除。對(duì)多次自由落體數(shù)據(jù)進(jìn)行統(tǒng)計(jì)分析，大部分?jǐn)?shù)據(jù)在120點(diǎn)后產(chǎn)生異常跳躍，所以暫選1～100數(shù)據(jù)段進(jìn)行參數(shù)擬合。

圖1　等距離間隔模式下時(shí)間測(cè)量曲線Fig.1　Time curve for equally spaced in time

圖2　傳統(tǒng)最小二乘與抗差最小二乘的重力誤差對(duì)比Fig.2　Gravity error comparison between ordinary least square and robust trend estimation

在相同參數(shù)設(shè)置下，將傳統(tǒng)最小二乘擬合(OLS)得到的重力值與FG5軟件處理結(jié)果進(jìn)行對(duì)比。圖2中虛線為傳統(tǒng)最小二乘擬合與FG5軟件計(jì)算結(jié)果的差異。在300次自由落體數(shù)據(jù)中，大部分傳統(tǒng)最小二乘擬合計(jì)算結(jié)果與FG5軟件計(jì)算結(jié)果的差異小于0.2μGal，僅有6次自由落體顯示巨大差異。對(duì)差異較大的6次落體時(shí)間序列進(jìn)行詳細(xì)分析(圖3)，可以看出，在1～100擬合數(shù)據(jù)中出現(xiàn)異常跳躍，這些異常跳躍直接導(dǎo)致最小二乘擬合的重力值與FG5軟件計(jì)算結(jié)果產(chǎn)生巨大偏離。

圖3　6次差異較大的落體時(shí)間間隔曲線Fig.3　Time duration curves of 6 abnormal drops

2　抗差最小二乘擬合

抗差參數(shù)估計(jì)通過合理設(shè)置各數(shù)據(jù)權(quán)重來抑制異常跳躍對(duì)重力值的影響：對(duì)正常觀測(cè)值采取保權(quán)處理；對(duì)非正常但可利用的觀測(cè)值采取降權(quán)處理；對(duì)粗差，使其權(quán)重為0，在參數(shù)擬合中不使用[9-10]。抗差參數(shù)估計(jì)主要基于迭代的方法來設(shè)置各數(shù)據(jù)點(diǎn)的權(quán)重因子，且?guī)Ъ訖?quán)估計(jì)參數(shù)矩陣的計(jì)算公式可表示為[8]：

圖5　傳統(tǒng)最小二乘擬合與抗差最小二乘擬合曲線對(duì)比(a)和抗差最小二乘中各測(cè)量數(shù)據(jù)的權(quán)重(b)Fig.5　Curve fitting comparison between ordinary least square and robust trend estimation(a) and weights of the data for robust trend estimation(b)

采用抗差最小二乘處理自由落體重力儀數(shù)據(jù)，并與傳統(tǒng)最小二乘擬合曲線進(jìn)行對(duì)比(圖5(a))。圖5(b)給出了抗差估計(jì)中各數(shù)據(jù)點(diǎn)在參數(shù)估計(jì)中的權(quán)重曲線。該權(quán)重曲線顯示，自由落體測(cè)量數(shù)據(jù)中最后部分?jǐn)?shù)據(jù)的權(quán)重幾乎為0，為異常數(shù)據(jù)，與圖5(a)中曲線相符。對(duì)比傳統(tǒng)最小二乘與抗差最小二乘的殘差位移直方圖(圖6)顯示，采用抗差參數(shù)估計(jì)使殘差位移由10-4m量級(jí)減小到約10-9m量級(jí)，相對(duì)傳統(tǒng)最小二乘而言，抗差參數(shù)估計(jì)的殘差位移更接近正態(tài)分布。

圖6　傳統(tǒng)最小二乘擬合與抗差最小二乘擬合的殘差位移直方圖Fig.6　Distance residuals histogram for both ordinary least square and robust trend estimation

若進(jìn)行參數(shù)擬合的數(shù)據(jù)段中沒有異常點(diǎn)，傳統(tǒng)最小二乘法(OLS)和抗差最小二乘法(Robust)與g9軟件處理結(jié)果的差異都在0.2μGal以內(nèi)。若進(jìn)行參數(shù)擬合的數(shù)據(jù)段含有異常跳躍，利用傳統(tǒng)最小二乘法與抗差最小二乘法以及g9軟件的數(shù)據(jù)處理結(jié)果見表1。表1顯示，用傳統(tǒng)最小二乘法處理數(shù)據(jù)得到的重力值與g9軟件計(jì)算結(jié)果差異巨大；采用抗差最小二乘法時(shí)，大部分異常數(shù)據(jù)對(duì)重力加速度擬合參數(shù)的影響得到抑制，重力值的精度提高約6～8個(gè)量級(jí)。

表1　傳統(tǒng)最小二乘法(OLS)、抗差最小二乘法分別與g9軟件處理結(jié)果的對(duì)比

3　結(jié)　語

抗差參數(shù)估計(jì)應(yīng)用到自由落體式絕對(duì)重力儀的數(shù)據(jù)處理中，極大地抑制了異常點(diǎn)對(duì)重力加速度計(jì)算結(jié)果的影響，提高了重力值的精度。相對(duì)傳統(tǒng)最小二乘法而言，抗差參數(shù)估計(jì)的計(jì)算結(jié)果與g9軟件處理結(jié)果差異更小，且能容忍部分異常數(shù)據(jù)，魯棒性能良好。同時(shí)，該方法也給出了一種鑒別自由落體中異常數(shù)據(jù)的途徑。

[1]Micro-g LaCoste. G9 User’s Manual[Z]. 2012

[2]Nagornyi V D. A New Approach to Absolute Gravimeter Analysis[J]. Metrologia, 1995, 32(3): 201-208

[3]Niebauer T M, Sasagawa G S, Faller J E, et al. A New Generation of Absolute Gravimeters[J]. Metrologia, 1995, 32(3): 159-180

[4]Orlob M. Spectral Analysis of Synthetically Affected FG5 Absolute Gravimeter Residuals[D]. Dallas:University of Texas, 2011

[5]Orlob M. Absolute Gravimetry Instrumentation of FG5 and an Overview of Applications[R]. ENGO 615-Advanced Physical Geodesy Project Report, 2008

[6]Orlob M, Braun A. Impact Estimation and Filtering of Disturbances in FG5 Absolute Gravimeter Observations[J]. International Journal of Geosciences, 2013, 4(2), 302-308

[7]Nagornyi V D, Svitlov S, Araya A. Improving Absolute Gravity Estimates by theLp-Norm Approximation of the Ballistic Trajectory[EB/OL].http://arxiv.org/abs/1509.00750,2015[8]費(fèi)業(yè)泰. 誤差理論與數(shù)據(jù)處理[M]. 北京：機(jī)械工業(yè)出版社，2000 (Fei Yetai. Error Theory and Data Processing.[M]. Beijing: China Machine Press, 2000)

[9]Holland P W, Welsch R E. Robust Regression Using Iteratively Reweighted Least-Squares[J]. Communications in Statistics-Theory and Methods, 1977, 6(9): 813-827

[10] 楊元喜. 抗差估計(jì)理論及其應(yīng)用[M]. 北京：八一出版社，1993 (Yang Yuanxi. Robust Estimation Theory and Its Application.[M]. Beijing: Bayi Press, 1993)

Foundation support:National Key Scientific Instrument and Equipment Development Project, No.ZDYZ2012-1-04; National Natural Science Foundation of China, No.41504140.

About the first author:HU Ming, PhD, assistant researcher, majors in inertial sensor and gravimeter, E-mail: huming@whigg.ac.cn.

Robust Trend Estimation and Its Application in the Free Fall Absolute Gravimeter

HUMing1ZHANGWeimin1YANGMeng1TIANWei1ZHONGMin1

1State Key Laboratory of Geodesy and Earth’s Dynamics, Institute of Geodesy and Geophysics,CAS, 340 Xudong Street, Wuhan 430077, China

Robust trend estimation is employed for some outliers in the free fall absolute gravimeter. We compare gravity as obtained by ordinary least square estimation and by robust trend estimation. The results show that we can reject outliers by adjusting the weighting factors of the time-distance pairs, thus improving measurement accuracy.

ordinary least square estimation; robust trend estimation; absolute gravimeter; free fall

2015-11-13

胡明，博士，助理研究員，主要從事慣性傳感器與重力儀相關(guān)研究，E-mail：huming@whigg.ac.cn。

10.14075/j.jgg.2016.09.018

1671-5942(2016)09-0833-04

P223

A

項(xiàng)目來源：國(guó)家重大科研裝備研制項(xiàng)目(ZDYZ2012-1-04)；國(guó)家自然科學(xué)基金(41504140)。