崔芳姿,吳 斌,史學磊,范振林,任 濤
(1.中國地質調查局水文地質環境地質調查中心,河北 保定 071051)
對于地質災害防治來說,監測并不是目的,通過大量監測數據預測災害變形趨勢乃至預報災害發生的時間,提前做出決策,做到防患于未然才是監測的最終目標。換句話說,即使是建立了先進有效的監測系統,取得了豐富的監測數據,而不具備對數據的綜合分析和對災害的預測預報能力,災難同樣可能發生[1-3]。因此,如何結合災體分析監測數據,預測災害變形趨勢,預報其失穩時間并提出合理的對策措施,是監測預警的核心內容[4-6]。
由于地質災害的復雜性,其預測預報技術是一項世界性的技術難題。盡管各種各樣的預報理論、模型不斷涌現,但應用成功的實例極為少見。預報成功的例證基本上是在連續監測,并對各種地質現象進行綜合分析的基礎上做出的,其成功率取決于技術人員對地質災害的認識水平、經驗和態度,具有很大的不確定性。另外,地質災害監測數據存在誤差以及由數據傳輸或設備問題導致的數據缺失,而國內外大多地質災害監測預警系統沒有對監測數據進行預處理[7-9]。因此,研究和提高地質災害監測數據預處理、監測數據綜合分析(特別是與致災地質作用的結合分析)和預測預報能力是十分必要的,也是極為緊迫的[10-11]。
本文提出了一套地質災害監測數據綜合處理與分析技術,并設計了地質災害數據信息化管理系統?;跀祿焯峁┑慕y一數據模型和數據服務,結合巫山縣地質災害監測預警示范站地質災害業務信息和地質災害監測數據,本文建立了一套統一的數據標準規范,實現了災害點、監測點、監測設備管理的數據綜合處理與分析服務,并通過該服務建立和完善了監測數據預處理方法、監測數據綜合分析方法和預警預報方法,進一步提高了地質災害監測預警信息的綜合處理能力,為地質監測管理決策人員提供輔助支持,保障地質災害監測預警工作順利進行。
數據綜合處理與分析服務構建于統一的地質災害監測數據和業務應用服務之上,基于統一的地質災害監測數據,首先將數據傳輸到數據庫,并對數據和統計量進行數據質量分析、相關性分析、異常值剔除、數據等距化、數據平滑等一系列處理;再繪制位移矢量圖、位移/速率-時間曲線,并利用切線角法進行預警預報,如圖1所示。數據綜合處理與分析服務能進一步提高地質災害監測預警信息的綜合處理能力,為地質監測管理決策人員提供輔助支持。

圖1 數據綜合處理與分析服務流程圖
基于數據綜合處理與分析服務,本文設計了一套地質災害數據信息化管理系統,如圖2所示,采用B/S技術架構,支持IE、360、Google Chrome等瀏覽器,總體包括數據采集層、地質災害監測數據庫、基礎管理、地質災害監測數據綜合處理與分析服務以及數據交換等部分,還包括安全保障體系、維護支撐體系、數據標準規范體系等外圍支持。系統基于巫山縣地質災害監測預警示范站的地質災害信息,通過采集巫山縣地質災害監測預警示范站地質災害業務信息和地質災害監測數據,構建統一的地質災害數據庫;利用數據庫提供的統一數據模型和數據服務,構建地質災害數據信息化管理系統,為地質災害監測預警業務提供應用服務;通過數據接口和數據交換,為其他政務系統提供各類信息服務。

圖2 系統技術架構圖
監測數據預處理包括異常數據剔除、差補與平滑等內容。
2.1.1 異常數據剔除
在科學研究實驗中,經常會出現因儀器故障或讀數錯誤而出現的誤差。對于大量數據可疑值的判別,常用的統計判別方法為拉依達檢驗(3σ準則)。該準則是以99.73%的置信度,根據隨機誤差正態分布理論判斷是否存在可疑值。若|xi-x|>3σ,則應舍去xi,否則應保留xi;可用max(abs(xi-mean(xi)))≥3std(xi)條件是否滿足來實現,其中abs()是絕對值函數,max()是最大值函數,std()是標準差函數。
2.1.2 數據等距化(差補)
由于儀器數據有時存在缺失或兩端數據的數據間隔不同,不方便進行數據分析,需對數據進行等距化,將數據轉化為間隔相等的數據。
本文采用三次樣條插值法進行等距化處理,即對數據進行樣條函數曲線擬合得到一條曲線,再由擬合后的曲線計算得出其他任意時刻的值。設有兩個傳感器A和B,分別以不同的頻率對目標進行采樣測量,各傳感器可均勻采樣,也可非均勻采樣。每個傳感器在采樣時刻ti的測量值記為(ti,yi),其中yi為向量,這樣每個傳感器均可得到一組測量值。由于傳感器的采樣周期不同,各傳感器獲得數據的時間值不盡相同,若直接進行融合,則可能會因時間差而使融合結果失去意義,反而不如單個傳感器的融合精度高,因此在進行數據融合前,必須將不同時刻的測量數據對準到同一時間點上。
設傳感器A在某一時間段[a,b]內對目標進行了n+1次測量,將整個時間區間按采樣時刻劃分為a=t0<t1<…<tn=b,給定的時刻點ti對應的觀測值為f(ti)=yi,i=0,1,…,n,構造一個三次樣條插值函數s(x),使其滿足下列條件:①s(ti)=yi,i=0,1,…,n;②s(x)在每個小區間[ti,ti+1]上是一個三次多項式,i=0,1,…,n-1;③s(x)在[a,b]上具有二階連續導數。
基于最小二乘法的樣條函數擬合是在樣條函數空間內找出f(t)關于范數的最佳逼近,即找到s*(t),使
三次樣條插值函數的構造過程為:記mi=s'(ti),i=0,1,2,…,n在每個小區間[ti,ti+1],i=0,1,2,…,n-1上,利用Hermite插值公式寫出三次樣條插值函數的計算公式,即

