孫偉紅，常直楊

(1.江蘇省測繪工程院，江蘇 南京 210013； 2.南京旅游職業(yè)學(xué)院，江蘇 南京 211100)

基于GIS的西秦嶺地區(qū)降水量空間插值方法探討

孫偉紅1，常直楊2

(1.江蘇省測繪工程院，江蘇 南京 210013； 2.南京旅游職業(yè)學(xué)院，江蘇 南京 211100)

反距離加權(quán)法；綜合插值法；克里金法；降水量；西秦嶺地區(qū)

以西秦嶺地區(qū)周邊77個氣象站1971—2000年平均降水量為數(shù)據(jù)源，選取反距離加權(quán)法(IDW法)、綜合插值法、泛克里金法(UK法)、回歸克里金法(RK法)等4種內(nèi)插方法，探討了西秦嶺地區(qū)降水量最佳插值方法的選取。交叉驗證結(jié)果表明：4種插值方法的精度,IDW法<綜合插值法

研究某一地區(qū)降水量年際變化特征及其空間分布，對于指導(dǎo)當(dāng)?shù)剞r(nóng)業(yè)生產(chǎn)、預(yù)報地質(zhì)災(zāi)害等具有重要意義[1-2]，但是在地形復(fù)雜的偏遠地區(qū)，受經(jīng)濟、技術(shù)和人力條件限制，氣象站點數(shù)量有限，獲取氣象資料難度較大。近些年隨著相關(guān)技術(shù)的發(fā)展，通過空間插值生成區(qū)域降水要素空間分布圖[3-6]成為一種有效的解決辦法。在眾多的空間插值方法中，最近距離法、反距離加權(quán)法等[3]體現(xiàn)了空間幾何距離的影響；趨勢面法、多元回歸法[4]體現(xiàn)了屬性數(shù)據(jù)之間的統(tǒng)計關(guān)系；克里金法[5-6]體現(xiàn)了氣象要素的原始數(shù)據(jù)及變異函數(shù)的結(jié)構(gòu)特點。每種方法都有其使用范圍、算法和優(yōu)缺點，沒有絕對最優(yōu)的空間內(nèi)插方法，只有特定條件下的最優(yōu)方法[7-9]。

西秦嶺地區(qū)地質(zhì)構(gòu)造活躍，滑坡和泥石流等地質(zhì)災(zāi)害多發(fā)，研究降水量空間分布對于該地區(qū)水土流失治理和地質(zhì)災(zāi)害防治尤為重要。本研究選用4種空間插值方法——反距離加權(quán)法(IDW法)、綜合插值法、泛克里金法(UK法)、回歸克里金法(RK法)，通過對比分析，選出適合該地區(qū)的最佳的空間插值方法，以更好地描述西秦嶺地區(qū)的降水量空間分布。

1 材料和方法

1.1 研究區(qū)概況

西秦嶺地區(qū)位于陜西省秦嶺以西、甘南山地高原以東，甘肅、陜西、四川三省交界處(圖1)，主體位于甘肅省隴南地區(qū)，面積約46 394 km2。全區(qū)海拔450～4 591 m，地勢起伏大，近20年來曾多次出現(xiàn)大暴雨和特大暴雨，年降水量的80%～87%集中在4—9月，年平均暴雨1.5次，豐水年可達3～5次。

圖1 西秦嶺地區(qū)地理位置示意

1.2 研究數(shù)據(jù)

