張鳴之 王鑫宇 趙文祎 吳明魁
1 清華大學工程物理系,北京市雙清路30號,100084 2 中國地質環境監測院,北京市大慧寺20號,100081 3 自然資源部地質災害智能監測與風險預警工程技術創新中心,北京市大慧寺20號,100081 4 武漢大學測繪學院,武漢市珞喻路129號,430079 5 中國地質大學(北京)信息工程學院,北京市學院路29號,100083 6 中國地質大學(武漢)地理與信息工程學院,武漢市魯磨路388號,430074
全球導航衛星系統(global navigation satellite system,GNSS)具有連續高精度、全天候監測、可直接提供三維變形信息等優勢,已廣泛應用于橋梁、水庫大壩、尾礦庫以及滑坡、沉降等地質災害的變形監測[1-2]。然而,受數據采集和通訊、觀測值粗差、多路徑效應等GNSS觀測值質量因素以及GNSS定位軟件性能的影響,GNSS精密定位獲取的三維位移監測序列中不可避免地存在粗差。因此,識別并剔除GNSS監測位移數據中的異常值,對提升GNSS監測的準確性,降低自動化監測預警的誤報、漏報概率具有重要意義。
異常檢測是數據挖掘的基礎問題之一,常用的方法有基于數理統計、聚類分析、數據挖掘、機器學習的方法等。針對GNSS變形監測時間序列中的粗差處理問題,吳浩等[3]提出一種基于小波分析的改進型3σ(本文簡稱小波-3σ)粗差探測方法,對GNSS位移序列中3σ~5σ粗差的探測效果較好,但難以剔除序列中的較小粗差;文獻[4-5]分別采用奇異譜分析、小波分析方法識別變形監測序列中的粗差。聚類方法中最常用的算法為K-means算法,其是利用異常值離群的特性,但聚類結果易受初始隨機聚類中心和樣本離群點的影響[6]。Hochreiter等[7]提出長短期記憶網絡算法(long short-term memory,LSTM),因其在處理單變量時序數據方面的優勢,近年來諸多學者提出基于LSTM的異常檢測算法。然而,上述大多數異常檢測算法為事后處理模式,難以滿足GNSS位移序列粗差的實時性、可靠性檢測需求。因此,需要一種初始化簡易、能實時適應數據變化、時間復雜度小的GNSS位移監測數據異常值實時檢測算法。
穩健隨機分割森林(RRCF)算法是一種面向動態數據流的異常值檢測算法[8]。該算法在數據流上動態維護隨機分割樹的分布,并根據新樣本點的插入或刪除對分割樹復雜度的影響來評判其異常程度。RRCF算法能通過一次遍歷樣本集獲得每個樣本點的異常評分,具有較好的運行效率,已在城市交通、水質監測、互聯網性能指標分析等領域得到廣泛應用并取得較好的效果[9-11]。
本文旨在研究基于RRCF的GNSS監測位移數據異常實時檢測方法。首先介紹RRCF算法的原理,然后介紹基于RRCF的GNSS監測位移數據異常實時檢測算法流程,最后通過仿真和實測數據評估算法性能。
RRCF算法的基本原理是在數據流上動態維護二叉樹,即穩健隨機分割樹(robust random cut tree,RRCT)結構的分布。給定樣本集S生成的樹T(S),樣本集中每個元素對應樹T(S)中的葉子節點,非葉子節點為其左右子節點的集合,根節點對應樣本集S,S上的RRCT生成如下[8]:
1)假設樣本數據為{x1,…,xi};
2)選取Xi,其中Xi符合[minx∈Sxi,maxx∈Sxi]上的均勻分布;
3)將S分割為左右子樹S1={x|x∈S,x≤Xi}、S2=S-S1,并對S1、S2進行遞歸分割,直到所有樣本點均在葉子節點上。
RRCF就是一組相互獨立的穩健隨機分割樹(RRCT),其核心問題是異常定義以及RRCF在數據流上的動態維護[8]。
給定樣本集S構建的樹Tree1,將x點從Tree1中插入或刪除得到的樹Tree2與直接由樣本集S±{x}構建的樹Tree3的概率分布一致。因此,要計算插入或刪除某個節點引起的樹復雜度變化,只需將節點直接從樹上插入或刪除,而無需使用新的數據集來構建樹[8]。這是RRCF適用于數據實時處理的理論依據。
若RRCF復雜度隨著點的加入而顯著增加,則該點為異常點。將某個樣本點刪除或插入后,模型的復雜度變化量(如所有葉子節點深度之和的變化量)為該點的異常程度。
給定樣本點集S和樣本點y(y∈S),定義f(y,S,T)為y在隨機分割樹T上的路徑長度。刪除點x后產生的樹為T(S-{x}),此時,y在樹T(S-{x})上的路徑長度為f(y,S-{x},T(S-{x}))。將樹中所有點的路徑長度之和定義為隨機分割樹的復雜度|M(T)|=∑f(y,S,T),y∈S,則將點x刪除后,RRCF復雜度變化可表示為[8]:
ET(S)[|M(T)|]-ET(S-{x})[|M(T(S-{x}))|]=
f(y,S-{x},T′))+
(1)
式中,ET(S)[|M(T)|]表示樣本為i時隨機分割樹復雜度的期望;Pr(T)表示每棵樹的概率密度;T′表示基于樣本集S-{x}構建的樹T(S-{x})。
為有效避免異常屏蔽現象所引起的漏檢問題,在執行刪除操作時并不只是刪除點x,而是刪除一個含有x的集合C(x∈C?S),此時集合C的異常得分將歸因于C中所有點。基于此,使用共謀異常得分衡量點x的異常程度[8]:
CODISP(x,S,|S|)=
f(y,S-C,T(S-C)))]
(2)
式中,集合C稱為共謀集合。C中點的個數,即共謀集合的大小n將影響算法對異常屏蔽現象的識別能力。若n過大,異常點所在共謀集合中正常點的異常評分偏高,易將正常值誤判為異常值;若n過小,則單獨刪除某個含有異常點的共謀集合時,模型復雜度不會發生顯著變化,可能會將異常值誤判為正常值,從而造成漏檢。因此,合理設置共謀集合的大小可以有效提高算法的穩健性和準確率。
RRCF算法利用新增數據引起的整體模型復雜度變化衡量其異常程度,實現對實時數據流的異常探測;通過對模型數據重復性隨機分割獲取統計特性,確保異常值識別的準確度;使用數據替換的方式代替逐次建模,實現對模型的動態維護,可顯著提高處理效率,滿足海量數據的快速處理與分析需求。
基于RRCF的GNSS監測位移數據異常實時檢測算法主要分為模型構建和異常評估兩個階段[12],算法架構如圖1所示。

