東北地區是中國最大的商品糧基地和農業生產最具發展潛力的地區之一，同時也是中國重要的工業和能源基地。 東北地區位于北半球的中高緯度，是我國緯度最高的地區，是世界著名的溫帶季風氣候區，是典型的氣候脆弱區和受氣候變暖影響最為敏感的地區。 近年來針對東北地區氣候變化已展開不少研究。 已有研究較多地關注氣溫、降水等方面，李莎等將時空 Kriging 方法應用于東北地區氣溫空間插值研究中；賀偉等將 Morlet 小波分析應用于東北地區氣溫和降水變化趨勢研究中；孫力等針對東北地區夏季降水變化的時空分布及變化規律進行研究。 這些研究對于正確認識東北地區氣候變化規律具有重要意義。 但是以往的研究大多都是基于氣溫或降水的單一氣象要素，而對東北地區空氣相對濕度的研究卻較少。 在氣候資源中，氣溫、降水和濕度與農作物及生物的生長發育有密切關系。 氣溫和降水固然是氣候資源中最主要的兩個因子，但也不能忽視空氣濕度這個氣候因子。 空氣濕度的大小是形成天氣演變主要因素之一， 與生產、生活有密切關系，如農作物病蟲害的發生直接取決于濕度的大??；空氣濕度過大或過小對人體健康也有很大影響等。

本研究對東北三?。ê邶埥?、吉林、遼寧）的空氣相對濕度進行空間插值，利用有限的空氣濕度站點數據來估算未知點的空氣濕度。 目前已發展了較多空間內插方法，如泰森多邊形、克里金（Kriging）內插法、反距離加權平均法、趨勢面分析法、多項式回歸法等。 這些方法只是局限于觀測站點的濕度值，沒有考慮到降水量、氣溫、地形地貌、大氣環流等因素對空氣濕度的影響。 而空氣濕度與降水量具有直接的關系，我國一般通過降水量來劃分濕潤區和干旱區。800 mm 的年降水量為濕潤區與半濕潤區界線，400 mm 的年降水量為半濕潤和半干旱區界線，200 mm 的年降水量為半干旱與干旱區界線。

因此，在進行空氣濕度插值時，考慮降水量的影響可以提高插值精度。 為此，本研究針對空氣相對濕度，選取不同的樣本數，結合降水量的觀測值，采用協同克里金方法對空氣相對濕度進行插值，并對插值結果的精度與普通克里金插值方法得到的結果進行對比分析與探討，同時對不同樣本數的插值結果進行對比。

1、研究區域概況與數據處理

1．1研究區域

研究的是中國東北三省 87 個氣象站點， 站點地處 38．90°-52．97°N，119．70°-132．97°E 之間。 地面觀測臺站所觀測的資料來自 1970－2009 年東北三省日空氣相對濕度和日降水量。 由于氣象臺站觀測資料不可避免存在缺測，本研究剔除了缺測較多的臺站，共選擇 77 個資料較好的觀測臺站。 空氣相對濕度選用不同的樣本數，選用兩種不同數量的觀測值分別是 77 個和 46 個， 而降水量共選用 77 個臺站的觀測值。 具體情況如圖 1 所示。

1．2數據處理

結合降水量對東北三省空氣相對濕度進行插值，以得到更加豐富和準確的空氣濕度信息。 利用數據相對較全的 1970－2009 年的空氣相對濕度觀測數據， 也就是對觀測站的空氣濕度的 40 年平均值分別進行普通克里金插值和協同克里金插值，并分別選取不同的樣本數量，對不同樣本數量的兩種不同空間插值方法進行比較。

40 年的平均空氣相對濕度分別選取 46 個樣本數和 77 個樣本數時的統計結果（表 1）。 從表 1 中可以看到，通過 Kolmogorov－Smirnov 法對兩種不同樣本數的平均空氣相對濕度分別進行正態檢驗的p－value 分別為 0．236 0 和 0．287 9， 均大于 0．05，因此認為符合正態分布。

