東北三省及內(nèi)蒙古自治區(qū)降水量空間插值方法比較

李 洪，郭 琳

（1.國家林業(yè)和草原局西南調(diào)查規(guī)劃院，云南昆明 650216；2.水資源安全北京實驗室 首都師范大學(xué)，北京 100048）

降水量對氣候變化監(jiān)測意義重大，直接關(guān)乎人類的生產(chǎn)和生活。隨著社會經(jīng)濟發(fā)展和生態(tài)環(huán)境變化，水資源的保護與開發(fā)利用等日益受到關(guān)注[1-2]。降水空間分布信息在水資源管理、災(zāi)害預(yù)測預(yù)警和區(qū)域可持續(xù)發(fā)展等方面的作用也越來越重要[3]。由于人力和經(jīng)濟等因素，目前氣象站點的布設(shè)存在較大局限性，不同空間尺度下降水?dāng)?shù)據(jù)的精度要求難以得到滿足[4]，定點觀測到的數(shù)據(jù)大多不能直接用于其他地點，更不能替代某一較大面積區(qū)域的平均值[4-5]。探索可靠和有效的空間插值技術(shù)，通過已知站點的降水?dāng)?shù)據(jù)進行空間內(nèi)插，轉(zhuǎn)化生成區(qū)域內(nèi)高精度的降水量空間分布圖，尤為重要[4，6]。

降水量具有間斷性和空間不連續(xù)性的特點，且小時降水量零值較多[7]，受多個地理要素的綜合影響，其空間插值表現(xiàn)出較多的不確定性[1，6]。常用的插值方法有多項式法、趨勢面法、反距離加權(quán)法、徑向基函數(shù)法和克里格法等[8]。已有研究表明，當(dāng)數(shù)據(jù)密度足夠大時，各插值方法用于同一組氣象數(shù)據(jù)的插值結(jié)果差異不大，但在數(shù)據(jù)分布稀疏和差異顯著時，插值方法的選擇及其計算參數(shù)的設(shè)定尤為關(guān)鍵[9]。應(yīng)根據(jù)插值數(shù)據(jù)的特征、精度及計算效果等，選擇合適的插值方法[10]。

本研究以東北三省及內(nèi)蒙古自治區(qū)123個氣象站臺的月平均降水?dāng)?shù)據(jù)為基礎(chǔ)，對協(xié)同克里格法、普通克里格法和反距離加權(quán)法3種空間插值方法進行比較分析，尋找適宜東北三省及內(nèi)蒙古自治區(qū)降水量的最佳插值方法，為這些地區(qū)的水資源監(jiān)測及應(yīng)用提供可靠的數(shù)據(jù)支撐。

1 材料與方法

1.1 研究區(qū)概況

研究區(qū)為東北三省及內(nèi)蒙古自治區(qū)（106°40′～132°92′E，38°90′ ～ 52°97′N）。東北三省包括黑龍江、吉林和遼寧3省，是我國的老工業(yè)基地和糧食主產(chǎn)區(qū)，具有綜合的工業(yè)體系、完備的基礎(chǔ)設(shè)施、豐富的農(nóng)產(chǎn)品資源和優(yōu)良的生態(tài)環(huán)境，是極具潛力的富饒之地。內(nèi)蒙古自治區(qū)土地廣闊，位于我國北部邊疆，是中國第3大省區(qū)，也是我國最大的草原牧區(qū)。研究區(qū)跨度大，地形和氣候差異較大，降水分布不均。研究區(qū)地理位置及氣象站點的分布情況如圖1所示。

圖1 研究區(qū)地理位置及氣象站點分布示意圖Fig.1 Geographical locations and distributions of meteorological stations in study area

1.2 數(shù)據(jù)分析與預(yù)處理

1.2.1 數(shù)據(jù)源

選取東北三省及內(nèi)蒙古自治區(qū)123個氣象站的月平均降水?dāng)?shù)據(jù)；其中，黑龍江省32個氣象站，遼寧省25 個氣象站，吉林省25 個氣象站，內(nèi)蒙古自治區(qū)41 個氣象站。用于空間插值的資料，還包括各站的經(jīng)緯度、海拔高度和平均氣溫，各氣象站的空間分布情況如圖2所示，各數(shù)據(jù)的統(tǒng)計見表1。氣象站分布不均勻，內(nèi)蒙古自治區(qū)較稀疏，吉林和遼寧省較密集。降水量、海拔和平均氣溫均波動較大。

表1 空間內(nèi)插數(shù)據(jù)統(tǒng)計Tab.1 Statistics of spatial interpolation data

圖2 研究區(qū)氣象數(shù)據(jù)的空間分布特征Fig.2 Spatial distribution characteristics of meteorological data in study area

1.2.2 降水量數(shù)據(jù)空間分布特征