圖1 RRCF算法總體架構
若已有GNSS位移監測歷史數據,則可使用歷史數據對模型進行初始化,即熱啟動建模:
1)輸入RRCF異常值檢測算法參數:隨機分割樹中樣本個數n、隨機分割森林中樹的數量N、異常評分閾值T等;
2)初始化創建N棵空樹,結合訓練數據集進行模型構建,并計算每個訓練集樣本點的共謀異常得分,將共謀異常得分與閾值T進行比較,判斷該點是否為異常點。
若無可作為訓練數據集的歷史數據,則依靠實時數據流完成模型的初始化構建,即冷啟動建模。與熱啟動不同,冷啟動將實時數據流中待檢測的點插入到N棵空樹中,并計算共謀異常得分。若樹中樣本個數達到預設的樣本個數n,隨機選擇一棵樹,用待測點替換樹中最舊的點,否則直接將待測點插入到樹中。
1)根據實時GNSS位移監測數據流的RRCF異常值檢測模型,用實時輸入數據流中每個待測樣本點替換樹中最舊的點;
2)使用RRCF異常值檢測模型計算待測樣本點的共謀異常評分,并與異常評分閾值進行比較,若大于閾值,則認為是異常數據,否則認為是正常數據。
3.1.1 實驗數據
為驗證算法的效果,采用GNSS位移分析中常用的坐標時間序列模型[13],模擬生成不含粗差的原始數據,函數表達式如下:
bmcos(2πfmti)]+rti
(3)
式中,i為坐標歷元時刻標識;x(ti)為GNSS測站某分量ti時刻的坐標;b為截距;v0為線性速度;m0為諧波個數;am和bm表示頻率為fm時周期項的振幅;rti為隨機噪聲,即rti∈N(0,σ2)。上述先驗理論模型并不足以模擬受多種因素影響的GNSS監測位移序列的變化特征,因此分別基于理想先驗模型和非理想先驗模型生成2組不含粗差的原始數據。其中,理想先驗模型生成的模擬數據參數為:b=5,v0=2,m0=2,a1=5,b1=5,a2=3,b2=3,σ=3,f1和f2分別表示周期為1 a 和半年對應的頻率;非理想先驗模型生成的模擬數據中,另外加入高頻周期項,其參數為a3=6,對應頻率f3=7。
在模擬粗差數據時,首先采用標準差為3σ(σ為隨機噪聲標準差)的正態分布模擬得到一組隨機誤差序列,然后將大于3σ的誤差序列隨機加入到不含粗差的原始數據中,最終得到受粗差污染的GNSS位移時間序列[3]。實驗數據由MATLAB平臺生成,歷元間隔4 h,共2 223個歷元。模擬的粗差總數為140,占總觀測歷元的6.30%,實驗數據如圖2所示。

