中國東北地區(qū)氣象數(shù)據(jù)的空間平穩(wěn)性檢驗

許家琦，李顏伶，舒 紅＊

（1.武漢大學(xué)測繪遙感信息工程國家重點實驗室，武漢430079；2.青島勘察測繪研究院，山東青島266032）

在對空間數(shù)據(jù)的分析中總是假設(shè)空間協(xié)方差結(jié)構(gòu)平穩(wěn)，然而這一假設(shè)在實際中常常不成立.例如降水數(shù)據(jù)就顯示出空間非平穩(wěn)性［1－2］，月降水中硫酸濃度也呈現(xiàn)出了非平穩(wěn)結(jié)構(gòu)［3］.空間數(shù)據(jù)的平穩(wěn)性檢驗有多種方法.地理加權(quán)回歸法將空間的概念引入到一般線性回歸的回歸參數(shù)中，回歸系數(shù)隨著空間位置的變化而變化［4］.地理加權(quán)回歸法檢驗的是變量之間相關(guān)關(guān)系的空間異質(zhì)性，地理加權(quán)回歸法（GWR）在已知數(shù)據(jù)非平穩(wěn)的情況下分析這種非平穩(wěn)性，而不是用于數(shù)據(jù)檢驗是否平穩(wěn).本文采用的3種方法（三維散點圖法，趨勢線法和鄰域圓法）都是用來檢驗屬性數(shù)據(jù)本身的空間異質(zhì)性.譜方法用局部平穩(wěn)過程的卷積來表示整個非平穩(wěn)過程，并用這個非平穩(wěn)模型來進行譜密度的參數(shù)和非參數(shù)估計［5－6］.由于每個被視為平穩(wěn)的局部區(qū)域中都要有足夠多的數(shù)據(jù)來計算空間變異函數(shù)然或者空間協(xié)方差函數(shù)［7］，因此這種方法并不適用于我國較為稀疏的氣象站點數(shù)據(jù).三維散點圖法擬合的是整個研究區(qū)域內(nèi)屬性數(shù)據(jù)的趨勢面；趨勢線法擬合的是投影面上所有屬性數(shù)據(jù)的趨勢；雖然鄰域圓法也是把空間分割為一些小的區(qū)域來研究，但鄰域圓法關(guān)注的是局部區(qū)域內(nèi)的二階矩.均值和方差的計算與經(jīng)驗變異函數(shù)的計算不同，對于參與計算的氣象站點沒有太多數(shù)量上的要求.在非平穩(wěn)的空間協(xié)方差結(jié)構(gòu)分析中，把不平穩(wěn)各向異性的地理空間映射到平穩(wěn)各向同性的頻譜空間中進行處理，然后將處理結(jié)果逆變換到地理空間中［8］.Sampson＆Guttorp法（SG模型）計算復(fù)雜，物理含義不夠直觀明確.而本文采用3種方法：（1）三維散點圖法有著簡單直觀的優(yōu)點，通過查看氣象數(shù)據(jù)在三維空間中的數(shù)值大小來初步判斷數(shù)據(jù)是否平穩(wěn)；（2）在趨勢線法中，把氣象數(shù)據(jù)投影到東西和南北兩個面上，在這兩個面上分別用多項式法擬合出兩條趨勢線，從而判斷出空間數(shù)據(jù)在東西方向和南北方向上是否平穩(wěn)［9］.（3）由于隨機變量的一階矩和二階矩已經(jīng)足夠為大多數(shù)實際問題提供可以接受的近似解決方法［10］，在鄰域圓法中通過統(tǒng)計圓內(nèi)樣本數(shù)據(jù)的均值和方差來檢驗我國東北地區(qū)氣象數(shù)據(jù)的空間平穩(wěn)性［11］.

1 研究區(qū)域概況

本文采用我國東北地區(qū)的氣象觀測站點數(shù)據(jù)，包括2006年1月1日的日平均氣溫數(shù)據(jù)，2006年1月1日的日平均風速數(shù)據(jù)以及2006年5月的月降水數(shù)據(jù)作為實驗數(shù)據(jù).氣象站點分布如圖1所示.

