宋敬茹， 胡 尹， 仝梓晨， 劉明月,2,3,4，高均海， 張永彬， 滿衛(wèi)東,2,3,4*

(1.華北理工大學(xué)礦業(yè)工程學(xué)院;2.唐山市資源與環(huán)境遙感重點實驗室;3.河北省礦區(qū)生態(tài)修復(fù)產(chǎn)業(yè)技術(shù)研究院;4.河北省礦業(yè)開發(fā)與安全技術(shù)重點實驗室;5.中煤科工生態(tài)環(huán)境科技有限公司唐山分公司：河北 唐山 063012)

降水是生態(tài)環(huán)境循環(huán)的必要環(huán)節(jié)和農(nóng)業(yè)用水的主要來源[1]。水資源管理、旱澇災(zāi)預(yù)防及生態(tài)治理等都離不開降水信息[2]?？臻g插值法是分析地理要素空間特征的強大工具[3]，主要有克里格插值法(Kriging)、反距離加權(quán)法(Inverse Distance Weight，IDW)和徑向基函數(shù)法(Radial Basis Functions，RBF)等[4]。

基于不同的插值方法對區(qū)域降水量進行分析，一直是當(dāng)前研究的熱點[5-6]。插值方法各自具有優(yōu)勢，插值精度存在很大的差異，一些學(xué)者根據(jù)研究領(lǐng)域特點采用了不同的研究方法。何紅艷等[7]根據(jù)研究區(qū)域的特點，采用混合插值法繪制降水圖。孟慶香等[8]將地統(tǒng)計學(xué)方法和普通插值方法進行對比研究，得到地統(tǒng)計學(xué)方法優(yōu)于普通插值方法。馬友平[9]利用遙感影像數(shù)據(jù)研究巖溶區(qū)石漠化，發(fā)現(xiàn)采用克里格插值法對數(shù)據(jù)的“缺失”區(qū)進行修復(fù)是值得借鑒的。朱珊珊[10]用地理信息系統(tǒng)的統(tǒng)計分析功能對流域年均降水量序列進行空間插值，認(rèn)為克里格法是最優(yōu)的插值。Sergio和Christopher Daly等[11-12]利用地統(tǒng)計學(xué)分析進行插值方法對比研究，獲得了高質(zhì)量的氣候圖。

如何依據(jù)已知降水?dāng)?shù)據(jù)和研究區(qū)域的特點選擇最優(yōu)的插值方法一直是需要探討的難點[13]。現(xiàn)有的研究表明，新疆地區(qū)就普通插值方法而言克里格插值法精度最高，能更好地反映新疆區(qū)域的降水空間分布特征[14]。進一步研究發(fā)現(xiàn)，運用地理加權(quán)回歸克里格法(geographically weighted regression Kriging，GWRK)建立插值變量與其影響因子之間的回歸模型，將克里格插值和回歸分析結(jié)合并引入與插值變量具有良好相關(guān)性的影響因子作為輔助變量，能夠有效提高插值方法的精度[15-16]。

該文采用克里格法(Kriging)、反距離加權(quán)法(IDW)、徑向基函數(shù)法(RBF)和地理加權(quán)回歸克里格法(GWRK)的對比來研究新疆的降水空間變化，為該地區(qū)的生態(tài)環(huán)境治理等工作提供最優(yōu)數(shù)據(jù)支持。

1 研究區(qū)與研究方法

1.1 研究區(qū)概況

新疆在中國西北部，面積約占全國總面積的1/6(圖1)。中部天山橫亙把新疆分為2部分，山脈盆地交替排列被喻“三山夾兩盆”[17]。新疆地處北溫帶，日照時間長，遠離海洋，深居內(nèi)陸又有高山阻隔，造成其日氣溫差大、降水量稀少[18]。新疆年均降水量僅有105 mm，南部降水量只有10 mm左右，差異較大。新疆特殊的地理環(huán)境使其具有獨特的氣候特征，常發(fā)生干旱災(zāi)害[19]，因而以新疆作為研究區(qū)域，揭示新疆區(qū)域降水分布特點與規(guī)律具有重要的意義。

1.2 數(shù)據(jù)來源和站點分布