圖2 實驗數據
3.1.2 評價指標
異常值檢測效果通過準確率A、精確率P、召回率R以及綜合評價指標F1進行評價:
A=(nTP+nTN)/(nTP+nTN+nFP+nFN)
(4)
P=nTP/(nTP+nFP)
(5)
R=nTP/(nTP+nFN)
(6)
F1=2P·R/(P+R)
(7)
式中,nTP、nTN、nFP、nFN分別表示正常點檢測為正常、異常點檢測為異常、異常點檢測為正常、正常點檢測為異常的樣本點個數。準確率可衡量算法預測正確的結果占總樣本的百分比;精確率可描述所有被預測為正常的樣本中實際為正常樣本的概率,可反映算法對異常值的誤判率;召回率可描述實際為正常的樣本被預測為正常樣本的概率,可反映算法對正常值的誤判率;綜合評價指標F1為綜合精確率和召回率對算法效果的整體評價。
3.1.3 結果分析
利用小波-3σ方法[3]和RRCF實時異常檢測方法對2組實驗數據進行處理與分析。其中,小波-3σ方法依據經驗選取變形分析中常用的sym7小波基函數,通過篩選確定分解層次為6;RRCF方法采用模擬實時數據處理的策略,使用前300組數據進行建模,剩余數據遞次送入模型進行檢驗,判斷是否存在異常。
理想先驗模型實驗數據處理結果如圖3所示,由圖可知,小波-3σ方法可以探測出GNSS變形序列中的大部分粗差,RRCF實時粗差探測效果與事后小波-3σ方法相當。非理想先驗模型實驗數據處理結果如圖4所示,由圖可知,相較于理想先驗模型實驗數據處理結果,2種方法均可探測出大部分粗差,但漏檢個數明顯增多。同時,小波-3σ方法略優于RRCF方法。

圖3 理想先驗模型實驗數據粗差探測結果

圖4 非理想先驗模型實驗數據粗差探測結果
表1為2種方法粗差探測結果統計,由表可知,在理想先驗模型仿真數據條件下,相較于小波-3σ方法,RRCF方法的檢測異常數和漏檢個數較多,準確率和精確率略低,但召回率和總體評價指標略高。由召回率可知,RRCF方法在理想先驗模型中將正常值誤判為異常值的概率較低。在非理想先驗模型仿真數據條件下,相較于小波-3σ方法,RRCF方法的檢測異常數和漏檢個數較多,準確率、精確率、召回率和總體評價指標均略低。與小波-3σ方法相比,RRCF方法更適用于緩變的GNSS位移監測序列的異常值檢測。小波-3σ方法在兩種條件下的總體評價指標F1分別為0.967 3和0.994 0,RRCF方法分別為0.975 7和0.976 8。對理想先驗模型仿真數據而言,RRCF方法略優于小波-3σ方法;而在非理想先驗模型仿真數據條件下,RRCF方法略低于小波-3σ方法。