利用s"(ti-)=s"(ti+),i=0,1,2,…,n-1,并附加邊界條件s"(t0)=s"(tn)=0,可得方程組:

式中,a0=1;

方程組系數矩陣為三角矩陣,其行列式不為0,所有方程組的解存在且唯一, 對方程組求解,即可得到遞推公式:

利用公式求出ai、bi,i=1,2,…,n,令mn+1=0,求出mn、mn-1、…、m0,將參數ti、yi、mi,i=1,2,…,n代入s(x)中即可求得三次樣條插值函數。
經過樣條插值擬合,可以得到一條平滑曲線,從而由其求得傳感器在任意時刻的值;再與其他傳感器進行時間對準,根據采樣時刻,從該條曲線中取出相應時刻的測量值。
2.1.3 數據平滑
儀器觀測數據存在隨機誤差(波動),數據趨勢分析時需摒棄誤差,找到數據真正的變化趨勢,因此需進行數據平滑處理。對于數據長度較小的數據,采用整體平滑處理,即三階一般多項式擬合平滑方法;對于長期監測數據,則采用窗口的移動方法,分時段進行平滑處理。
對于給定的一組數據(xi,yi),i=0,1,…,m要求在函數類ψ={ψ0,ψ1,…,ψn}中找一個函數y=S*(x),使誤差平方和最小,即

式中,S(x)=a0ψ0+a1ψ1+…+anψn,n<m。
這就是一般的最小二乘逼近,即曲線擬合的最小二乘法。為了使問題的提法更有一般性,通常把最小二乘中的考慮為加權平方和,即

這里w(x)≥0為[a,b]上的權函數,表示不同點(xi,f(xi))處的數據比重。利用最小二乘法求擬合曲線的問題,就是在S(x)中求函數y=S*(xi),使式(6)取得最小。將其轉化為求多元函數I(a0,a1,…,an)=問題。根據多元函數極值的必要條件,則有:

若記

則可改寫為:

將其寫成矩陣Ga=d,其中a=(a0,a1,…,an)T,d=(d0,d1,…,dn)T,則有:

由于ψ0,ψ1,…,ψn線性無關,因此|G|≠0,式(6)存在唯一解,ak=ak*,k=0,1,…,n,從而得到函數f(x)的最小二乘解為:

可以證明這樣得到的S*(x),對任何形如式(1)的S(x),均有:

2.2.1 相關性分析
相關性分析可用來驗證兩個變量間的線性關系,通過相關系數可反映兩個變量是否呈線性關系、線性關系強弱以及是正相關還是負相關。最簡單最常用的相關系數計算方法為積差法,即利用兩個變量的協方差與二者標準差的乘積之比來表示,計算公式為:

相關性系數介于-l~1之間,若r接近0,則兩個變量沒有線性關系;若r接近-l或1,則說明兩個變量線性相關性很強;r為正,表示正相關,r為負,表示負相關。
2.2.2 變形趨勢分析
數據綜合處理與分析服務不僅針對原始監測數據,還針對原始數據的統計量。各個統計量的計算規則為:累計變形量,即當前值-初始值;平均變形速率,即(當前值-初始值)/月;上半年變形增量,即當前值-上半年初始值;上半年變形速率,即上半年變形增量/6;上半年數據中誤差,即sqrt(((dx12+dx22…)-(dx1+dx2)2/n)/(n-1));變形方向,即dx累計變形量和dy累計變形量計算。
2.3.1 變形趨勢預測方法
位移矢量圖(平面圖)反映了變形監測點隨時間在空間上的變化情況,主要表示變形量和變形方向,即二維空間上的x方向累計變形和y方向累計變形,能直觀顯示滑坡整體的變形態勢。
2.3.2 切線角法
對于某個滑坡來說,等速變形階段的位移速率v為一個恒定值,則可通過累計位移S/v的辦法將S-t曲線的縱坐標變換為與橫坐標相同的時間量綱。定義S(i)為某一單位時間段(一般采用一個監測周期,如1 天、1 周等)內的斜坡累計位移量,v為等速變形階段的位移速率,T(i)為變換后與時間相同量綱的縱坐標值。
由T-t曲線(圖3)可得到改進的切線角的表達式為:

圖3 滑坡位移/速率—時間曲線
本文基于數據綜合處理與分析技術設計了一套地質災害數據信息化管理系統,主要包括基礎地理信息模塊、地災監測信息“一張圖”模塊、災害信息管理模塊、地災專業監測模塊、數據綜合處理與分析模塊和系統維護管理模塊,如圖4所示。

圖4 系統功能結構圖
通過綜合處理與分析服務可以判斷監測數據是否可靠,調整監測重心,預測邊坡的失穩時間,掌握宏觀變形特征,判斷所處變形階段和穩定性現狀,進而研判未來一段時間的地質災害變形發展趨勢。
1)監測數據預處理方法應用。將斜坡變形演化階段應變(位移)—時間曲線的三階段過程視為滑坡預測預報的基本標準:第I階段為初始變形階段,又稱減速變形階段,巖土體變形以減速發展,變形曲線斜率逐漸減??;第Ⅱ階段為等速變形階段,又稱穩定變形階段,巖土體變形大致以等速發展,變形曲線近似一條傾斜直線,應變速率大體不變;第Ⅲ階段為加速變形階段,巖土體變形速率迅速增加,至巖土體破壞。
通過異常數據剔除、數據等距化與數據平滑等預處理,提高了地質災害監測數據的可用性,形成了多參數歷史曲線圖。當坡體已進入加速蠕變階段和臨滑階段時,可根據多參數歷史曲線圖上位移隨時間增長速率的變化來預測邊坡的失穩時間。
2)監測數據綜合分析方法應用。通過變形趨勢分析可統計監測運行以來各監測點的累計水平位移量、累計垂直位移量、累計水平變形速率、累計垂直變形速率、上半年度水平位移增量、上半年度垂直位移增量、上半年度水平變形速率、上半年度垂直變形速率、上半年度水平位移中誤差、上半年度垂直位移中誤差和變形方向。
利用滑坡變形量和變形速度能調整滑坡監測頻次。對于多年監測緩慢變形或基本不變形的災害點應降低專業監測頻次甚至改為群測群防監測,以減少監測工作投入和數據處理與分析的工作量。通過變形方向可以判斷各監測點的變形方向與主滑方向是否一致。將中誤差與技術規范中要求的監測數據誤差進行對比分析,實測誤差是否在容許范圍之內,以判斷數據的可靠性。
3)預警預報方法應用。位移矢量圖可反映各監測點的變形量、變形方向、變形發生部位,各監測點上變形差異,屬于整體變形還是局部變形以及宏觀變形特征等。通過曲線可判斷各監測點變形的主要發生時間;綜合變形特征與宏觀巡查分析,可判斷所處變形階段和穩定性現狀,從而研判未來一段時間的地質災害變形發展趨勢。
當切線角α≤45°時,斜坡變形處于等速變形階段,進行藍色預警;當切線角45°<α<80°時,斜坡變形進入初加速變形階段,進行黃色預警;當切線角80°≤α<85°時,斜坡變形進入中加速變形階段,進行橙色預警;當切線角α≥85°時,斜坡變形進入加加速變形(臨滑)階段,進行紅色預警;當切線角α≈89°時,滑坡進入臨滑狀態,應發布臨滑警報,如表1所示。

表1 滑坡變形預警級別劃分
地質災害監測數據綜合處理與分析服務是地質災害數據信息化管理系統的關鍵,通過異常值剔除、數據等距化、平滑、繪制位移—時間曲線等處理,能大大提高地質災害監測數據的可用性,提高地質災害監測預警信息的綜合處理能力。通過數據綜合處理與分析服務,建立和完善監測數據預處理方法、監測數據綜合分析方法以及預警預報方法,結合變形特征判斷所處的變形階段和穩定性現狀,進而研判未來一段時間的地質災害變形發展趨勢。利用數據預處理方法、綜合分析方法和預警模型解決了地質災害預警工作中的數據缺失、數據質量差、預警模型缺失等問題,為地質災害預警工作和防災減災工作提供了技術支撐,對于防治決策服務成功實施和地質調查工作順利進行具有重要意義。