Wilcoxon秩和檢驗在熱帶氣旋強度預報方法評定中的應用

周聰，余暉，傅剛

(1．中國海洋大學海洋環(huán)境學院，山東青島266100;2．中國氣象局上海臺風研究所，上海200030)

0 引言

熱帶氣旋(tropical cyclone，簡稱TC)強度變化的研究和預報歷來受到重視(于潤玲等，2013)。目前國內外常見的熱帶氣旋強度預報方法主要有數值預報方法和統(tǒng)計預報方法。研究表明，數值預報對TC路徑預報效果較好，而統(tǒng)計預報對TC強度預報的優(yōu)勢較明顯(湯杰等，2011)。在對TC強度業(yè)務預報的精度進行評定時，一般考慮TC強度的平均絕對誤差、預報趨勢一致率、均方根誤差、箱線圖和技巧水平等指標(占瑞芬等，2010;湯杰等，2011)。

在TC強度的業(yè)務預報上，除了關注上述指標外，往往還需要知道TC強度預報的顯著性檢驗結果。從統(tǒng)計的角度看，由于對TC強度預報精度評定時所考慮的樣本總體分布是非正態(tài)分布，所以在顯著性檢驗時不能直接利用常規(guī)的t、F等參數檢驗方法，而應先對樣本采取正態(tài)化轉換，或者使用簡單有效的非參數檢驗法。事實上，氣象統(tǒng)計樣本總體分布往往是非正態(tài)分布，或僅滿足局部正態(tài)分布(Gunn and Marshall，1958;White，1980;Ulbrich，1983;Birmili et al．，2001)，因此研究非參數檢驗法在氣象上的應用具有十分重要的科學意義和實用價值。氣象上常見的非參數檢驗方法是蒙特卡羅檢驗，在使用蒙特卡羅檢驗時，為使其統(tǒng)計檢驗有較高的模擬精度，樣本量N需要很大(超過十萬甚至上百萬)。本文將介紹一種不依賴樣本總體分布，且不需要滿足大樣本條件的非參數Wilcoxon檢驗法。

目前中國TC強度的業(yè)務客觀預報方法(中國氣象局，2012)有4個:氣候持續(xù)法模型(余暉等，2006)、統(tǒng)計預報模型(Chen et al．，2011)、遺傳神經網絡模型(姚才等，2007)、偏最小二乘模型(宋金杰等，2011)。為便于討論，本文選擇兩種TC強度業(yè)務預報方法進行比較，即:氣候持續(xù)模型(Climatology and Persistence Scheme，簡稱CLIPER)和偏最小二乘法(Partial Least Squares Scheme，簡稱 PLS)。CLIPER于2005年正式投入業(yè)務使用，PLS于2010年投入準業(yè)務使用。

1 熱帶氣旋強度誤差分布

將某一時次的TC強度誤差(Error，簡稱E)定義為

式中:i指第 i個時次樣本，i=1，2，…，n;n指預報時次總數;I代表臺風強度(以TC近中心最大風速表征);f、r分別代表強度預報模型和實際情形。

計算2005—2012年6月12 h預報時刻 CLIPER模型的預報誤差，記為Ei。

將 Ei(i=1，2，…，n)標準化，記為，繪制頻率分布，如圖1所示，其中x軸表示E(1)i值，y軸表示對應值出現的頻數。與標準正態(tài)分布曲線(虛線)比較，初步判斷的分布是非正態(tài)分布。進一步對E(1)i做描述性統(tǒng)計分析，發(fā)現其偏度與峰度超過正態(tài)分布臨界標準(Mardia，1970;Srivastava，1984)，可以判定該數據的總體分布是顯著非正態(tài)分布。如果對此類對象進行顯著性檢驗，常用的t檢驗、F檢驗等參數檢驗無法使用。

記 Ei(i=1，2，…，n)中位數為 θ，記θ(i=1，2，…，n)，并繪制頻率分布，如圖 2 所示，其中x軸表示值，y軸表示對應值出現的頻數。若將看作近似的對稱分布，此時可采取適用于對稱分布總體的非參數檢驗方法。

2 非參數檢驗法

參數檢驗指用某些參數(如正態(tài)分布中的平均值和方差)來描述總體的特征并對樣本所屬總體的性質做出若干假設檢驗，是近代氣候學研究中最常用的統(tǒng)計方法。由這類方法導出的結論會受到參數假設的限制。而對參數不作嚴格的假設，相對于參數檢驗而言，對資料的分布形態(tài)沒有特殊的要求或只作一點諸如對稱性之類的簡單假設，一般稱作非參數檢驗。

圖1 E的頻數分布Fig．1 Frequency distribution of E

