基于Kriging方法的負偏態分布數據的空間插值研究

黃 靜，王愛倩，翟世龍

（1.新疆農業大學 水利與土木工程學院，新疆 烏魯木齊 830052；2.新疆克孜爾水庫管理局，新疆 庫車 842000；3.新疆地震局，新疆 烏魯木齊830011）

Kriging插值方法又稱空間局部插值法，是以變異函數理論和結構分析為基礎，在有限區域內對區域化變量進行無偏最優估計的一種方法［1-2］。該理論和方法建立在二階平穩或內蘊假設的基礎上，所以要求進行分析計算的數據服從正態分布，但在實際應用中，數據常具有異常值，高偏度以及非正態分布性質的特點，這對于變異函數擬合及插值的穩健性有著極大的影響［3］。目前在對實際數據進行穩健處理方面做了很多研究，最為有效的方法是將非正態的原始數據轉化為正態或近似正態分布的變換［4］。變換的方法有直接變換，如取倒數、取平方根、取對數，以及變換能夠較強的冪變換和Johnson變換。但實際中，正態變換并不是對所有的數據都奏效，比如存在趨勢的拖尾負偏態數據，無法進行正態轉換時，就要先對數據進行趨勢分析，采用全局多項式擬合配合Kriging插值方法來達到穩健性插值的目的。

1 數據基本情況

本文選取的樣本數據是某水庫庫底高程數據，庫區范圍：南北方向約長13km，東西方向約長8km，樣本數共計466個，采樣間隔不等，間隔最小約150m，最大約1 000m左右。

2 數據基本統計分析

對466個樣本數據進行直方圖統計分析，如圖1所示，在直方圖中現實樣本數據是存在拖尾的負偏態數據，峰值為2.55，偏度為-0.74，對其進行KS正態檢驗，P值小于0.01，說明該樣本數據不符合正態分布，見圖2。

3 數據趨勢分析和全局趨勢擬合

對該樣本數據進行正態變換，仍無法通過K-S正態檢驗。對數據進行分析，采用趨勢分析來判斷樣本數據在特定方向是否存在趨勢，若存在趨勢，先選擇一個擬合效果最好的多項式對散點進行內插，然后再對殘差隨機短程變異進行Kriging插值建立殘差模型，見圖3。

在圖3中，X軸代表實際地理位置的南北方向，Y軸代表實際地理位置的東西方向，因此從該圖中可以看出，該樣本數據南北方向呈現南高北低趨勢，在東西方向呈現兩頭高中間低的趨勢，因此，采用二次多項式對該數據進行全局趨勢擬合，并得到殘差數據。

圖1 樣本高程數據頻率直方圖

圖2 樣本高程數據正態分布概率圖

圖3 樣本高程數據南北-東西方向投影趨勢

4 殘差數據的隨機短程變異函數模型參數估計

對殘差數據計算變異函數值，通過變異函數模型進行模擬，目的在于確定一個最佳的擬合模型，最終計算出變異函數模型的幾個主要參數：塊金值、基臺值和變程。要想獲得精度較高的Kriging插值結果，變異函數模型的參數估計至關重要，首先要確定離散數據計算變異函數值的步長、最大步長和步數，其次要根據不同方向上的變異函數值分布情況，判斷數據是否存在各項異性，將此影響加入到模型中。

最大步長的設定一般是研究區域沿某一方向的最大尺度［5］，最大步長等于步長×步長組數。本文采用平均最鄰近法確定步長大小，步長為203m，步長組數為30，最大步長是6 090m，約是最大方向長度的一半。

在圖4中，變異函數表面圖上每個柵格單元用顏色進行編碼，從圖4中可以看出，在西北-東南方向的變異性要比南西-北東方向增加的略快。此時，選用的變異函數模型為tetraspherical，變程為2.5 km，基臺值為17.682，不存在塊金效應。下面對數據進一步進行方向效應分析，見圖5。

圖4 未考慮方向效應的變異函數散點圖及表面

圖5 考慮方向效應的變異函數圖及表面

在圖5中，散點圖顯示的是在某一方向的變異函數分布情況，在變異函數表面圖中橢圓的長軸表示變異函數在該方向有最大變程，其值為5.2km，超過這個變程殘差數據空間不相關，長軸方向表示變異函數在此方向有最小變程，其值為1.9km。當變異函數值達到穩定水平時的值相同時，這就是函數模型的基臺值為18.507，塊金值為0.023 57。

5 Kriging插值結果及交叉驗證

得到變異函數擬合模型，并且該模型中考慮了趨勢面的剔除和方向效應，從而計算Kriging插值系數，對數據進行插值。

為了選擇最佳的變異函數模型來提高精度，分別采用不同的變異函數模型進行計算，并通過交叉驗證［6］結果進行比較，本文把最佳插值精度結果列出，選用的模型是Tetraspherical，見圖6。

圖6 Kriging插值交叉驗證結果

從圖6可以看出，選用Tetraspherical作為變異函數擬合模型進行Kriging插值，其交叉驗證預測精度結果是：平均預測誤差為-0.029 53，該值應接近于0；均方根標準誤差為1.001，其值應接近于1；均方根誤差為1.479，平均標準差為1.417，其值越小，說明精度越高，是選取最佳變異函數模型的依據。

6 結 論

本文對拖尾的負偏態數據進行深入分析，將數據剝離為趨勢項和殘差兩部分，分別采用二次多項式曲面擬合及Kriging插值的方法進行插值，在變異函數模型擬合時考慮方向效應，對模型參數計算進行優化，其插值結果采用交叉驗證進行精度評定，最終選定Tetraspherical作為最終的變異函數模型，其插值精度最高。

通過對Kriging插值結果進行交叉驗證發現，預測精度較低的點主要分布在測區的邊緣地帶，其原因有二，一是樣本數據選用的是某水庫庫底的高程數據，其地形變化在庫邊變化較快；二是在采樣時，庫邊數據采樣間隔大，數據量少。在對部分空間數據進行插值時，邊緣數據的處理至關重要，采用什么模型和方法來對空間數據進行插值來提高邊緣數據的插值精度，是以后研究的重點。

［1］曾懷恩，黃聲享.基于Kriging方法的空間數據插值研究［J］.測繪工程，2007，16（5）：5-8.

［2］張梁，林韜，汪慶鋒.使用克里金插值法進行CGCS2000到海南地方坐標系的轉換［J］.測繪與空間地理信息，2014，37（9）：219-222.

［3］李曉輝，袁峰，白曉宇，等.典型礦區非正態分布土壤元素數據的正態變換方法對比研究［J］.地理與地理信息科學，2010，26（6）：102-105.

［4］陳道貴，胡乃聯，李國清.區域化變量非正態分布的穩健性［J］.北京科技大學學報，2009，31（4）：412-417.

［5］吳學文，晏路明.普通Kriging法的參數設置及變異函數模型選擇方法［J］.地球信息科學，2007，9（3）：104-108.

［6］陳華鑫，姜藝，李碩，等.基于ArcGIS的交通量預測模型［J］.同濟大學學報，2010，38（8）：104-108.