應用時空復合趨勢面Kriging 方法插值地表沉降

黃丙湖，李欣芮，范芷睿，潘海燕，廖一蘭

（1. 中國石油大學（華東） 海洋與空間信息學院，山東 青島 266580；2. 新疆大學建筑工程學院，新疆 烏魯木齊 830046；3. 中國科學院地理科學與資源研究所資源與環境國家重點實驗室，北京 100101）

地鐵施工中土質、荷載等因素將導致地面沉降，易引發坍塌、開裂等事故，因此必須建立相應的觀測體系進行風險評估，以便采取相應措施防止事故災害的發生。目前的常用方法為在地鐵施工線上設立地表沉降監測點，并定期觀測監測點沉降數據，但所得的采樣數據較離散，無法直觀體現整體沉降趨勢，必須對采樣數據進行插值[1-2]。空間插值考慮了沉降監測點在空間上的相關性，能對整體沉降區域的監測點進行預測。常用的空間插值方法包括幾何方法、統計方法、函數方法和綜合方法[3]，其中以Kriging 為代表的地學統計方法可定量描述空間隨機分布變量的相關性，已應用于測繪、地質學、GIS 等學科，具有高度可行性和有效性[4-8]。陳雷[9]等利用Kriging 插值方法對地面沉降數據進行內插，但選取樣本量小，只能反映局部趨勢；劉志平[10]等利用Kriging插值方法研究了邊坡變形速率，但其將變形速率演化為與時間相關的函數，弱化了時空關聯性。在對地表沉降量等時空變量進行插值前，除分析變量空間相關性外，還應考慮時間相關性，王建民[11]等利用時空Kriging插值方法，將時間域擴展到空間域，建立了時空半變異函數，該方法計算量大，選擇變異函數模型時需人工干預，會造成一定誤差；Kyriakidis P C[12]等利用時空Kriging 插值方法研究區域變化量的時空分布特征，但未與傳統地統計插值方法進行比較，無法體現其優勢。時空復合趨勢面（CSTM）兼顧了沉降點的時間相關性和空間相關性，利用已知沉降點監測值進行趨勢面擬合得到區域內每個時空變量的全局平均趨勢。因此，本文在傳統地統計插值模型的基礎上引入了CSTM，使模型能更準確地反映時空變量特征，從而提高模型預測精度；并與傳統方法進行了精度對比，驗證了時空插值的優勢。

1 基于CSTM的Kriging插值方法

1.1 CSTM

CSTM 首先通過區域內已知時空變量的屬性值進行趨勢面擬合，再利用指數平滑法得到每個時空變量的屬性值。該方法將區域內每個屬性值分為趨勢值和殘差值兩個部分，以此反映時空變量的變異情況。Lee[13]等提出CSTM 是基于時空隨機場等理論定義的。時空隨機場[14]即z(p)=z(s't)，包括時間尺度和空間尺度兩個維度，s為二維空間位置（本文中為沉降監測點空間位置），t為時間（本文中為沉降監測點觀測時間），將Y(p)=log[z(p)]作為z(p)的對數變換。全局平均趨勢mY(p)是由Y(p) 構建的確定性函數，反映Y(p)的全局平均趨勢。CSTM是基于復合時間/空間指標獲得mY(p)的有效方法[13]。傳統理論獲得mY(p)（即mY(s't)）的方法是將其分為時間和空間兩部分，而CSTM的優勢在于兼顧了時間分量與空間分量[15]。