當選擇多個變量信息進行插值時，第一信息有限時，合理有效地選擇第二信息輔助插值是提高預測精度的關鍵。 輔助因子的選擇需要符合較易獲取、與第一信息關聯度強等特點。 根據研究目的，將空氣濕度作為主變量， 降水量作為輔助變量，將空氣濕度和降水量這兩種變量進行 Person 相關系數分析。 從表 1 可以看到，當空氣濕度選取 77 個觀測值時，這兩種變量的相關系數為 0．548，同時達到了 0．01 顯著水平；當空氣濕度選取 46 個觀測值時，這兩種變量的相關系數為 0．596， 并且達到了 0．01的顯著水平， 因而空氣濕度與降水量密切相關，可以將降水量作為提高空氣濕度插值精度的輔助變量。

2、克里金插值及降水量信息利用

空間插值是通過已知樣本值估計未知點的屬性值或由不規則觀測點數據生成連續表面或等值線。 它的本質是根據給定的一組或多組離散點的屬性值與空間位置擬合一個函數，該函數能充分反映樣本點屬性值與空間位置間的數學關系，從而推斷研究區域范圍內未知空間點的屬性值。 空間插值結果的精度受到插值模型、模型算法、尺度范圍和樣本代表性等多種因素的影響。

2．1區域化變量

當一個變量呈空間分布時，就稱之為區域化變量。 這種變量反映了空間某種屬性的分布特征。 地質、海洋、土壤、氣象、水文、生態等領域都具有某種空間屬性。 區域化變量具有雙重性，在觀測前區域化變量 Z（x）是一個隨機場，觀測后是一個確定的空間結構函數值。

區域化變量具有兩個最顯著的特征，即隨機性和結構性。 一是區域化變量 Z（x）是一個隨機函數，多數時候它具有局部的、隨機的、異常的特征；其次是區域化變量具有一般的或平均的結構性質，即變量在點 X 與偏離空間距離為 h 的點 X＋h 處的隨機量 Z（x）與 Z（x＋h）具有某種程度的自相關性，而且這種自相關性依賴于兩點間的距離 h 與變量其他結構特征。 在某種意義上說這就是區域化變量的結構性特征。

2．2普通克里金方法

克里金法是地統計學的主要內容之一，從統計意義上說，是從變量相關性和變異性出發，在有限區域內對區域化變量的取值進行無偏、最優估計的一種方法；從插值角度講是對空間分布的數據求線性最優、無偏內插估計的一種方法。 首先假設區域化變量 Z（x）滿足二階平穩假設和本征假設，其數學期望為 m，協方差函數 c（h）及變異函數 γ（h）存在。 即

假設在待估計點 x0的臨域內共有 n 個實測點，即 x1，x2， … ，xn。 那 么 ， 普通克里 金法 （OrdinaryKriging）的插值公式：

式中，Z*（x0）為待估計空氣濕度柵格值 ；n 為用于空氣濕度插值站點數；Z（xi）為第 i 個站點平均濕度值；λi為賦予第 i 個站點平均濕度的一組權重系數。 引入拉格朗日系數 u 進行推導可得：

因此，權重系數（λi，i＝1，…，n）和拉格朗日乘數u 可由以上式（5）求得，從而由式（4）可得研究區域內任意點的插值估計。

2．3引入降水量信息的協同克里金方法

協同克里金法（CoKriging）是普通克里金法的擴展，把區域化變量的最佳估值方法從單一屬性發展到兩個或兩個以上的協同區域化屬性，用一個或多個輔助變量對所感興趣的變量進行插值估算，這些輔助變量與主要變量都有相關關系，并且假設變量之間的相關關系能用于提高主要預測值的精度。 當研究區內的輔助信息（如降水量）較容易獲取且變化平穩時，可將這類信息作為輔助影響因素引入協同克里金方法， 用協同克里金法對東北三省年均濕度進行插值時， 將年均降水量作為一個協變量有利于區域估值的結果， 在計算中要用到兩屬性各自的變異函數和協變異函數。 在二階平穩假設下，其期望為：