新疆地區(qū)降水?dāng)?shù)據(jù)來自國家氣象數(shù)據(jù)網(wǎng)(http://data.cma.cn/)的66個國家基礎(chǔ)氣象站，另外收集各氣象站的經(jīng)緯度和高程等信息。各氣象站點分散均勻，具備一定的獨立性、實用性和真實性，并選取20%的氣象站點作為驗證點站點，氣象站點分布如圖1所示。同時還結(jié)合了地形因子，獲取了來自地理空間數(shù)據(jù)云的16幅包含新疆區(qū)域的DEM影像數(shù)據(jù)。

1.3 插值方法

從數(shù)量有限的氣象站點獲得了降水?dāng)?shù)據(jù)的實測值，那么怎樣根據(jù)站點的數(shù)據(jù)指標(biāo)信息推算到整個研究區(qū)域，國內(nèi)外已經(jīng)基于這一問題開展了大量的研究，并總結(jié)了多種方法[20-21]。其中就存在應(yīng)用地統(tǒng)計學(xué)結(jié)合GIS平臺闡明數(shù)據(jù)的空間變異性，并對比差異。在實現(xiàn)空間插值的過程中，該文主要采用Kriging、IDW、RBF和GWRK對空間降水進行模擬，并對這4種方法進行對比分析。

1.3.1 降水?dāng)?shù)據(jù)序列的正態(tài)性檢驗

各氣象站點的年降水量數(shù)據(jù)變異程度中等，經(jīng)擬合優(yōu)度檢驗(Jarque-Bera)氣象站點的年降水量數(shù)據(jù)序列(圖2)符合正態(tài)分布。

1.3.2 克里格插值法(Kriging)

克里格插值法(Kriging)是一種經(jīng)典的無偏內(nèi)插估計算法[22]，在地理、環(huán)境與大氣等科學(xué)領(lǐng)域都被應(yīng)用[23]，一般公式：

(1)

式中，Yxi為觀測值；x0為待估觀測點；αi為加權(quán)因子；n為觀測點總數(shù)。加權(quán)因子αi由無偏估計和最小方差來確定[24]，其公式：

(2)

1.3.3 反距離加權(quán)法(IDW)

IDW插值方法的原理是樣本點周圍數(shù)值隨著其到樣本點的距離而發(fā)生變化，會以距離為權(quán)重進行加權(quán)平均，也就是距離采樣點越近的點被賦予權(quán)重越大[26-27]。通過對插值數(shù)據(jù)的模擬和比較，選擇球形模型作為降水的空間預(yù)測模型，其表示公式為：

(3)

式中，Y為降水預(yù)測值；Yi為第i個樣本的數(shù)值(i=1，…，n)；Di為距采樣點距離；p為Di的冪，準(zhǔn)則為最小平均絕對誤差。

1.3.4 徑向基函數(shù)法(RBF)

徑向基函數(shù)法(Radial Basis Functions，RBF)是指一種隨距離變化且徑向?qū)ΨQ的標(biāo)量函數(shù)[28]。此方法利用不同的基函數(shù)來擬合成過已知樣本點的平滑表面[29-30]。其計算公式為：

(4)

式中，Y為估算點降水預(yù)測值；Yi為第i(i=1，2，3，……n)個站點實測值；di為估算點到第i站點的距離；n為參與插值的實測站點數(shù)；G為歐氏距離的單調(diào)函數(shù)。

1.3.5 地理加權(quán)回歸克里格(GWRK)

GWRK是結(jié)合地理加權(quán)回歸(geographically weighted regression，GWR)和克里格插值法(Kriging)的一種混合插值方法[31]。其中，GWR是一種基于距離觀測點的地理加權(quán)回歸系數(shù)來探索空間非平穩(wěn)性的局部地理回歸方法[32]。相較于普通最小乘法(ordinary least squares，OLS)的全局回歸技術(shù)，GWR能客觀地量化空間上的變化關(guān)系。

該文考慮到新疆區(qū)域地形特色，采用GWRK方法，選取經(jīng)度、緯度和高程為建立模型的自由變量，年降水量為模型的因變量進行分析。

1) 根據(jù)氣象站樣本點實測降水量值和GWR模型預(yù)測值得到樣本點殘差值，GWR模型[16]如下：

(5)

式中,YGWR(s0)為點s0處的預(yù)測值；k為自變量的數(shù)目；Xk(s0)為點s0處第k個自變量的測量值；βi(s0)為局部估量的系數(shù)。

2) 利用Kriging法對GWR回歸的殘差進行插值，其插值結(jié)果與GWR結(jié)果進行空間加運算[33]，GWRK法表達式[16]為：

(6)

式中，ε(s0)為插值后得到的殘差結(jié)果。

模型建立前對自由變量進行分析，篩選符合建模要求的因子，發(fā)現(xiàn)緯度和高程的回歸系數(shù)均達到了顯著水平(P<0.05)，但經(jīng)度的P>0.05,故不再將經(jīng)度作為自由變量(表1)。方差膨脹因子(Variance Inflation Factor，VIF)均遠小于7.5，說明自由變量間不存在多重共線性問題,結(jié)果可靠。

表1 模型回歸計算結(jié)果

