周聰,余暉,傅剛
(1.中國海洋大學(xué)海洋環(huán)境學(xué)院,山東青島266100;2.中國氣象局上海臺風(fēng)研究所,上海200030)
熱帶氣旋(tropical cyclone,簡稱TC)強(qiáng)度變化的研究和預(yù)報(bào)歷來受到重視(于潤玲等,2013)。目前國內(nèi)外常見的熱帶氣旋強(qiáng)度預(yù)報(bào)方法主要有數(shù)值預(yù)報(bào)方法和統(tǒng)計(jì)預(yù)報(bào)方法。研究表明,數(shù)值預(yù)報(bào)對TC路徑預(yù)報(bào)效果較好,而統(tǒng)計(jì)預(yù)報(bào)對TC強(qiáng)度預(yù)報(bào)的優(yōu)勢較明顯(湯杰等,2011)。在對TC強(qiáng)度業(yè)務(wù)預(yù)報(bào)的精度進(jìn)行評定時,一般考慮TC強(qiáng)度的平均絕對誤差、預(yù)報(bào)趨勢一致率、均方根誤差、箱線圖和技巧水平等指標(biāo)(占瑞芬等,2010;湯杰等,2011)。
在TC強(qiáng)度的業(yè)務(wù)預(yù)報(bào)上,除了關(guān)注上述指標(biāo)外,往往還需要知道TC強(qiáng)度預(yù)報(bào)的顯著性檢驗(yàn)結(jié)果。從統(tǒng)計(jì)的角度看,由于對TC強(qiáng)度預(yù)報(bào)精度評定時所考慮的樣本總體分布是非正態(tài)分布,所以在顯著性檢驗(yàn)時不能直接利用常規(guī)的t、F等參數(shù)檢驗(yàn)方法,而應(yīng)先對樣本采取正態(tài)化轉(zhuǎn)換,或者使用簡單有效的非參數(shù)檢驗(yàn)法。事實(shí)上,氣象統(tǒng)計(jì)樣本總體分布往往是非正態(tài)分布,或僅滿足局部正態(tài)分布(Gunn and Marshall,1958;White,1980;Ulbrich,1983;Birmili et al.,2001),因此研究非參數(shù)檢驗(yàn)法在氣象上的應(yīng)用具有十分重要的科學(xué)意義和實(shí)用價(jià)值。氣象上常見的非參數(shù)檢驗(yàn)方法是蒙特卡羅檢驗(yàn),在使用蒙特卡羅檢驗(yàn)時,為使其統(tǒng)計(jì)檢驗(yàn)有較高的模擬精度,樣本量N需要很大(超過十萬甚至上百萬)。本文將介紹一種不依賴樣本總體分布,且不需要滿足大樣本條件的非參數(shù)Wilcoxon檢驗(yàn)法。
目前中國TC強(qiáng)度的業(yè)務(wù)客觀預(yù)報(bào)方法(中國氣象局,2012)有4個:氣候持續(xù)法模型(余暉等,2006)、統(tǒng)計(jì)預(yù)報(bào)模型(Chen et al.,2011)、遺傳神經(jīng)網(wǎng)絡(luò)模型(姚才等,2007)、偏最小二乘模型(宋金杰等,2011)。為便于討論,本文選擇兩種TC強(qiáng)度業(yè)務(wù)預(yù)報(bào)方法進(jìn)行比較,即:氣候持續(xù)模型(Climatology and Persistence Scheme,簡稱CLIPER)和偏最小二乘法(Partial Least Squares Scheme,簡稱 PLS)。CLIPER于2005年正式投入業(yè)務(wù)使用,PLS于2010年投入準(zhǔn)業(yè)務(wù)使用。
將某一時次的TC強(qiáng)度誤差(Error,簡稱E)定義為
式中:i指第 i個時次樣本,i=1,2,…,n;n指預(yù)報(bào)時次總數(shù);I代表臺風(fēng)強(qiáng)度(以TC近中心最大風(fēng)速表征);f、r分別代表強(qiáng)度預(yù)報(bào)模型和實(shí)際情形。
計(jì)算2005—2012年6月12 h預(yù)報(bào)時刻 CLIPER模型的預(yù)報(bào)誤差,記為Ei。
將 Ei(i=1,2,…,n)標(biāo)準(zhǔn)化,記為,繪制頻率分布,如圖1所示,其中x軸表示E(1)i值,y軸表示對應(yīng)值出現(xiàn)的頻數(shù)。與標(biāo)準(zhǔn)正態(tài)分布曲線(虛線)比較,初步判斷的分布是非正態(tài)分布。進(jìn)一步對E(1)i做描述性統(tǒng)計(jì)分析,發(fā)現(xiàn)其偏度與峰度超過正態(tài)分布臨界標(biāo)準(zhǔn)(Mardia,1970;Srivastava,1984),可以判定該數(shù)據(jù)的總體分布是顯著非正態(tài)分布。如果對此類對象進(jìn)行顯著性檢驗(yàn),常用的t檢驗(yàn)、F檢驗(yàn)等參數(shù)檢驗(yàn)無法使用。
記 Ei(i=1,2,…,n)中位數(shù)為 θ,記θ(i=1,2,…,n),并繪制頻率分布,如圖 2 所示,其中x軸表示值,y軸表示對應(yīng)值出現(xiàn)的頻數(shù)。若將看作近似的對稱分布,此時可采取適用于對稱分布總體的非參數(shù)檢驗(yàn)方法。
參數(shù)檢驗(yàn)指用某些參數(shù)(如正態(tài)分布中的平均值和方差)來描述總體的特征并對樣本所屬總體的性質(zhì)做出若干假設(shè)檢驗(yàn),是近代氣候?qū)W研究中最常用的統(tǒng)計(jì)方法。由這類方法導(dǎo)出的結(jié)論會受到參數(shù)假設(shè)的限制。而對參數(shù)不作嚴(yán)格的假設(shè),相對于參數(shù)檢驗(yàn)而言,對資料的分布形態(tài)沒有特殊的要求或只作一點(diǎn)諸如對稱性之類的簡單假設(shè),一般稱作非參數(shù)檢驗(yàn)。
圖1 E的頻數(shù)分布Fig.1 Frequency distribution of E
圖2 E的頻數(shù)分布Fig.2 Frequency distribution of E
非參數(shù)檢驗(yàn)把數(shù)據(jù)按大小排隊(duì),使每個數(shù)據(jù)都有自己的“地位”,稱為秩(Rank)。雖然不知道原始數(shù)據(jù)的分布形式,但是這些秩及由其產(chǎn)生的一些統(tǒng)計(jì)量的性質(zhì)和分布是可以得到的。吳喜之(2006)研究認(rèn)為,由這種方法得出的結(jié)論更具普遍意義。
常用的非參數(shù)統(tǒng)計(jì)方法有符號檢驗(yàn)、秩和檢驗(yàn)、等級相關(guān)檢驗(yàn)、Mann-Whitney檢驗(yàn),超大樣本如氣候?qū)W研究中常見的蒙特卡羅檢驗(yàn)等。在氣象學(xué)研究中,很多情況下參數(shù)檢驗(yàn)不適用,這是因?yàn)樵趨?shù)檢驗(yàn)中,數(shù)據(jù)的分布一般默認(rèn)為正態(tài)分布、方差齊性,然而大多數(shù)氣象數(shù)據(jù)并不滿足此種特征,在這種情況下假如堅(jiān)持使用參數(shù)檢驗(yàn),勢必會歪曲數(shù)據(jù)的本來面目,從而使人們對檢驗(yàn)的有效性產(chǎn)生懷疑。
Wilcoxon(1945)針對對稱分布提出了一種簡單有效的非參數(shù)檢驗(yàn)方法,在總體分布未知、或者不考慮研究對象具體分布的情況下,非參數(shù)檢驗(yàn)可以用于研究目標(biāo)總體與理論總體分布是否相同,或者各樣本所在總體的分布位置是否相同(張文彤和閆潔,2004)。Wilcoxon(1945)巧妙地將非參數(shù)過程轉(zhuǎn)化為參數(shù)過程(Conover and Iman,1981),該方法在醫(yī)學(xué)統(tǒng)計(jì)、生物、行為科學(xué)中已有廣泛應(yīng)用(Gehan,1965;Tarone and Ware,1977;蔣知儉,1997;Roth et al.,1998)。
在 X1,X2,…,Xn中,Xi如果是第 Ri個最小值,即Xi=X(Ri)(第 Ri個順序統(tǒng)計(jì)量),稱 Ri為 Xi的秩。R統(tǒng)計(jì)量只考慮了樣本點(diǎn)的大小而未考慮其絕對值的大小,然而其絕對值的大小有時也必須考慮。
Wilcoxon符號秩統(tǒng)計(jì)量用W來表示。樣本的絕對值|X1|,|X2|,…,|Xn|排序,其順序統(tǒng)計(jì)量為表示|Xj|在絕對值樣本中的秩,即,還用S(x)表示符號函數(shù)I(x>0),它在x>0時為1,否則為0。Wilcoxon符號秩統(tǒng)計(jì)量定義為
它是正的樣本點(diǎn)按絕對值所得的秩的和(Wilcoxon,1945)。
假設(shè)F(x-θ)為連續(xù)分布的總體,通常零假設(shè)為H0:θ=0,按照以上定義,有以下定理:
如果零假設(shè) H0:θ=0成立,則 W1,W2,…,Wn是獨(dú)立同分布的(吳喜之,2006),其分布為P(Wi=
用正態(tài)近似求W+的分布,則
由中心極限定理,在n大時,可近似地認(rèn)為
對于上面的正態(tài)近似,也可以如下表示,當(dāng)n很大時,有
以上檢驗(yàn)過程不僅適用于單樣本,而且可推廣到配對樣本上(Wilcoxon,1945;吳喜之,2006)。假設(shè)(X1,Y1),(X2,Y2),…,(Xn,Yn)代表隨機(jī)抽樣的一組配對數(shù)據(jù)。
令 Pi=Yi- Xi,Qi=Xi- Yi(i=1,2,…,n),根據(jù)中位數(shù)性質(zhì),Pi、Qi的總體也具有連續(xù)對稱性,設(shè)其總體的中位數(shù)為 θP、θQ,則原假設(shè) H0:θ1= θ2可替換為H'0:θP=θQ=0,即回歸到單樣本情形。
在TC強(qiáng)度預(yù)報(bào)中實(shí)現(xiàn)Wilcoxon秩和檢驗(yàn),首先定義誤差E1i、E2i:
其中:i指第 i個時次樣本,i=1,2,…,n;n代表獨(dú)立樣本個數(shù);代表兩種臺風(fēng)強(qiáng)度預(yù)報(bào)方法的預(yù)報(bào)結(jié)果;f、r分別代表強(qiáng)度預(yù)報(bào)模型和實(shí)際情形。
經(jīng)過整理簡化,Wilcoxon符號秩檢驗(yàn)在TC強(qiáng)度預(yù)報(bào)檢驗(yàn)中的流程為:
1)令 Wi=E2i- E1i(i=1,2,…,n),如果 Wi是正數(shù),記為正號;Wi為負(fù)數(shù),記為負(fù)號。
2)計(jì)算Wi的絕對值。
3)將Wi的絕對值按升序排列,求出相應(yīng)的秩。
4)分別計(jì)算正、負(fù)號Wi的秩和,分別記為W+和W-。
若正、負(fù)秩總和大致相當(dāng),可以認(rèn)為兩配對樣本數(shù)據(jù)正負(fù)變化程度基本相當(dāng),分布差距較小。
兩配對樣本的Wilcoxon符號秩和檢驗(yàn),按照下面的公式計(jì)算Z統(tǒng)計(jì)量,當(dāng)樣本容量足夠大時,Z近似服從正態(tài)分布。
對于單邊檢驗(yàn)H0:θ=θ0,Hα:θ>θ0及水平α,W的臨界值為
其中:Zα通過查正態(tài)分布表獲得。若W的計(jì)算值超過臨界值c,則否定原假設(shè),認(rèn)為兩個樣本的分布不同。
分別計(jì)算2005—2012年12 h預(yù)報(bào)時刻下CLIPER模型、PLS模型的預(yù)報(bào)誤差,記 ECLIPERi、EPLSi。根據(jù)Wilcoxon秩和檢驗(yàn)方法,得到Z統(tǒng)計(jì)量超過臨界值,拒絕原假設(shè),即認(rèn)為2005—2012年6月的CLIPER模型和PLS模型在12 h的預(yù)報(bào)誤差存在顯著性差異。
TC強(qiáng)度的預(yù)報(bào)結(jié)果與觀測結(jié)果的誤差分布是非正態(tài)分布,此時若進(jìn)行顯著性檢驗(yàn),可使用非參數(shù)檢驗(yàn)方法。作為最常用的非參數(shù)檢驗(yàn)方法之一,Wilcoxon秩和檢驗(yàn)不要求確保樣本所屬的總體符合某種理論分布,僅需總體具有對稱性,具有十分廣泛的適用性(吳喜之,2006)。
Wilcoxon秩和檢驗(yàn)也可以與其他TC強(qiáng)度預(yù)報(bào)評定參數(shù)相結(jié)合,有利于進(jìn)一步提高當(dāng)前精度預(yù)報(bào)水平。然而,該方法也有一定的局限性。在總體分布已知時,相比參數(shù)方法,非參數(shù)方法無法充分利用樣本中信息,因此,在大多情況下參數(shù)檢驗(yàn)方法仍然是最有效的顯著性檢驗(yàn)方法。
蔣知儉.1997.醫(yī)學(xué)統(tǒng)計(jì)學(xué)[M].北京:人民衛(wèi)生出版社:168-174.
宋金杰,王元,陳佩燕,等.2011.基于偏最小二乘回歸理論的西北太平洋熱帶氣旋強(qiáng)度統(tǒng)計(jì)預(yù)報(bào)方法[J].氣象學(xué)報(bào),69(5):745-756.
湯杰,陳國民,余暉.2011.2010年西北太平洋臺風(fēng)預(yù)報(bào)精度評定及分析[J].氣象,37(10):1320-1328.
吳喜之.2006.非參數(shù)統(tǒng)計(jì)[M].北京:中國統(tǒng)計(jì)出版社:23-25;32-41.
姚才,金龍,黃明策.2007.遺傳算法與神經(jīng)網(wǎng)絡(luò)相結(jié)合的熱帶氣旋強(qiáng)度預(yù)報(bào)方法試驗(yàn)[J].海洋學(xué)報(bào),29(4):11-19.
余暉,胡春梅,蔣樂貽.2006.熱帶氣旋強(qiáng)度資料的差異性分析[J].氣象學(xué)報(bào),64(3):357-363.
于潤玲,余暉,端義宏.2013.登陸華南熱帶氣旋強(qiáng)度變化與大尺度環(huán)流的關(guān)系[J].大氣科學(xué)學(xué)報(bào),36(5):619-625.
占瑞芬,湯杰,余暉.2010.2009年西北太平洋熱帶氣旋定位和業(yè)務(wù)預(yù)報(bào)精度評定[J].氣象,36(10):114-121.
張文彤,閆潔.2004.SPSS統(tǒng)計(jì)分析基礎(chǔ)教程[M].北京:高等教育出版社:279-301.
中國氣象局.2012.臺風(fēng)業(yè)務(wù)和服務(wù)規(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,F(xiàn)ossella 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.