式中，Z2,CK*（x0）為 x0處空氣濕度估算值；Z2（x2j）為各點空氣濕度；λ2j為賦予各個站點空氣濕度的一組權重系數；Z1（x1i）為各點降水量；λ1i為賦予各個站點降水量的一組權重系數；N1、N2分別為降水量和空氣濕度用于濕度插值站點數，其中 N1＞N2。 引入兩個拉格朗日系數 u1和 u2進行推導可得：

式中 γ11和 γ22分別是 Z1和 Z2的變異函數模型，γ12和 γ21是這兩個變量的變異函數模型， 其中，γ12（h）＝γ21（h）求解線性方程組式（9）即可獲得權重系數（λ1i，i＝1，2，…，N1；λ2j，j＝1，2，…，N2）以及兩個拉格朗日乘數 u1和 u2可由以上兩式求得，從而由（8）式可得研究區域內任意點的插值估計。

2．4估值精度評價方法

采用交叉驗證（Cross－validation）比較不同模型和尋找最佳的插值結果。 在交叉驗證過程中， 首先將觀測值 Z（xi）暫時去除， 然后通過其他觀察值Z（x1），…，Z（xn）和變異函數模型來預測 Z（xi）值， 最后將預測值放回到原始數據中去， 重復以上過程直到對所有觀測點進行估值， 所以在所有的觀測點處既有實際觀測值又有估計值，用統計方法一一比較實際值和預測值。 均方根誤差（RMSE）可以用來評價預測值與觀測值的接近程度。 利用協同克里金（COK）相對于普通克里金（OK）的均方根誤差（RMSE）減少的百分數（RRMSE）表示預測精度的提高程度。

式中，Z（xi）和 Z*（xi）分別為在 xi處的測量值和預測值；n 為樣點數。

3、結果與分析

從表 2 可以看出，當空氣濕度觀測值個數為 46時， 協同克里金的 RMSE 值比采用了全部觀測值（77 個） 的協同克里金的 RMSE 值降低了 20．58％；比普通克里金的 RMSE 值降低了 7．66％。試驗證明，引入降水量作為協變量的協同克里金與普通克里金相比可以提高插值精度。 同時，也可以發現，當空氣濕度樣本數為 46 個時，比樣本數為 77 個的插值結果要好。 因此，在樣本采樣存在困難時，可以采用協同克里金插值，減少采樣點個數，從而節省人力、物力和財力。

圖 2a 是選取 77 個站點， 經普通克里金內插的東北三省平均濕度分布圖， 圖 2b 是選取 77 個站點，由協同克里金插值后的東北三省平均濕度分布圖。圖 3a 是選取 46 個站點，經普通克里金內插的東北三省平均濕度分布圖， 圖 3b 是選取 46 個站點，由協同克里金插值后的東北三省平均濕度分布圖。

由圖中可以看出，本試驗區內最大空氣濕度值集中分布在東北平原，并沿東、西、北逐漸降低。 最小值則集中在西邊的大興安嶺一帶。

降水量作為協變量來進行 CoKriging 的插值，從表 2 可以看到，與 Ordinary Kriging 相比可以提高插值精度。 當空氣濕度的樣本數為 77 個觀測值時，CoKriging 方法比 Ordinary Kriging 的插值方法的均方根誤差 RMSE 值降低了 2．02％； 當減少主變量的樣本數， 空氣濕度的樣本數從中選取較完整的 46個站點時，CoKriging 方法比 Ordinary Kriging 的插值方法的均方根誤差 RMSE 值降低了 7．66％。 從中可以發現，引入了高程作為協變量的協同克里金插值方法效果更佳。 另一方面，減少了主變量的（只選取 46 個站點的觀測值）的 CoKriging 方法比 77 個站點的 CoKriging 方法均方根誤差 RMSE 降低了20．58％，插值結果更好。 因此 ，試驗結果表明，CoKriging 不但可以提高插值精度 ， 而且可以在一定程度上減少實地觀測樣本數，從而可以省時省力節省資金。