研究所需SRTM3-DEM數(shù)據(jù)下載自計算機網(wǎng)絡(luò)信息中心國際科學(xué)數(shù)據(jù)鏡像網(wǎng)站(http://datamirror.csdb.cn/)，經(jīng)緯度柵格圖層均由DEM提取，所有圖層都采用相同的投影坐標(biāo)系統(tǒng)。

氣象站臺經(jīng)緯度、海拔及歷年年平均降水量數(shù)據(jù)均來自中國氣象科學(xué)數(shù)據(jù)網(wǎng)，時間段為1971—2000年，氣象站臺中分布于西秦嶺地區(qū)的有23個，分布于鄰近區(qū)域的有54個，總計77個(圖2)。其中，考慮插值精度檢驗，隨機選取氣象站臺總數(shù)的20%(14個氣象站臺)的降水量作為精度評估檢驗站點數(shù)據(jù)，其余80%的降水量作為實驗站點數(shù)據(jù)。

圖2 西秦嶺地區(qū)氣象站臺空間分布

1.3 插值方法

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

反距離加權(quán)法是基于相近相似的原理[3]，即兩個物體離得越遠性質(zhì)就越不相似，反之，則相似性越高。計算時，以樣本點與插值點間的距離作為權(quán)重進行加權(quán)平均，離插值點越近的樣本點被賦予的權(quán)重越大。

1.3.2 綜合插值法(改進后的多元線性回歸法)

大尺度的降水與海拔、海陸相對位置有關(guān)[2]，因此海拔、經(jīng)緯度等都可能影響降水量的分布。利用多元線性回歸法進行區(qū)域整體插值，多元線性回歸方程為

Y=a0+a1λ+a2Φ+a3h

(1)

式中：Y為多元線性回歸模擬的降水量，mm；λ為經(jīng)度；Φ為緯度；h為海拔；a0、a1、a2、a3為參數(shù)。

由于采用多元線性回歸進行降水量插值模擬誤差較大，不能區(qū)分不同區(qū)域的降水特征，因此將回歸統(tǒng)計模型與插值模型相結(jié)合，即采用綜合插值方法[10]進一步提高降水量模擬精度。選擇IDW法插值結(jié)果與回歸模擬結(jié)果分別作為兩個影響因子，做二者與降水量實測值的多元線性回歸，建立多元線性回歸方程為

Y′=a0′+a1′I+a2′Y

(2)

式中：Y′為改進后的多元線性回歸模擬的降水量，mm；I為用IDW法插值的降水量結(jié)果，mm；a0′、a1′、a2′為參數(shù)。

1.3.3 泛克里金法(UK法)

泛克里金法是以地理坐標(biāo)X和Y的趨勢面方程分離趨勢項，再對分離漂移趨勢項后的殘差進行普通克里金預(yù)測[6]，最后將趨勢面方程的預(yù)測結(jié)果與殘差的普通克里金預(yù)測結(jié)果相加，泛克里金方法中的漂移趨勢項是地理坐標(biāo)X和Y的函數(shù)。

1.3.4 回歸克里金法(RK法)

影響氣候要素空間分布的地理因子可分為兩類[11]：一類是大氣候因素或宏觀地理因子，例如大的山脈走向、地勢高低和地理位置(經(jīng)度、緯度和距海遠近)等；另一類是小氣候因素或微觀地形因子，如地形遮蔽度、坡度、坡向和地形起伏量等。前者影響氣候要素的空間分布趨勢(稱為宏觀本底值)，可用經(jīng)度、緯度和海拔3個變量來表述；后者影響氣候要素的局地變化(稱為小地形訂正，即殘差項)，可借助ArcGIS軟件計算生成。使用公式可以表示為

P=F(λ，Φ，h)+ξ

(3)

式中：P為實測降水量，mm；λ為經(jīng)度；Φ為緯度；h為海拔，m；F(λ，Φ，h)為降水的趨勢值；ξ為殘差項。

將F(λ，Φ，h)展開成三維二次趨勢面方程為

F(λ，Φ，h)=b0+b1λ+b2Φ+b3h+b4λΦ+b5λh+

b6Φh+b7λ2+b8Φ2+b9h2

(4)

式中：b0，b1，…，b9為待定系數(shù)。

利用逐步回歸優(yōu)化模型，這種插值方法稱為回歸克里金法。

2 精度檢驗

采用交叉驗證方法對插值結(jié)果進行對比分析，準(zhǔn)確驗證不同方法的相對精度。交叉驗證法是先假定每一站點的降水量數(shù)值未知，用周圍站點的數(shù)值來估算，再計算估計值與實測值的差值。本研究對14個檢驗站點(圖1)的插值結(jié)果進行了交叉驗證，采用平均絕對誤差(MAE)、平均相對誤差(MRE)、均方根誤差(RMSE)作為評價指標(biāo)，計算公式分別為

(5)

(6)

(7)

式中：n為檢驗站點數(shù)量；Pbi、Pai分別為第i個檢驗站點降水量的估算值和實測值，mm。

MAE值定量地給出誤差，反映了估計值的實測誤差范圍。MRE能通過反映不同數(shù)據(jù)量或不同要素的誤差相對值，定性地給出誤差范圍，效果更直觀。RMSE值反映了利用樣點數(shù)據(jù)估值的極值效應(yīng)和靈敏度，其值越小模型越準(zhǔn)確。

3 結(jié)果分析

3.1 降水量空間插值結(jié)果

本研究中將IDW法的站點搜索范圍設(shè)為12，距離權(quán)重指數(shù)設(shè)為2，利用ArcGIS和SPSS軟件得出基于IDW法、綜合插值法的西秦嶺地區(qū)多年平均降水量空間分布結(jié)果見圖3(a)、(b)。

利用UK法進行空間插值時，西秦嶺地區(qū)多年平均降水量三維二次趨勢面方程為

F(λ,Φ,h)=-15 372.029+852.16Φ+3.644λ2-

20.541λΦ+2.249E-5h2+17.369Φ2

(8)

在對殘差項ξ進行普通克里金插值時，模型和參數(shù)的選擇可以影響到插值的精度，為了選擇合適的模型、基臺值、變程值等參數(shù)，利用GS+軟件對殘差數(shù)據(jù)進行了半方差分析。半方差分析顯示球狀模型效果最好，R2=0.958，塊金值為770，基臺值為7 111，變程值為242 300 m，偏基臺值與基臺值之比為89.2%，空間相關(guān)性強。最后將三維二次趨勢面方程的預(yù)測結(jié)果與殘差的普通克里金預(yù)測結(jié)果相加，得出多年平均降水量空間分布，見圖3(c)。

利用RK法進行空間插值時， 相關(guān)性分析結(jié)果顯示，經(jīng)緯度與氣象站臺多年平均降水量、海拔之間均為極顯著相關(guān)，表明西秦嶺地區(qū)多年平均降水量的空間分布與空間位置(坐標(biāo))極顯著相關(guān)。趨勢面方程為

圖3 西秦嶺地區(qū)多年平均降水量空間分布

F=2 253.699-121.403Φ+23.936λ

(9)

式中：F為多年平均降水量趨勢項；λ為經(jīng)度；Φ為緯度。

由于研究區(qū)經(jīng)緯度和海拔之間也存在極顯著相關(guān)性，因此趨勢面方程在一定程度上也考慮了海拔的影響。對殘差進行分析得出最佳模型和參數(shù)，然后將趨勢面方程的預(yù)測結(jié)果與殘差的普通克里金預(yù)測結(jié)果相加，得出多年平均降水量的空間分布，見圖3(d)。

3.2 誤差分析

多年平均降水量空間插值交叉檢驗誤差分析結(jié)果見表1。從表1可知，考慮了未知樣本點的空間方位與已知樣本點的空間分布關(guān)系的UK法、RK法的精度比IDW法、綜合插值法高；RK法較UK法精度略高，但提高幅度并不大，RK法插值計算精度參數(shù)MAE=43.44 mm，MRE=7.4%，RMSE=57.35 mm。預(yù)測值與實測值的散點關(guān)系圖如圖4所示，其中RK法預(yù)測值與實測值之間的相關(guān)系數(shù)最大，為0.807 9。綜合MAE、MRE、RMSE及相關(guān)系數(shù)的分析結(jié)果(表1、圖4)，整體來看，各插值方法精度RK法>UK法>綜合插值法>IDW法。這表明，在西秦嶺地區(qū)利用RK法對多年平均降水量進行空間插值時，直接將海拔作為影響因子考慮進去，能夠提高插值精度；在西秦嶺地區(qū)，多年平均降水量的空間分布主要與空間位置(地理坐標(biāo))有關(guān)。

表1 不同插值方法精度比較

4 結(jié) 論

不同的插值方法對于不同地區(qū)、不同研究數(shù)據(jù)有其特定的應(yīng)用范圍。對比分析了常用的4種空間插值方法，分析得出西秦嶺地區(qū)多年平均降水量空間插值方法的精度RK法>UK法>綜合插值法>IDW法。

圖4 4種插值方法實測值與預(yù)測值關(guān)系散點圖

[1] 尚宗波,高瓊,楊奠安.利用中國氣候信息系統(tǒng)研究年降水量空間分布規(guī)律[J].生態(tài)學(xué)報,2001,21(5):689-694.

[2] 朱會義,賈紹鳳.降水信息空間插值的不確定性分析[J].地理科學(xué)進展,2004,23(2):34-42.

[3] TSANIS I K,GAD M A.A GIS precipitation method for analysis of storm kinematics[J]. Environmental Modelling & Software,2001,16(3):273-281.

[4] WACKEERNAGEL H.Multivariate Geostatistics[M].Berlin:Springer-Verlag,1995:256.

[5] 邵曉梅,嚴(yán)昌榮,魏紅兵.基于Kriging插值的黃河流域降水時空分布格局[J].中國農(nóng)業(yè)氣象,2006,27(2):65-69.

[6] 蘇曉燕,趙永存,楊浩,等.不同采樣點數(shù)量下土壤有機質(zhì)含量空間預(yù)測方法對比[J].地學(xué)前緣,2011,18(6):34-40.

[7] SINGH V P,CHOWDHURY P K.Comparing some methods of estimating mean areal rainfall[J].Water Resources Bull,1986,22(2):275-282.

[8] 李新,程國棟,盧玲.空間內(nèi)插方法比較[J].地球科學(xué)進展,2000,15(3):260-265.

[9] 趙傳燕,馮兆東,南忠仁.隴西祖厲河流域降水插值方法的對比分析[J].高原氣象,2008,27(1):208-214.

[10] 劉勁松,陳輝,楊彬云,等.河北省年均降水量插值方法研究[J].生態(tài)學(xué)報,2009,29(7):3493-3500.

[11] 林忠輝,莫興國,李宏軒,等.中國陸地區(qū)域氣象要素的空間插值[J].地理學(xué)報,2002,5(1):47-56.

(責(zé)任編輯 李楊楊)

江蘇省高校自然科學(xué)研究面上項目(16KJB170014)

TP79

A

1000-0941(2017)06-0054-03

孫偉紅(1986—)，女，山西運城市人，工程師，碩士，主要從事地理信息制圖研究工作。

2016-11-15