吳 瓊,滕云田,王曉美
(中國地震局地球物理研究所,北京 100081)
高精度絕對重力測量技術可以實現(xiàn)地球重力場及其隨時間變化和空間分布的μGal 量級測量(1 μGal= 10-8m/s2),是重力匹配導航[1]、深空探測[2]等領域最重要、最基本的測量技術之一。雖然現(xiàn)階段國內的絕對重力測量儀器還依賴進口,但隨著碘穩(wěn)頻激光器研制、高精度時頻測量、真空技術以及精密機械加工等領域的不斷進步,自主研發(fā)的高精度絕對重力儀正逐漸成熟并向實際應用領域拓展[3]。
絕對重力儀一般采用分組多次測量求平均的方式獲得最終測量結果。在實際測量過程中,儀器自身的系統(tǒng)振動和測點的背景振動對測量結果的影響比較明顯,是引起測量結果離散度變大和導致測量結果中出現(xiàn)異常點的主要原因[4]。對單次測量的時間-位移坐標進行最小二乘擬合時,由于測量誤差的存在,同樣會導致測量結果出現(xiàn)異常點[5]。此外,激光器和銣原子頻標的電子系統(tǒng)長期連續(xù)工作時可能存在的瞬時擾動,使得提供的長度基準和時間基準會發(fā)生小概率的跳變,從而造成單個測量結果異常甚至突變。
目前針對測量過程中的異常點(離群程度較?。┖屯蛔凕c(離群程度較大),一般采用的處理方式是在組測量結果計算時,采用基于正態(tài)分布的一元異常點檢測算法進行異常點的檢測和舍棄。但實際測量結果表明,在單組100 次下落的絕對重力測量過程中,如果異常點出現(xiàn)的個數(shù)較多(大于5 個)或離群程度不一樣的時候,該算法會造成某些異常點的漏檢,降低組測量結果的測量精度并可能在最終測量結果中引入測量偏差。
隨著人工智能(AI)技術的不斷進步,異常點檢測作為一個特定的研究方向,在過程工業(yè)、金融業(yè)、通信領域得到普遍關注和廣泛應用[6]。異常點檢測算法一般基于某個存在的模型,定義一個區(qū)間,處于區(qū)間外的點被定義為異常點。根據(jù)模型的定義可以將異常值檢測算法分為基于概率分布模型[7]、基于距離和聚類[8]、基于數(shù)據(jù)分類[9]等。
考慮到在自主研發(fā)的絕對重力儀測量過程中,單組測量數(shù)據(jù)一般會被預先設定為50 或100 個,但要求快速在線完成異常值檢測,故本文利用基于歐式距離的局部異常因子LOF 算法(Local Outlier Factor)[10]對組內測量數(shù)據(jù)進行異常值檢測。通過對LOF 算法關鍵唯一參數(shù)——“鄰域寬度”的選取及模擬仿真,并與傳統(tǒng)的基于正態(tài)分布的一元異常點檢測算法進行對比確定,基于LOF 算法的異常值檢測可以實現(xiàn)對組內測量異常點在線、快速、精準檢測,同時還可以消除漏檢和誤檢現(xiàn)象,從而提高自主研發(fā)絕對重力儀測量的精度和準確度。
目前絕對重力測量中,針對測量過程中出現(xiàn)的異常值一般采用一元正態(tài)分布異常值檢測算法進行判斷和排除,基本算法原理如下:
對于一次測量獲得的n個測量結果中,定義第i個測量結果xi的異常值Zi為:
其中,為本組測量結果的平均值,σ為本組測量的標準差。
如果Zi≥ 3,則認為該值為異常值,在本次測量結果中應予以舍棄。但是在自主研發(fā)的絕對重力儀測量過程中,這種檢測算法很容易造成異常值的漏檢。為了說明這一現(xiàn)象,首先構建一組測試用的數(shù)據(jù)集,包含異常點和突變點,具體如下:
1)構建數(shù)據(jù)集,包含100 個數(shù)據(jù),分布區(qū)間為400~600,單位為μGal,服從標準正態(tài)分布;
2)在第10、30、49、51、90 位置的測量結果上分別疊加600、-600、-400、-400、400,構建異常點;
3)在50 和70 次測量結果上疊加-2500 和1500,構建突變點。
如圖1所示,當不考慮50、70 位置的突變點,對構建的數(shù)據(jù)集使用一元異常值檢測算法時,設置的異常點被正常檢出。但是當數(shù)據(jù)集的第50 和70 位置被設置為突變點后,第10、30、49、51、90 位置的異常點均被認為是正常值保留??梢哉f當測量的數(shù)據(jù)中出現(xiàn)突變點時,離群程度較小的異常點被漏檢并參與最終測量結果的計算。
這種現(xiàn)象出現(xiàn)在項目組自主研發(fā)絕對重力的測量 結果中,直接影響了最終測量結果的準確度和精度。故需要引入新的異常值檢測算法,不僅可以完成對離群程度較大的突變點進行檢測,同時還要避免對離群程度較小、連續(xù)出現(xiàn)的異常點的漏檢。本文選用計算速度快,檢測效果好的局部異常因子LOF 算法。
局部異常因子LOF 算法的計算流程如下:
Step 1:定義測試數(shù)據(jù)集內的 p(x1,y1)點和 q(x2,y2)點之間的距離:
Step 2:定義點p 的鄰域 D { p} 寬度k:
對于點o ∈ D:
1)至少有k個點 o' ∈ D { p},滿足:
2)至多有k-1 個點 o' ∈ D { p},滿足:
Step 3:計算點p 的可達距離:
Step 4:計算點p 局部可達密度:
Step 5:計算點p 的局部可達因子LOF:
從LOF 的定義可以看出,在點p 的鄰域內,其他點的局部可達密度大,而點p 的局部可達密度小,則點p的局部異常因子就大,即點p是一個局部異常值。
Breunig 等[10]還證明了在測量結果組成的數(shù)據(jù)集 合中,正常值的LOF 一般在1 左右并小于2,因此可 以將LOF 值大于2 的作為異常值或突變值。
從1.2 節(jié)LOF 的定義可以看出,對測量數(shù)據(jù)的每個測點計算其LOF 時,只有1 個參數(shù),即鄰域寬度k影響最終檢測結果。因此本部分通過數(shù)值模擬確定k在利用LOF 算法完成絕對重力測量數(shù)據(jù)檢測時的取值。
為了更好地測試k對最終檢測結果的影響程度,基于1.1 節(jié)構建的測試數(shù)據(jù),主要針對以下幾方面進行LOF 算法的檢測效果測試:
1)存在突跳點的情況下離群程度較小的異常點檢測;
2)存在連續(xù)異常、突跳情況的異常點檢測;
3)原始數(shù)據(jù)標準差較大情況下異常點檢測。
不同k值對應的各個測點的LOF 計算結果如圖2所示。LOF 的大小以不同的顏色標識,可以看出隨著鄰域寬度k取值從0 到100 的變化,LOF 算法的檢測效果的變化:
k=1~4:漏檢率和誤檢率較高;
k=5~20:對1)、2)、3)項內容均可以正確檢測,各個異常點的LOF 值明顯可辨(均大于20);
k=21~100:出現(xiàn)不同程度的漏檢現(xiàn)象。
圖3顯示了不同的參數(shù)k對應的異常點檢測效果。其中第1 和第2 條曲線分別顯示了將檢測出的異常點刪除后得到的新數(shù)據(jù)集的均值和標準差,第3 條曲線顯示的是檢出的異常點個數(shù)。
測試數(shù)據(jù)集共設置7 個異常點,從圖2和圖3的第3 條曲線可以明顯看出,在5≤k≤60 的區(qū)間范圍內,7 個異常點可被準確檢出,包含第49、50 和51三個連續(xù)異常和突變點,異常點檢出后的均值和方差 與模擬數(shù)據(jù)預設的均值和方差相比均小于5 μGal。
圖2 不同k 值對應的各個測點的LOF Fig.2 Local outlier factor of each measurement point corresponding to different k values
圖3 不同k 值時異常點檢測效果 Fig.3 Effect of outlier detection under different k values
在實際的絕對重力測試過程中,每組測量次數(shù)一般設置為50~100,24 h 完成48 組測量。
綜合圖2、圖3,考慮運算的速度、異常點LOF 值的可辨識性和算法檢測的準確性,鄰域寬度k的取值設定為10,并要求測量時每組測量次數(shù)不少于20 次。
自主研發(fā)的激光干涉絕對重力儀于2018年2月9日至4月1日,在中國計量科學研究院昌平園區(qū)2017國際絕對重力比對基地的8 號測點進行連續(xù)觀測。測試期間設置每組測量50 次,組間時間間隔是4 h,共完成14 650 次測量。測量結果的準確度和測量精度依據(jù)《中國大陸構造環(huán)境檢測網(wǎng)絡技術規(guī)程—絕對重力測量部分》之“成果計算”部分完成計算。
圖4所示為連續(xù)2 組,共100 個原始測量結果、一元異常值檢測結果和本文討論的LOF 算法完成的異常點檢測結果(k=10,LOF 大于2 認為是異常值[10])。
從圖4可以看出,LOF 算法將負向離群程度較大的第27 點作為異常值準確檢出,但是一元異常值檢測算法則將其作為有效值保留,這樣在計算這100 個測量結果的平均值時,基于一元異常值檢測算法的結果肯定小于基于LOF 算法檢測的結果,影響了測量結果的準確度。
圖4 實測數(shù)據(jù)的異常值檢測算法對比 Fig.4 Comparison on outlier detection algorithms based on measured data
圖5是分別利用一元異常值檢測算法和LOF 算法,按照50 個測量結果為一組計算絕對重力組測量結果的誤差棒圖對比??紤]到在測量過程中,2018年2月 12日河北廊坊永清縣Ms4.3 地震對儀器參考棱鏡的影響造成的測量結果異常,圖5中顯示的是本次連續(xù)測量期間震后2018年2月22日至4月13日共計50 天的連續(xù)測量數(shù)據(jù)。可以看到,LOF 算法完備地完成測量結果中異常值的檢測,組測量結果的誤差棒圖沒有出現(xiàn)異常結果,200 組的測量結果為980XXX783.0 μGal(儀器測量高度位置),測量精度為1.94 μGal。
基于一元異常值檢測的組測量結果中,第138 組 和152 組出現(xiàn)異常結果,表明這兩組內有較為明顯的離群點未被正確檢出,且參與了組結果的計算,200組的測量結果為980XXX782.39 μGal(儀器測量高度位置),測量精度僅為427.15 μGal,遠大于利用LOF算法進行異常值檢測后的計算結果。
圖6為200 組測量結果,利用一元異常值檢測算法和LOF 算法進行異常值檢測后單組測量精度,第138 組和152 組內由于存在未被檢出的異常點,造成 單組測量精度偏差。而對其他的單組測量精度對比,基于LOF 算法的組測量結果精度與一元異常值檢測算法的結果相比,平均小9.37 μGal。
圖5 兩種異常值檢測算法完成異常值檢測后組測量結果的誤差棒圖對比 Fig.5 Comparison on error bars of the two outlier detection algorithms
圖6 兩種異常值檢測算法完成異常值檢測后單組測量精度對比 Fig.6 Comparison on precisions of the two outlier detection algorithms
本文首先描述了目前自主研發(fā)絕對重力儀中普遍使用的基于正態(tài)分布的一元異常值檢測算法的漏檢現(xiàn) 象。然后引入LOF 算法,通過對同一組測試數(shù)據(jù)集進行計算的結果表明,LOF 算法規(guī)避了存在離群程度不同的異常點時的漏檢現(xiàn)象,同時對連續(xù)出現(xiàn)的異常點和突跳點均可以準確檢出。該算法非常適合自主研發(fā)絕對重力儀組測量結果的異常值檢測。
從LOF 算法的定義可以看出,影響該算法計算結果的唯一參數(shù)是鄰域寬度k。對此利用構建的數(shù)據(jù)樣本進行數(shù)值模擬計算的結果表明,對于100 個數(shù)據(jù)樣本組成的數(shù)據(jù)集,5≤k≤60 的區(qū)間范圍內可以檢出離群程度不同的全部異常點,無漏檢現(xiàn)象??紤]實際測量中單組測量次數(shù)一般為50~100 次,從計算時間和異常值檢測效果的綜合考慮,取k=10。
最后利用LOF 算法,取鄰域寬度k=10,對自主研發(fā)絕對重力儀長期連續(xù)觀測數(shù)據(jù)進行處理,對比了一元異常值檢測算法和LOF 算法的檢驗效果。結果表明,利用LOF 算法可以在線、快速、精準地完成異常值檢測,明顯提高總體計算結果的精度,同時不會影響測量結果的準確度。在下一步工作中考慮將LOF 算法代替一元異常值檢測,作為自主研發(fā)絕對重力儀單組異常值檢測的通用算法。