4、小結與討論

以東北三省東經 119．7°～132．97°，北緯 38．9°~52．97°區域為試驗區。 利用該地區 1970－2009 年全年氣象站的空氣相對濕度觀測數據的 40 年平均值， 分別采用普通克里金和協同克里金空間內插方法， 對兩種不同方法的實驗結果進行對比分析。

而且選用不同的樣本數量來估計該區的空氣濕度，對不同樣本數量的空氣濕度的協同克里金的插值結果進行對比分析。 試驗結果表明，CoKriging 不但可以提高插值精度，而且可以在一定程度上減少實地觀測樣本數，從而可以節省人力財力和物力。

由于空氣濕度不僅與降水量有關， 還與氣溫、日照時間、經度緯度、海拔等都有關系，要想提高空氣相對濕度的插值精度，還要進一步考慮這些因素的影響。 同時，由于氣象臺站分布稀疏，且空間分布不均勻，而且不可避免的存在缺測，導致局部地區誤差較大，對空間插值而言，要想提高精度，能夠反映數據空間變化趨勢和周期的合理的采樣設計是必要的前提。 在后續研究中，對缺測數據，可以先進行填補再進行插值，或者采用專門針對這種復雜數據進行處理的統計模型來進行研究。 而且東北地區冬天干燥夏天濕潤，在下一步研究中還可以考慮時間的周期性， 在空間插值過程中將時間結合起來，以解決空氣濕度的時空分布問題。

參考文獻：
［1］ 姜曉燕，劉樹華，馬明敏，等．東北地區近百年降水時間序列變化規律的小波分析［J］．地理研究，2009，28．（2）：354－362．
［2］ 吳正方，靳英華，劉吉平，等．東北地區植被分布全球氣候變化區域響應［J］．地理科學，2003，23（5）：564－570．
［3］ 李 莎 ， 舒 紅 ， 徐正全 ． 利 用時空 Kriging 進行氣溫插值研究［J］．武漢大學學報（信息科學版），2012，37（2）：237－241．
［4］ 賀 偉，布仁倉，熊在平，等．1961－2005 年東北地區氣溫和降水變化趨勢［J］．生態學報，2013，33（2）：519－531．
［5］ 孫 力，安 剛，丁 立，等．中國東北地區夏季降水異常的氣候分析［J］．氣象學報，2000，58（1）：70－82．
［6］ 閆慶武，卞正富，王 紅．利用泰森多邊形和格網平滑的人口密度空間化研究———以徐州市為例［J］.武漢大學學報（信息科學版），2011，36（8）：987－990．
［7］ OLIVER M A，WEBSTER R． Kriging： A method of interpola-tion for geographical information systems ［J］. International Jour-nal of Geographical Information Systems，1990，4（3）：313－332．
［8］ LOUIS D E． Pollution models and inverse distance weighting：Some critical remarks ［J］. Computers ＆ Geosciences，2013，52：459－469．
［9］ 胡緒福，唐 瑋，王 博，等．趨勢面分析法應用效果研究———以YL 油田某區塊為例 ［J］. 長 江大學學報 （ 自然科學版 ），2010，7（3）：536－538．
［10］ 劉小飛，張寄陽，劉祖貴，等．噴灌大田土壤水分不同空間插值方法對比分析［J］. 灌溉排水學報，2008，27（4）：116－118．
［11］ 薛 毅，陳立萍．R 統計建模與 R 軟件［M］. 北京：清華大學出版社，2007．