表1 粗差探測結果統計
RRCF算法作為一種實時異常檢測方法,僅依據當前歷元收集的數據對下個歷元數據進行異常判斷,對異常點識別的準確率、精確率與召回率接近事后方法,可以滿足海量實時數據流的異常探測需求。
黑臺滑坡位于甘肅省永靖縣鹽鍋峽鎮,使用該滑坡05GP01測站2019-11-01~2021-10-25實測GNSS位移監測數據,樣本個數為9 589個。原始數據包含水平位移、垂直位移與合位移,本文將以垂直位移處理結果為例進行詳細分析。
圖5(a)為垂直方向GNSS位移監測原始數據,紅色方框中數據點為明顯的離群異常點。對數據進行劃分,取前1 000組數據為訓練集對隨機分割森林進行初始化,剩余數據作為測試集。數據測試結果表明,隨機分割樹的最大容量取256、512和1 024時影響不大,本文設置為256;共謀集合大小依據經驗取為5,既不會將粗差值湮沒,也可避免因捕捉噪聲導致誤警。對大量實測GNSS位移監測數據進行處理,結果表明,異常評分閾值設置為50時符合大多數GNSS累積位移監測數據粗差檢測的實際需要。

圖5 RRCF方法異常評分
圖5(b)為實時數據的異常得分,將異常得分與閾值進行比較,可篩選出異常值。最終在9 589個樣本中檢測出60個異常值,GNSS位移數據異常值識別效果與刪除異常值之后的正常數據如圖6(a)和6(b)所示。由圖可知,GNSS位移監測數據發生異常突變時,RRCF方法可以很好地識別異常值,其檢測結果與實際異常值情況吻合且誤判率較低。上述結果表明,基于RRCF的GNSS位移監測數據異常實時檢測算法具有可行性,且對實時動態變化的數據流具有很好的適應性。

圖6 原始GNSS位移監測數據及RRCF方法粗差探測結果
為進一步驗證RRCF實時算法粗差探測的準確性,使用小波-3σ方法對該數據進行事后異常檢測。圖7(a)為異常點探測情況,將異常點刪除后的正常數據如圖7(b)所示。小波-3σ方法在9 589個樣本中檢測出77個異常值,由圖7可知,該方法可識別并剔除大部分異常點。

圖7 原始GNSS位移監測數據及小波-3σ方法粗差探測結果
小波-3σ方法利用小波分解出的高頻系數計算得到觀測數據的中誤差估值σ[3],容易受到小波分解的趨勢項與實際變形特征的符合程度、粗差大小及分布的影響和干擾,導致求取的中誤差估值不準確,進而影響粗差探測效果。因此,小波-3σ方法容易因小波基和分解層數設置不合理而出現誤判或漏判。實時RRCF方法根據數據插入或刪除后模型的復雜度變化來衡量數據異常程度,可根據數據實時變化特征進行靈活調整,且評分標準統一,不易受到數據特征的影響,可適用于GNSS位移數據的實時監測。
在GNSS位移粗差實時檢測過程中,若水平方向位移、垂直方向位移與合位移均發生異常,則應注意該測站其他傳感器或該滑坡體的其他測站是否也同時出現警報,綜合判斷是否發布地質災害預警。
本文提出基于RRCF的GNSS位移監測數據異常實時檢測方法,同時采用仿真數據與實測數據對算法進行驗證與分析,得到以下結論:RRCF方法在理想先驗模型條件下略優于小波-3σ方法,準確率和精確率略低,召回率略高,異常檢測的總體評價指標F1為0.975 7;而在非理想先驗模型仿真數據條件下,RRCF方法略差于小波-3σ方法,準確率、精確率、召回率均略低,總體評價指標為0.976 8。針對數據異常實時檢測場景,采用地質災害位移監測數據進行驗證,結果表明GNSS位移監測數據發生異常突變時,RRCF方法可以很好地識別異常值,其檢測結果與實際異常值情況吻合且誤判率較低,具有較高的檢測準確率和可用性,可為GNSS位移實時監測預警提供技術支撐。