冬季，受極地大陸氣團控制，東北區(qū)氣候嚴寒.1月是全年最冷月，全區(qū)的月平均氣溫均在－4℃以下.溫度梯度之大，為全年之冠.大小興安嶺北段的平均氣溫在－30℃左右，是我國冬季最冷的地方.冬季（以1月為代表）在蒙古冷高壓影響下，空氣穩(wěn)定，風速較小.特別是松花江流域以北地區(qū)，因位置最北，常在高壓內(nèi)部控制之下，冬季為全年風速最小季節(jié).1月份全區(qū)平均風速在1～5 m／s之間；山區(qū)13 m／s；平原多在3～4 m／s.本區(qū)的降水主要集中于夏季.各地降水一般從5月份明顯增加，9月份又陡降.降水量自東南向西北有遞減的趨勢.降水日多于20 d的區(qū)域由東北部的山區(qū)延伸到了松遼平原的東部，降水日數(shù)的高值區(qū)仍是小興安嶺南部和長白山一帶［12］.

圖1 中國東北地區(qū)氣象觀測站點分布圖Fig.1 Location map of meteorological stations in the Northeast area of China

2 研究方法

設(shè)一隨機函數(shù)Z，其空間分布律不因平移而改變，即若對任一向量h，兩個k維向量的隨機變量｛Z（x1），Z（x2），…，Z（xk）｝和｛Z（x1＋h），Z（x2＋h），…，Z（xk＋h）｝.

有相同的分布律.這種平穩(wěn)假設(shè)至少要求Z（x）的各階矩均存在，且平穩(wěn)，而在實際工作中卻很難滿足.在線性地質(zhì)統(tǒng)計學(xué)研究中，只需假設(shè)其1、2階矩存在且平穩(wěn)就夠了，因而提出二階平穩(wěn)假設(shè).即

假定：1）在整個研究區(qū)內(nèi)，區(qū)域化變量Z（x）的期望存在且等于常數(shù)：

2）在整個研究區(qū)內(nèi)，區(qū)域化變量的空間協(xié)方差函數(shù)存在且平穩(wěn)［13］：

2.1 三維散點圖法

空間現(xiàn)象可以分解為兩部分：區(qū)域性趨勢（平均值）和殘差（變化性）.二階平穩(wěn)意味著區(qū)域范圍內(nèi)隨機變量的平均值和方差均沒有趨勢變化.平穩(wěn)性就是觀測值圍繞固定常數(shù)上下波動，且離開均值上下波動的幅度大體相似（方差為常量）.反之，如果趨勢面不是一個水平面或者數(shù)據(jù)圍繞趨勢面波動的幅度不均衡就是非平穩(wěn).所以可以對觀測值作散點圖，同時顯示趨勢面，通過人眼初步觀察一下是否滿足平穩(wěn)性.

2.2 趨勢線法

將三維散點圖分別向東西方向和南北方向上投影，繪制氣象數(shù)據(jù)與x軸，y軸的二維散點圖.用一條紅色的趨勢線（既可以是直線也可以使曲線）來擬合這些數(shù)據(jù)，通過兩個方向上的剖面圖可以清楚看到數(shù)據(jù)趨勢在哪個方向上.

2.3 鄰域圓法

所謂鄰域圓法就是在一個觀測點的周圍以該點為圓心畫一個圓，根據(jù)圓內(nèi)的數(shù)據(jù)計算均值和相對標準差，然后對每個觀測點都做同樣的計算，通過圓內(nèi)均值和相對標準差在上是否平穩(wěn).由于這種方法可以根據(jù)氣象站點之間的距離來計算出合適的鄰域圓半徑并且可以準確地遍歷整個研究區(qū)域，比起移動窗口法［11］，這種方法更適合我國稀疏的氣象站點數(shù)據(jù).

標準差描述的是區(qū)域化變量觀測值間的離散情況，因此計算其標準差就相當于求得離散情況，但是由于均值不同時變異程度就不同，即區(qū)域內(nèi)的標準差依賴于區(qū)域內(nèi)的均值，均值不同不能直接采用標準差，因此這里采用相對標準差，也稱變異系數(shù)，是一個無量綱的統(tǒng)計量，能反映單位均值上的離散程度，可以消除單位或平均數(shù)不同對兩個或多個資料變異程度比較的影響［14］.用公式表示為：

式中s表示圓內(nèi)樣點標準差，m表示圓內(nèi)樣點均值.

如果數(shù)據(jù)平穩(wěn)，圓內(nèi)的均值應(yīng)該圍繞某個固定常數(shù)上下隨機波動，沒有明顯的趨勢，相對標準差應(yīng)該基本在0值附近波動.

3 結(jié)果分析