降水量數(shù)據(jù)的空間分布特征直接影響降水量內(nèi)插方法的選擇。通過ArcGIS平臺構(gòu)建123個氣象站降水量數(shù)據(jù)的泰森多邊形，獲得研究區(qū)平均降水量空間分布信息（圖3）；通過趨勢分析工具，以降水量為高度構(gòu)建三維透視圖（圖4）；通過直方圖、QQ?plot 圖等，對降水量數(shù)據(jù)分布的整體趨勢、分布均衡性、正態(tài)性、聚集性和隨機性進行分析，掌握降水量空間分布的基本特征。

圖3 降水量數(shù)據(jù)泰森多邊形Fig.3 Tyson polygons of precipitation data

圖4 降水量數(shù)據(jù)整體分布趨勢Fig.4 Overall distribution trends of precipitation data

降水量從西至東呈增加趨勢，在南-北方向上呈U 型分布，存在明顯趨勢性。在數(shù)據(jù)的自相關(guān)分析中，考慮數(shù)據(jù)的Moran′s ⅠIndex = 0.4，信度Z 為9.05，參照統(tǒng)計先驗信息，數(shù)據(jù)有高度的自相關(guān)性。氣象站的分布均衡性是判斷內(nèi)插數(shù)據(jù)質(zhì)量的一個重要指標(biāo)；正態(tài)分布是多種隨機數(shù)據(jù)處理的基礎(chǔ)假設(shè)，數(shù)據(jù)如果能符合這種假設(shè)，則能應(yīng)用多種基于隨機數(shù)據(jù)思想的內(nèi)插方法進行內(nèi)插。降水量數(shù)據(jù)本身不呈正態(tài)分布，需對數(shù)據(jù)進行變換處理，使數(shù)據(jù)更趨于正態(tài)分布。有學(xué)者運用對數(shù)、立方根和平方根等變換方式對數(shù)據(jù)進行預(yù)處理，提高空間插值預(yù)測的精度[11]。本研究中，首先對樣本數(shù)據(jù)進行整體分析，揭示降水量的空間分布特征，并選擇Log 變換方式對降水量數(shù)據(jù)進行預(yù)處理。

1.2.3 降水量數(shù)據(jù)與其他數(shù)據(jù)相關(guān)性分析

從機理上考慮，降水量受經(jīng)度、緯度、海拔和氣溫等多方面因素的影響，對降水量數(shù)據(jù)與這些數(shù)據(jù)進行相關(guān)分析（表2）。降水量與經(jīng)度呈極顯著正相關(guān)（P< 0.01）；與平均氣溫呈極顯著負(fù)相關(guān)（P<0.01）；與海拔和緯度無明顯相關(guān)性。通過繪制散點圖直觀分析降水量與海拔、經(jīng)緯度和平均氣溫的變化趨勢（圖5）。隨經(jīng)度增大，降水量逐漸增多，可能是因為經(jīng)度高的區(qū)域更靠近海洋，更濕潤。緯度上整體變化不大，可能是因為本研究區(qū)域基本屬于中溫帶，緯度跨度不大。當(dāng)月平均氣溫為20～27 ℃時，降水量與氣溫有明顯趨勢關(guān)系，氣溫越高，降水量越低。綜合分析，降水量數(shù)據(jù)與其他數(shù)據(jù)存在一定相關(guān)關(guān)系，且自身分布存在明顯特征，具有良好的正態(tài)性，適于用多種內(nèi)插方法進行空間插值。

表2 降水量數(shù)據(jù)與其他數(shù)據(jù)的相關(guān)性Tab.2 Correlations among precipitation data and other data

圖5 降水量與海拔、經(jīng)緯度和平均氣溫分布的散點圖Fig.5 Distribution scatter plots of precipitation with altitude,latitude,longitude and average temperature

1.3 離群值剔除

離群值產(chǎn)生的原因有很多，可能是真實且正常的數(shù)據(jù)，與其他觀測值屬同一總體；也可能是由于測量時的偶然性，或記錄和計算時的失誤，產(chǎn)生的非正常、錯誤的數(shù)據(jù)，與其他觀測值不屬于同一總體，在插值計算前需將其改正或剔除[12-14]。本研究通過直方圖、半變異/協(xié)方差函數(shù)云和Voronoi 圖查找出的離群值基本相同，分別分布于研究區(qū)東北方向和南部的兩個樣本點，將其剔除。

2 空間插值方法選取及依據(jù)

考慮降水量的特殊性，即影響降水的因素較多，包括氣象站的經(jīng)緯度、海拔高度、坡度、坡向及風(fēng)速等[3]。本研究選擇協(xié)同克里格法（Co-kriging）、普通克里格法（Ordinary kriging，OK）和反距離加權(quán)法（In?verse distance weighting，IDW）3 種方法對東北三省及內(nèi)蒙古自治區(qū)的降水量進行空間內(nèi)插（圖6）。將離散站點的觀測數(shù)據(jù)轉(zhuǎn)換為連續(xù)的空間數(shù)據(jù)曲面，以便更快捷和精確地獲取未知點的降水量數(shù)值。