由表1可知，方差膨脹因子(Variance Inflation Factor，VIF)均遠小于7.5，說明自由變量間不存在多重共線性問題,結(jié)果可靠。

1.3.6 評價方法

1) 模型評價 R2是模型的擬合優(yōu)度的檢驗量，介于0～1，越接近1模型擬合效果越佳，RSS是用于檢驗回歸模型的常用指標(biāo)，其值越小，擬合效果越佳[34]。

(7)

(8)

2) 方法評價 將獲取的降水量數(shù)據(jù)進行基本的預(yù)處理并求出其年均降水量。選取20%樣本點的氣象站點，即12個作為驗證點，對其進行準(zhǔn)確性驗證，并用于空間插值。通過計算準(zhǔn)確性(Accuracy)、均方根誤差(RMSE)和平均相對誤差(MRE)等指標(biāo)評價不同插值方法的預(yù)測精度，公式如下：

(9)

(10)

(11)

(12)

2 結(jié)果與分析

2.1 插值結(jié)果比較

以全新疆區(qū)域66站點的多年年均降水量為指標(biāo)進行降水量空間插值分析，對比發(fā)現(xiàn)Kriging、IDW、RBF和GWRK得出的結(jié)果與精度有所不同。根據(jù)空間數(shù)據(jù)的特征采用交叉驗證的方法分析結(jié)果，用R2、RSS、Accuracy、RMSE和MRE指標(biāo)評價不同插值方法的預(yù)測精度(表2)， GWRK的R2值最接近1,且其RSS值相對最小；從Accuracy看，AccuracyKriging>AccuracyGWRK>AccuracyRBF>AccuracyIDW，即Kriging準(zhǔn)確度最高；從RMSE看，RMSEGWRK

表2 不同插值方法的平均相對誤差比較

2.2 新疆降水量空間分布特征

基于上述4種插值方法生成了降水量的空間分布圖，直觀地分析新疆降水量的空間分布特征。依全局觀察，新疆的降水量北部明顯多于南部，并隨緯度從北向南大致呈帶狀遞減，有明顯的南北分界。從局部觀察來看，南北疆的山區(qū)部分會因地形因素有小范圍的波動，降水量發(fā)生增減，降水量空間分析分布圖見圖3。

3 討論

降水量屬于氣象數(shù)據(jù)，會受到多種驅(qū)動因素的影響，獲取降水量信息單靠從國家氣象觀測站點觀測的數(shù)據(jù)是不嚴(yán)謹(jǐn)?shù)?，存在片面性。隨著研究的不斷深入，已經(jīng)有大量學(xué)者從單純考慮降水?dāng)?shù)據(jù)本身到逐漸加入高程、經(jīng)緯度等影響因素[35]，且研究表明，區(qū)域降水會受到位置、地形以及大氣的制約。該文在研究時主要考慮位置和地形因素，進行多方面的區(qū)域降水探討，從而提高了降水插值的精度，其中，就新疆地區(qū)而言高程因素對降水影響相對較大，這與新疆地區(qū)三山脈夾兩盆地的獨特地形密切相關(guān)。但在實際生活中影響降水空間分布的還有許多隨機因素，加入隨機因素以及大氣因素會不會進一步提高降水插值精度還有待于進一步的探討。

4 結(jié)論

每種空間插值方法都有自己的適用范圍、研究目的、算法和優(yōu)缺點。采用1981—2010年新疆范圍內(nèi)的66個國家基本站點年平均降水?dāng)?shù)據(jù)，分別運用Kriging、IDW、RBF和GWRK探討了新疆降水的空間分布。新疆地區(qū)降水量差別大、氣象臺站少、高程變化大，不能僅僅考慮站點的降水?dāng)?shù)據(jù)和位置因素，故該文加入高程影響因子進一步研究。結(jié)果表明：GWRK得出的結(jié)果要優(yōu)于Kriging、IDW和RBF。采用GWRK插值模型的方法進行新疆區(qū)域的降水分布特征的研究以及降水量的預(yù)測，會比普通的插值方法得到更高的插值精度。

研究表明，有關(guān)降水插值的研究應(yīng)多加考慮位置、地形等相關(guān)因素的影響。針對不同研究區(qū)域的不同空間位置、地理和地形特點，分析研究出相關(guān)性較大的影響因子建模分析，插值分析獲取更高插值精度。

此外，該文在研究時考慮的影響因子較少，后期可以嘗試引入坡度、坡向、地形曲率、氣溫和濕度指數(shù)等進行分析研究。同時，還可以開發(fā)新的插值方法與遙感指數(shù)和雷達數(shù)據(jù)相結(jié)合，如引入歸一化植被指數(shù)(NDVI)、歸一化水體指數(shù)(NDWI)等以提高精度，這在未來還需要進一步研究。