本文用R語言實現(xiàn)了對我國東北地區(qū)氣象要素（氣溫，風速，降水）的空間平穩(wěn)性檢驗.R語言有著優(yōu)秀的統(tǒng)計分析能力和繪圖環(huán)境.R里面自帶的包中包含了統(tǒng)計學(xué)家們開發(fā)的前沿算法.通過直接調(diào)用這些包里面的函數(shù)（如persp3d（））可以避免通過大量編程實現(xiàn)統(tǒng)計算法.

3.1 三維散點圖法

利用rgl包中的persp3d（）函數(shù)繪制三維散點圖［15］，如圖2所示.

圖2 中國東北地區(qū)氣象觀測數(shù)據(jù)散點圖（氣溫，風速，降水）Fig.2 Scatterplots of meteorological data in the Northeast area of China（temperature，wind speed and precipitation）

從圖2中可以看出，日平均氣溫有一個由北向南線性增加的趨勢，且增幅明顯.東北地區(qū)2006年1月1日的日平均氣溫從高緯地區(qū)的－25℃到低緯地區(qū)的－5℃.但是氣溫變量的隨機波動很小，都是圍繞趨勢面小幅度的上下波動.風速有一個自西向東減少的趨勢，趨勢不明顯.但是風速變量的隨機波動很大，即數(shù)據(jù)偏離趨勢面的幅度很大.降水數(shù)據(jù)的趨勢面是個明顯的曲面，趨勢面不平直，變量的隨機性一般.

3.2 趨勢線法

利用r中的l m（）和predict（）函數(shù)繪制趨勢線，如圖3所示.

從圖3中可以看出，氣溫數(shù)據(jù)有一個極為明顯的由北向南增加的趨勢，由西向東方向則微弱減少.風速數(shù)據(jù)基本上可以看成是一個平穩(wěn)的數(shù)據(jù)集，不論是在南北還是東西方向上都圍繞平均風速2 m／s上下隨機波動.而5月的降水數(shù)據(jù)由西向東增加，由南向北減少，有著明顯的趨勢且隨機波動較大.

3.3 鄰域圓檢驗

圖3 東北地區(qū)氣象觀測數(shù)據(jù)趨勢線圖（氣溫，風速，降水）Fig.3 Trend line maps of meteorological data in the Northeast area of China（temperature，wind speed and precipitation）

用R實現(xiàn)鄰域圓法的思路為：定義圓半徑值的大小，計算第一個點與所有點（包括第一個點本身）之間的距離，存儲距離小于半徑值的點并計算其均值和相對標準差，得到第一個點鄰域圓內(nèi)的均值和相對標準差.同理可計算得所有點鄰域圓內(nèi)的均值和相對標準差，并作圖顯示以判斷該區(qū)域化變量在該空間上是否平穩(wěn).鄰域圓按照從西向東、從北向南的順序遍歷整個研究區(qū)域.其結(jié)果如圖4～圖6所示.

圖4 東北地區(qū)日平均氣溫數(shù)據(jù)鄰域圓法檢驗圖Fig.4 Neighborhood circle plots of daily average temperature data in the Northeast area of China

圖5 東北地區(qū)日平均風速數(shù)據(jù)鄰域圓法檢驗圖Fig.5 Neighborhood circle plots of daily average wind speed data in the Northeast area of China

圖6 東北地區(qū)月降水數(shù)據(jù)鄰域圓法檢驗圖Fig.6 Neighborhood circle plots of monthly precipitation data in the Northeast area of China

從圖4～圖6中可以看到2006年1月1日的風速數(shù)據(jù)可以看成是平穩(wěn)的數(shù)據(jù)，而氣溫數(shù)據(jù)和降水數(shù)據(jù)均不平穩(wěn).鄰域圓內(nèi)的氣溫均值逐步上升，由于氣溫圍繞均值波動的幅度始終較小，因此相對標準差隨著均值的增加而逐步減小.鄰域圓內(nèi)的降水均值呈螺旋式上升的趨勢，相對標準差從樣本點號1到7是平穩(wěn)的，從9號樣本點到81也是平穩(wěn)的.但是從樣本號7到9之間降水數(shù)據(jù)的相對標準差急劇下降.樣本點號為7，8，9的氣象站如圖所示，由一個三角形圍住.從圖中可以看出，由于1～7號樣本點處于東北地區(qū)的西北角上，分布稀疏且呈一個狹長的形狀，因此落在每個鄰域圓的樣本數(shù)相對較少，標準差較大.而1～7號樣本點周圍鄰域圓內(nèi)的樣本均值卻比較小.所以1～7號樣本點周圍鄰域圓內(nèi)的相對標準差較大.