圖6 空間插值方法Fig.6 Spatial interpolation methods

3 插值精度評價

在進行氣象數(shù)據(jù)空間插值時，交叉驗證法（Cross validation）是應(yīng)用最廣泛的一種精度評價方法[18]。本研究通過計算估計值和實測值間的差值，采用預(yù)測誤差（Prediction error）中的指標(biāo)作為交叉驗證法的評價指標(biāo)，進行插值結(jié)果精度評價和模型合理性分析[19-20]。理想模型的標(biāo)準(zhǔn)為標(biāo)準(zhǔn)平均值（Mean standardized）最接近0，均方根預(yù)測誤差（Root-mean-square）最小，平均標(biāo)準(zhǔn)誤差（Average standard error）最接近均方根預(yù)測誤差，標(biāo)準(zhǔn)均方根預(yù)測誤差（Root-mean-square standardized）最接近1，實測值與預(yù)測值更接近[18，21]。

3種插值方法中，反距離加權(quán)法的均方根預(yù)測誤差最?。?4.750 0 mm）（表3）。協(xié)同克里格法的標(biāo)準(zhǔn)平均值為0.000 2，最接近0；平均誤差最?。?.011 1）；標(biāo)準(zhǔn)均方根預(yù)測誤差為0.863 5 mm，最接近1；均方根預(yù)測誤差為46.320 0 mm，與平均標(biāo)準(zhǔn)誤差53.790 0 mm相差較小。

表3 空間插值結(jié)果的精度評價Tab.3 Accuracy evaluation of spatial interpolation results

4 插值結(jié)果分析

根據(jù)數(shù)據(jù)分析，原始降水量數(shù)據(jù)的空間分布存在二次曲線趨勢，且不呈正態(tài)分布；克里格插值過程中，選擇Log 變換方式，使其呈正態(tài)分布，并移除二次趨勢面；協(xié)同克里格插值過程中，數(shù)據(jù)中加入經(jīng)緯度、海拔和平均氣溫作為影響變量，其他參數(shù)選擇與普通克里格一致。3 種插值方法得到的研究區(qū)降水量分布圖如圖7所示。

圖7 空間插值結(jié)果（a:IDW；b:OK；c:Co-kriging）Fig.7 Spatial interpolation results(a:IDW;b:OK;c:Co-kriging)

3 種插值方法所得結(jié)果趨勢大致相同，且呈現(xiàn)兩個極值地區(qū)，黑龍江省北部和遼寧省南部地區(qū)為降水高值區(qū)，內(nèi)蒙古自治區(qū)西南部和黑龍江省東部地區(qū)為降水低值區(qū)?？傮w上，內(nèi)蒙古自治區(qū)的降水比其他省份少，與研究區(qū)的實際情況吻合。

5 結(jié)論與討論

將反距離加權(quán)法、普通克里格法和協(xié)同克里格法應(yīng)用于東北三省及內(nèi)蒙古自治區(qū)的降水量插值運算，并對其插值效果和預(yù)測精度進行分析。結(jié)果顯示，綜合3種方法的插值結(jié)果與精度評價，協(xié)同克里格法的插值效果最好，各項評價指標(biāo)均較符合要求，插值效果最平滑。

反距離加權(quán)法的均方根預(yù)測誤差最小，但此方法易受極值影響，插值空間分布圖中易出現(xiàn)“牛眼”現(xiàn)象，相鄰站點的空間依賴關(guān)系表現(xiàn)不佳。協(xié)同克里格法綜合考慮了經(jīng)緯度、海拔和平均氣溫等與降水量有關(guān)的影響因素，并將其納入數(shù)據(jù)集一起參與插值估算，在其他因素均對降水量有影響的情況下，協(xié)同克里格法的插值效果及插值精度更理想。此外，在插值前，對數(shù)據(jù)的分布趨勢及分布規(guī)律進行分析，對離群值進行查找和剔除均非常必要，有助于插值方法及參數(shù)的選擇。

無論哪種插值方法均需在特定的假設(shè)條件和適用范圍內(nèi)應(yīng)用，與數(shù)據(jù)自身特點和空間分布規(guī)律有關(guān)。本研究在進行降水量插值研究時，僅考慮經(jīng)緯度、海拔及平均氣溫的影響，未將其他影響因素（風(fēng)速、坡度和坡向等）納入考慮；除插值方法外，站點數(shù)量及分布均勻程度、降水量觀測時間間隔尺度和柵格單元大小等因素也會造成空間信息的不確定性。在以后的研究中，還需對這些因素進行深入分析。