圖2 E的頻數分布Fig．2 Frequency distribution of E

非參數檢驗把數據按大小排隊，使每個數據都有自己的“地位”，稱為秩(Rank)。雖然不知道原始數據的分布形式，但是這些秩及由其產生的一些統(tǒng)計量的性質和分布是可以得到的。吳喜之(2006)研究認為，由這種方法得出的結論更具普遍意義。

常用的非參數統(tǒng)計方法有符號檢驗、秩和檢驗、等級相關檢驗、Mann-Whitney檢驗，超大樣本如氣候學研究中常見的蒙特卡羅檢驗等。在氣象學研究中，很多情況下參數檢驗不適用，這是因為在參數檢驗中，數據的分布一般默認為正態(tài)分布、方差齊性，然而大多數氣象數據并不滿足此種特征，在這種情況下假如堅持使用參數檢驗，勢必會歪曲數據的本來面目，從而使人們對檢驗的有效性產生懷疑。

Wilcoxon(1945)針對對稱分布提出了一種簡單有效的非參數檢驗方法，在總體分布未知、或者不考慮研究對象具體分布的情況下，非參數檢驗可以用于研究目標總體與理論總體分布是否相同，或者各樣本所在總體的分布位置是否相同(張文彤和閆潔，2004)。Wilcoxon(1945)巧妙地將非參數過程轉化為參數過程(Conover and Iman，1981)，該方法在醫(yī)學統(tǒng)計、生物、行為科學中已有廣泛應用(Gehan，1965;Tarone and Ware，1977;蔣知儉，1997;Roth et al．，1998)。

3 Wilcoxon符號秩和檢驗

在 X1，X2，…，Xn中，Xi如果是第 Ri個最小值，即Xi=X(Ri)(第 Ri個順序統(tǒng)計量)，稱 Ri為 Xi的秩。R統(tǒng)計量只考慮了樣本點的大小而未考慮其絕對值的大小，然而其絕對值的大小有時也必須考慮。

Wilcoxon符號秩統(tǒng)計量用W來表示。樣本的絕對值|X1|，|X2|，…，|Xn|排序，其順序統(tǒng)計量為表示|Xj|在絕對值樣本中的秩，即，還用S(x)表示符號函數I(x＞0)，它在x＞0時為1，否則為0。Wilcoxon符號秩統(tǒng)計量定義為

它是正的樣本點按絕對值所得的秩的和(Wilcoxon，1945)。

假設F(x－θ)為連續(xù)分布的總體，通常零假設為H0:θ=0，按照以上定義，有以下定理:

如果零假設 H0:θ=0成立，則 W1，W2，…，Wn是獨立同分布的(吳喜之，2006)，其分布為P(Wi=

用正態(tài)近似求W+的分布，則

由中心極限定理，在n大時，可近似地認為

對于上面的正態(tài)近似，也可以如下表示，當n很大時，有

以上檢驗過程不僅適用于單樣本，而且可推廣到配對樣本上(Wilcoxon，1945;吳喜之，2006)。假設(X1，Y1)，(X2，Y2)，…，(Xn，Yn)代表隨機抽樣的一組配對數據。

令 Pi=Yi－ Xi，Qi=Xi－ Yi(i=1，2，…，n)，根據中位數性質，Pi、Qi的總體也具有連續(xù)對稱性，設其總體的中位數為 θP、θQ，則原假設 H0:θ1= θ2可替換為H'0:θP=θQ=0，即回歸到單樣本情形。

在TC強度預報中實現Wilcoxon秩和檢驗，首先定義誤差E1i、E2i:

其中:i指第 i個時次樣本，i=1，2，…，n;n代表獨立樣本個數;代表兩種臺風強度預報方法的預報結果;f、r分別代表強度預報模型和實際情形。

經過整理簡化，Wilcoxon符號秩檢驗在TC強度預報檢驗中的流程為:

1)令 Wi=E2i－ E1i(i=1，2，…，n)，如果 Wi是正數，記為正號;Wi為負數，記為負號。

2)計算Wi的絕對值。

3)將Wi的絕對值按升序排列，求出相應的秩。

4)分別計算正、負號Wi的秩和，分別記為W+和W－。

若正、負秩總和大致相當，可以認為兩配對樣本數據正負變化程度基本相當，分布差距較小。

兩配對樣本的Wilcoxon符號秩和檢驗，按照下面的公式計算Z統(tǒng)計量，當樣本容量足夠大時，Z近似服從正態(tài)分布。

對于單邊檢驗H0:θ=θ0，Hα:θ＞θ0及水平α，W的臨界值為