4 結(jié)論與展望

本文介紹的3種平穩(wěn)性檢驗方法各有優(yōu)缺點，三維散點圖法簡單直觀，可以迅速看出數(shù)據(jù)的均值是否是一個水平的平面.如果是一個平直的面那么數(shù)據(jù)均值平穩(wěn)，反之不平穩(wěn).三維散點的缺點在于不便于發(fā)現(xiàn)均值不平穩(wěn)的方向性，以及二階矩是否平穩(wěn).趨勢線法把空間數(shù)據(jù)投影到東西和南北兩個方向上，通過兩個方向上的剖面圖可以清晰地看出空間不平穩(wěn)的方向性.鄰域圓法同時顯示了局部區(qū)域內(nèi)樣本點的均值圖和相對標準差圖，它可以同時清楚地展現(xiàn)一階矩和二階矩是否平穩(wěn).3種方法同時使用，對同樣3組數(shù)據(jù)（我國東北地區(qū)日氣溫，日風速，月降水）進行空間平穩(wěn)性分析，得到相同的結(jié)論.即日氣溫數(shù)據(jù)有一個由北向南迅速遞增的線性趨勢，這與我國大部分地區(qū)處于北半球的溫帶，氣溫南高北低的氣候規(guī)律相吻合.東北地區(qū)日風速數(shù)據(jù)在空間上比較平穩(wěn)，始終在2 m／s的平均值上下波動，但是隨機波動幅度較大，從0.5 m／s到4 m／s不等.月降水量有著東高西低，南高北低的趨勢.我國東北地區(qū)南面臨近渤海、黃海，東面臨近日本海，臨海地區(qū)空氣濕潤，雨量相對充沛.降水的空間趨勢檢驗結(jié)果也與我國東北地區(qū)的地理環(huán)境相吻合.

本文有效檢驗了中國東北地區(qū)氣象要素的空間平穩(wěn)性.在統(tǒng)計建模中如何處理不平穩(wěn)的空間數(shù)據(jù)是一個有待進一步研究的問題.

［1］ Berndtsson R.Temporal variability in spatial correlation of daily rainfal［J］.Water Resources Research，1988（24）：1511－1517.

［2］ Obled C，Creutin J D.Some developments in the use of empirical orthogonal functions for mapping meteorological Field［J］.Journal of Climate and Applied Meteorolog，1986（25）：1189－1204.

［3］ Loader C，Switzer P.Spatial covariance estimation for monitoring data［R］.Stanford：Department of statistics，Stanford University，1989.

［4］ Brunsdon C，F(xiàn)otheringham S，Charlton M.Geographically weighted regression－modeling spatial non－stationarity［J］.The Statistician，1998（47）：431－443.

［5］ Fuentes M.Spectral method for nonstationary spatial processes［J］.Biometrika，2002，1（89）：197－210.

［6］ Fuentes M.A formal test for nonstationarity of spatial stochastic processes［J］.Journal of Multivariate Analysis，2005（96）：30－54.

［7］ Goovaerts P.Geostatistics for natural resources evaluation［M］.New York，Oxford：Oxford University Press，1997.

［8］ Sampson P D，Guttorp P.Nonparametric estimation of nonstationary spatial covariance structure［J］.Journal of the A－merican Statistical Association，1992，87（417）：108－119.

［9］ 湯國安，楊 昕.ArcGIS地理信息系統(tǒng)空間分析實驗教程［M］.北京：科學(xué)出版社，2006.

［10］ Journel A G，Huijbregts C J.Mining Geostatistics［M］.London，New York，San Francisco：Academic Press，1978：30－39.

［11］ Hass T C.Lognormal and moving window methods of estimating acid deposition［J］.Journal of the American Statistical Association，1990，85（412）：950－963.

［12］ 周 琳.東北氣候［M］.北京：氣象出版社，1991：25－28.

［13］ 嚴 冰，陽正熙，吳塹虹，等.地學(xué)數(shù)據(jù)分析教程［M］.北京：科學(xué)出版社，2008：91－96.

［14］ 常文淵，戴新剛，陳洪武.地質(zhì)統(tǒng)計學(xué)在氣象要素場插值的實例研究［J］.地球物理學(xué)報.2004，6（47）：982－989.

［15］ 湯銀才.R語言與統(tǒng)計分析［M］.北京：高等教育出版社，2008.