對于一個周期性變量，CSTM 可利用基于指數平滑函數的權重來計算擬合函數，以獲得每個時空坐標平滑的平均趨勢值，在任意時空坐標pj=(sj'tj)下，平均趨勢值的計算公式為[13]：

式中，Y(sj'tj)為沉降量在可測量監測點pj處的對數變換；wi為權重；d(si'sj)為沉降監測點pi、pj之間的歐氏距離，；as、at分別為指數平滑函數的空間范圍和時間范圍。

1.2 基于CSTM的Kriging插值

CSTM 將區域內每個屬性值分成趨勢值和殘差值，基于CSTM的Kriging插值是對殘差進行插值，再與趨勢值求和，最后得到時空變量的預測值。

1）數據空間分布的判斷。使用Kriging 插值方法的前提是待插值數據滿足正態分布且具有空間自相關性，因此可利用Moran’s I等空間自相關性統計量指標判定待插值數據的空間分布特征，滿足該前提才能使用該方法。

2）時空尺度的確定。式（2）中as的選取原則為小于等于兩個測量沉降監測點之間的最大歐氏距離，at的設置原則為小于等于沉降監測點多次觀測的最大時間差。根據式（2）計算得到指數平滑函數的權重。

3）mY(p)的計算。沉降過程中根據已知監測點沉降量來構建CSTM，利用式（1）計算得到每個時空監測點的mY(p)。

4）Kriging插值。Y(p)由確定性變換X(p)=Y(p)-mY(p) 產生，由局部穩定平均值mx(p)=E[X(p)]=E[Y(p)-mY(p)]獲得固有平穩時空隨機場，其中E[.]表示期望算子[13]，X(p) 為殘差。由于沉降監測點p=(s't)和p'=(s''t')之間空間滯后γ=‖s-s'‖ 和時間滯后τ=|t-t' |的影響產生誤差，二者之間的總體誤差用協方差計算，能衡量一對沉降監測點間的時空相關性，即

Kriging 插值方法考慮了樣本之間的空間相關性，通過最小化均方誤差權重來獲得未監測沉降點殘差的最佳估計值。X(p)在點k處的估計值⌒是測量值xd的線性組合，計算公式為：

式中，λ為采用傳統Kriging 插值估算權重的列向量，沉降點pk、pd的時空復合距離越近，權重越大；mx(pk)、mx(pd)為X(p)在預測點pk、pd的平均趨勢。

式中，cx(pk'pd)為估計點與數據點之間X(p)的協方差行向量；cx(pd'pd)為數據點之間X(p)的協方差矩陣[16]。

5）監測量插值結果的估計。將區域內每個沉降監測點的mY pj=(sj'tj) 與Kriging 插值得到的殘差對應求和，最后進行指數變換得到每個監測點的沉降量。

2 實驗與結果分析

2.1 數據來源

本文的實驗數據來源于烏魯木齊市軌道交通3 號線一期工程05標明園站—友好站區間，該工程處于施工階段，以2 號豎井和橫通道開挖。考慮到地層自穩能力差，受開挖、支護工程和地下水等因素影響，地表會出現不同程度沉降，因此豎井施工過程中針對可能沉降變化較大的區域布設41個監測點，構成了6條監測線（圖1），在井口頂部圈梁短邊周圍布設6個監測點，其余監測點分別布設在豎井橫通道附近和隧道線路附近。采用徠卡電子水準儀（±0.3 mm/km），按照二等水準測量方法，對地表監測點進行觀測，經平差合格后，獲得監測點處的高程；再統計觀測期間監測點的高程變化情況，獲得沉降變化量。實驗獲取了2017-06-24—2017-09-08 共15 期沉降數據，監測頻率為5 d/次，監測過程中個別監測點數據缺失，導致總數據量小于n=41×15。

圖1 地鐵地表沉降監測點布設圖

2.2 插值過程

首先，對獲取的數據進行空間自相關性檢驗，利用Moran’s I 進行判定，p值小于0.05 且z 得分超過1.65，認為通過95%置信度檢驗，結果見表1，可以看出，4個時段的數據未通過顯著性檢驗，且它們的Moran’s I 均接近于0，空間自相關性很弱，將其剔除。本文分別利用Kriging插值和基于CSTM的Kriging插值方法對研究區內符合條件的11個時段的監測點沉降量進行插值。實驗中選取每期監測數據利用Kriging方法進行獨立插值發現，預測沉降監測點與半徑為2 m 區域內的已知沉降點空間相關性最強，因此確定CSTM基于指數平滑函數權重的空間尺度為2 m，時間尺度選取待插值數據的前后1 個時間段。值得注意的是，由于剔除了7 月30 日—8 月4 日的數據，導致該時段周期變為10 d，即7月30日—8月9日，本文將這種時段相距天數不相等的間隔稱為非等時間間隔，反之稱為等時間間隔。然后，計算區域內每個沉降監測點的mY(p)。最后，對每個地表沉降監測點的殘差進行Kriging插值。

表1 實驗數據空間自相關性檢驗

兩種方法的插值結果見圖2，可以看出，圈梁、豎井橫通道附近的沉降程度相對明顯，因此施工過程中需加強監測這兩個部分；基于CSTM的Kriging插值方法考慮了時空相關性，7月15—25日圈梁附近地表緩慢下降，而7月25日—8月14日則開始抬升，后期逐步趨于穩定狀態；豎井和橫通道附近的點在不同時段內沉降波動較大，這可能是由施工推進造成的。

2.3 精度評價

實驗采用留一交叉驗證的方法，即對于N個樣本，每次選擇N-1 個樣本訓練數據，一個樣本驗證模型預測的好壞。該方法的優勢在于，每個回合中幾乎所有的沉降數據皆用于訓練模型，因此最接近原始沉降數據的分布，評估結果比較可靠。首先假定某一監測點沉降量為未知，然后利用插值方法計算該監測點沉降量，最后計算該監測點實際觀測值與插值得到的沉降值之間的誤差。這種誤差主要是由插值算法的準確性和監測點的分布情況造成的[17]。實驗選取均方根誤差（RMSE）和R2作為精度驗證指標（表2），RMSE 可反映監測點數據插值結果的離散程度，計算公式為：

表2 兩種方法的RMSE和R 2 比較

式中，vi為第i個監測點的實際沉降量；v'i為第i個監測點的差值估算值；n為參與交叉驗證的監測點數目。

R2可體現模型的擬合效果，由殘差平方和總體平方組成，其值越接近1，表明線性相關程度越高，計算公式為：

由表2 可知，空間插值和時空插值都表現出較強的線性相關性，本文方法比傳統方法的RMSE 降低約35%，說明本文方法的總體離散程度明顯小于傳統方法；傳統方法的R2大部分在0.5以下，而本文方法的擬合程度更好；綜合來看，本文方法提高了模型精度，更適用于地表沉降監測。

然而，CSTM 未給出標準的時間尺度確定方法，因此本文在前后一個以非等間隔時間段為指數平滑的時間尺度實驗之后，又將時間尺度設為等間隔時間段進行探索性實驗，結果見表3，可以看出，采用等間隔時間段的沉降數據構建CSTM 得到的精度更好。插值結果見圖3。

表3 不同時間尺度的本文方法的精度評價

圖3 同時期非等間隔與等間隔本文方法插值結果比較

CSTM 同時考慮了沉降監測點的時間相關性與空間相關性，因此選取研究區域內均勻時段進行預測將得到更好的效果，由于本次實驗數據量的局限性，等間隔沉降量差值的時間序列較短，未來可選取更長時間序列進行實驗。

3 結 語

本文將CSTM 引入Kriging 插值，同時考慮了空間相關性和時間相關性，使該插值方法能更真實地反映時空變量變化趨勢，提高了地表沉降監測精度。實驗表明，與傳統方法相比，本文方法的RMSE 降低約35%，具有更高的擬合精度；從插值結果來看，能更精細地反映沉降量的時空差異，以實現對地鐵地表沉降的有效監測。受限于數據觀測天數較少，本文具有一定局限性，未來可考慮對更長時間的沉降監測進行實驗，同時可引入周邊建筑物、土質、施工方式、人口密度、車流量等輔助變量來優化插值模型，實現長時間、大范圍的地表沉降趨勢分析。除此之外，本文方法還可與三維激光掃描數據的點云數據相結合，對格網化后的點云數據進行插值，實現更高效快捷的整體監測分析。