其中:Zα通過查正態(tài)分布表獲得。若W的計算值超過臨界值c，則否定原假設，認為兩個樣本的分布不同。

分別計算2005—2012年12 h預報時刻下CLIPER模型、PLS模型的預報誤差，記 ECLIPERi、EPLSi。根據Wilcoxon秩和檢驗方法，得到Z統(tǒng)計量超過臨界值，拒絕原假設，即認為2005—2012年6月的CLIPER模型和PLS模型在12 h的預報誤差存在顯著性差異。

4 結論與討論

TC強度的預報結果與觀測結果的誤差分布是非正態(tài)分布，此時若進行顯著性檢驗，可使用非參數檢驗方法。作為最常用的非參數檢驗方法之一，Wilcoxon秩和檢驗不要求確保樣本所屬的總體符合某種理論分布，僅需總體具有對稱性，具有十分廣泛的適用性(吳喜之，2006)。

Wilcoxon秩和檢驗也可以與其他TC強度預報評定參數相結合，有利于進一步提高當前精度預報水平。然而，該方法也有一定的局限性。在總體分布已知時，相比參數方法，非參數方法無法充分利用樣本中信息，因此，在大多情況下參數檢驗方法仍然是最有效的顯著性檢驗方法。

蔣知儉．1997．醫(yī)學統(tǒng)計學［M］．北京:人民衛(wèi)生出版社:168-174．

宋金杰，王元，陳佩燕，等．2011．基于偏最小二乘回歸理論的西北太平洋熱帶氣旋強度統(tǒng)計預報方法［J］．氣象學報，69(5):745-756．

湯杰，陳國民，余暉．2011.2010年西北太平洋臺風預報精度評定及分析［J］．氣象，37(10):1320-1328．

吳喜之．2006．非參數統(tǒng)計［M］．北京:中國統(tǒng)計出版社:23-25;32-41．

姚才，金龍，黃明策．2007．遺傳算法與神經網絡相結合的熱帶氣旋強度預報方法試驗［J］．海洋學報，29(4):11-19．

余暉，胡春梅，蔣樂貽．2006．熱帶氣旋強度資料的差異性分析［J］．氣象學報，64(3):357-363．

于潤玲，余暉，端義宏．2013．登陸華南熱帶氣旋強度變化與大尺度環(huán)流的關系［J］．大氣科學學報，36(5):619-625．

占瑞芬，湯杰，余暉．2010.2009年西北太平洋熱帶氣旋定位和業(yè)務預報精度評定［J］．氣象，36(10):114-121．

張文彤，閆潔．2004．SPSS統(tǒng)計分析基礎教程［M］．北京:高等教育出版社:279-301．

中國氣象局．2012．臺風業(yè)務和服務規(guī)定［S］．北京:氣象出版社．

Birmili W，Wiedensohler A，Heintzenberg J，et al．2001．Atmospheric particle number size distribution in central Europe:Statistical relations to air masses and meteorology［J］．J Geophys Res，106(D23):32005-32018．

Chen P Y，Yu H，Chan J C．2011．A western North Pacific tropical cyclone intensity prediction scheme［J］．Acta Meteor Sinica，25:611-624．

Conover W J，Iman R L．1981．Rank transformations as a bridge between parametric and nonparametric statistics［J］．Am Stat，35(3):124-129．

Gehan E A．1965．A generalized two-sample Wilcoxon test for doubly censored data［J］．Biom，52(3/4):650-653．

Gunn K L S，Marshall J S．1958．The distribution with size of aggregate snowflakes［J］．J Meteor，15(5):452-461．

Mardia K V．1970．Measures of multivariate skewness and kurtosis with applications［J］．Biom，57(3):519-530．

Roth J A，Atkinson E N，Fossella F，et al．1998．Long-term follow-up of patients enrolled in a randomized trial comparing perioperative chemotherapy and surgery with surgery alone in resectable stage III:A non-small-cell lung cancer［J］．Lung Cancer，21(1):1-6．

Srivastava M S．1984．A measure of skewness and kurtosis and a graphical method for assessing multivariate normality［J］．Stat Probabil Lett，2(5):263-267．

Tarone R E，Ware J．1977．On distribution-free tests for equality of survival distributions［J］．Biom，64(1):156-160．

Ulbrich C W．1983．Natural variations in the analytical form of the raindrop size distribution［J］．J Climate Appl Meteor，22(10):1764-1775．

White G H．1980．Skewness，kurtosis and extreme values of Northern Hemisphere geopotential heights［J］．Mon Wea Rev，108(9):1446-1455．

Wilcoxon F．1945．Individual comparisons by ranking methods［J］．Biometrics Bull，1(6):80-83．