Monte-Carlo模擬與經(jīng)驗(yàn)路徑模型預(yù)測(cè)臺(tái)風(fēng)極值風(fēng)速的對(duì)比

郭云霞，侯一筠，齊鵬

( 1. 中國(guó)科學(xué)院海洋研究所，山東 青島 266071；2. 中國(guó)科學(xué)院大學(xué)，北京 100049；3. 青島海洋科學(xué)與技術(shù)試點(diǎn)國(guó)家實(shí)驗(yàn)室海洋與氣候動(dòng)力學(xué)實(shí)驗(yàn)室，山東 青島 266237；4. 中國(guó)科學(xué)院海洋大科學(xué)研究中心，山東 青島 266071)

1 引言

中國(guó)東南沿海區(qū)域是世界上遭受臺(tái)風(fēng)災(zāi)害嚴(yán)重的地區(qū)之一。臺(tái)風(fēng)，在東太平洋和大西洋被稱為颶風(fēng)，常常引起狂風(fēng)、暴雨、巨浪和風(fēng)暴潮等自然災(zāi)害。利用臺(tái)風(fēng)風(fēng)場(chǎng)模型與數(shù)值模擬方法分析臺(tái)風(fēng)災(zāi)害風(fēng)險(xiǎn)具有十分重要的意義。

20世紀(jì)下半葉，為了解決臺(tái)風(fēng)風(fēng)速觀測(cè)數(shù)據(jù)匱乏的問(wèn)題，Monte-Carlo模擬作為可替代以及可行的方法被廣泛用來(lái)分析臺(tái)風(fēng)的危險(xiǎn)性。它利用臺(tái)風(fēng)風(fēng)場(chǎng)模型和臺(tái)風(fēng)關(guān)鍵參數(shù)資料，采用數(shù)值模擬方法構(gòu)造大量虛擬臺(tái)風(fēng)，進(jìn)而對(duì)臺(tái)風(fēng)影響區(qū)域進(jìn)行極值風(fēng)速的預(yù)測(cè)（又稱為臺(tái)風(fēng)危險(xiǎn)性分析）。采用Monte-Carlo隨機(jī)抽樣構(gòu)造虛擬臺(tái)風(fēng)的具體過(guò)程是：首先給定特定地點(diǎn)臺(tái)風(fēng)關(guān)鍵參數(shù)的統(tǒng)計(jì)概率模型，包括臺(tái)風(fēng)中心氣壓差、最大風(fēng)速半徑、移動(dòng)速度、移動(dòng)方向以及臺(tái)風(fēng)路徑與研究點(diǎn)的最小距離。給定這些參數(shù)的概率模型之后，用Monte-Carlo方法對(duì)這些模型進(jìn)行隨機(jī)抽樣，臺(tái)風(fēng)在滿足抽樣出來(lái)的關(guān)鍵參數(shù)的同時(shí)沿著直線路徑發(fā)展，同時(shí)模擬的風(fēng)速被記錄下來(lái)。在該過(guò)程中臺(tái)風(fēng)的強(qiáng)度保持不變直到登陸，當(dāng)臺(tái)風(fēng)登陸后則采用衰減模型描述中心壓差的變化，由此即可得到虛擬臺(tái)風(fēng)路徑[1]。

美國(guó)[2]和澳大利亞[3]均采用Monte-Carlo模擬方法編制設(shè)計(jì)風(fēng)速圖。該模擬方法最初是由Russell[4-5]提出并應(yīng)用于德克薩斯州海岸颶風(fēng)風(fēng)速的預(yù)測(cè)。自此以后，許多學(xué)者都對(duì)該方法進(jìn)行了擴(kuò)展與改進(jìn)。Batts等[6]和Shapiro[7]相繼提出了Batts颶風(fēng)風(fēng)場(chǎng)模型以及Shapiro模型。Georgiou[8]通過(guò)對(duì)臺(tái)風(fēng)歷史數(shù)據(jù)進(jìn)行統(tǒng)計(jì)分析，結(jié)合Georgiou臺(tái)風(fēng)模型，發(fā)展了基于Monte-Carlo方法的臺(tái)風(fēng)數(shù)值模擬方法。Vickery和Twisdale[9]研究表明Shapiro風(fēng)場(chǎng)模型要優(yōu)于Batts風(fēng)場(chǎng)模型。Meng等[10]研究了臺(tái)風(fēng)邊界層的風(fēng)場(chǎng)，建立了Yan Meng（YM）臺(tái)風(fēng)風(fēng)場(chǎng)模型。Simiu和Scanlan[11]回顧了在結(jié)構(gòu)設(shè)計(jì)中描述風(fēng)氣候的問(wèn)題，并對(duì)設(shè)計(jì)風(fēng)速的定義提出了相應(yīng)的標(biāo)準(zhǔn)。Thompson和Cardone[12]采用數(shù)值后報(bào)方法反演墨西哥灣歷史上發(fā)生的臺(tái)風(fēng)個(gè)例，從而獲得相應(yīng)臺(tái)風(fēng)引起的最大風(fēng)速序列。Vickery和Twisdale[13]指出雖然這些研究使用的方法是相似的，但在風(fēng)場(chǎng)模型、衰減模型、研究區(qū)域大小以及歷史臺(tái)風(fēng)提取方法等方面存在顯著差異。

自2000年以來(lái)，多位學(xué)者[14-20]逐步發(fā)展了全路徑模型方法。Vickery等[14]首次開(kāi)發(fā)了經(jīng)驗(yàn)路徑模型，該模型可以生成臺(tái)風(fēng)從產(chǎn)生到消亡的完整軌跡，允許臺(tái)風(fēng)在發(fā)展過(guò)程中改變移動(dòng)方向、移動(dòng)速度和強(qiáng)度。Vickery等[14]利用經(jīng)驗(yàn)路徑模型構(gòu)造虛擬臺(tái)風(fēng)的過(guò)程是：根據(jù)歷史臺(tái)風(fēng)資料，抽樣產(chǎn)生虛擬臺(tái)風(fēng)初始時(shí)刻的關(guān)鍵參數(shù)；基于大西洋海域5°×5°網(wǎng)格內(nèi)歷史臺(tái)風(fēng)數(shù)據(jù)，擬合得到經(jīng)驗(yàn)路徑模型的系數(shù)；結(jié)合初始時(shí)刻臺(tái)風(fēng)關(guān)鍵參數(shù)以及臺(tái)風(fēng)所處相應(yīng)網(wǎng)格的經(jīng)驗(yàn)路徑模型可以預(yù)測(cè)下一時(shí)刻臺(tái)風(fēng)的位置與強(qiáng)度，以此類推可以產(chǎn)生完整的虛擬臺(tái)風(fēng)路徑。Vickery等[14]指出，與Monte-Carlo模擬方法相比，經(jīng)驗(yàn)路徑模擬方法的優(yōu)點(diǎn)在于它不需要依賴子區(qū)域氣候一致性的假設(shè)。因此，即使是在臺(tái)風(fēng)氣候?qū)W變化較大的區(qū)域也可以采用經(jīng)驗(yàn)路徑模型進(jìn)行臺(tái)風(fēng)災(zāi)害分析，這有助于分析大型系統(tǒng)的臺(tái)風(fēng)風(fēng)險(xiǎn)。而Monte-Carlo模擬方法更加適合分析單個(gè)站點(diǎn)的臺(tái)風(fēng)危險(xiǎn)性。美國(guó)建筑結(jié)構(gòu)荷載規(guī)范[21]推薦的設(shè)計(jì)風(fēng)速也是基于Vickery等[14]的經(jīng)驗(yàn)路徑模型。

自Vickery等[14]以來(lái)，臺(tái)風(fēng)危險(xiǎn)性分析方法得到了一些重要改進(jìn)。Li和Hong[22]采用地理加權(quán)回歸方法對(duì)Vickery等[14]的經(jīng)驗(yàn)路徑模型進(jìn)行了簡(jiǎn)化，并對(duì)簡(jiǎn)化后模型進(jìn)行了充分的驗(yàn)證。Vickery和Wadhera[23]及Vickery等[24]利用颶風(fēng)偵察飛機(jī)測(cè)得的壓力數(shù)據(jù)和颶風(fēng)研究部門提供的颶風(fēng)風(fēng)場(chǎng)照片信息，更新了最大風(fēng)速半徑（Rmax）和Holland氣壓剖面參數(shù)（B）的統(tǒng)計(jì)模型。Vickery[25]開(kāi)發(fā)了一種新的颶風(fēng)登陸后的衰減模型。Mudd等[26]采用颶風(fēng)預(yù)測(cè)模型模擬了未來(lái)氣候變化情景下2100年的颶風(fēng)。Lee和Rosowsky[27]分別在2005年以及未來(lái)的氣候情景下為美國(guó)西北部沿岸區(qū)域建立了颶風(fēng)風(fēng)速數(shù)據(jù)庫(kù)。Rosowsky等[28]對(duì)Mudd等[26]的研究進(jìn)行了相應(yīng)的改善。此外，Powell等[29]、Legg 等[30]、Apivatanagul等[31]、Pei等[32]、Li和 Hong[33]都對(duì)經(jīng)驗(yàn)路徑模型進(jìn)行了其他一些衍生研究。

過(guò)去幾十年，我國(guó)在臺(tái)風(fēng)風(fēng)險(xiǎn)分析領(lǐng)域也開(kāi)展了許多研究。歐進(jìn)萍等[34]利用Batts風(fēng)場(chǎng)模型分析了我國(guó)東南沿海重點(diǎn)城市的臺(tái)風(fēng)危險(xiǎn)性。Xiao等[1]結(jié)合Thompson和Cardone[12]臺(tái)風(fēng)模型和Monte-Carlo模擬方法對(duì)我國(guó)東南沿海地區(qū)11個(gè)重點(diǎn)城市進(jìn)行臺(tái)風(fēng)風(fēng)災(zāi)分析。Li和Hong[22]對(duì)經(jīng)驗(yàn)路徑模型進(jìn)行了簡(jiǎn)化并使用該模型分析了我國(guó)東南沿海的臺(tái)風(fēng)災(zāi)害。這些研究為沿海重點(diǎn)城市的抗風(fēng)設(shè)計(jì)提供了有益的參考。

利用Monte-Carlo模擬與經(jīng)驗(yàn)路徑模型分析臺(tái)風(fēng)危險(xiǎn)性的思路是相似的。兩種方法都是首先構(gòu)造大量虛擬臺(tái)風(fēng)，然后采用適當(dāng)?shù)娘L(fēng)場(chǎng)模型計(jì)算這些臺(tái)風(fēng)的風(fēng)速，從中提取最大風(fēng)速構(gòu)成極值風(fēng)速序列，最后采用極值分布進(jìn)行擬合并預(yù)測(cè)不同重現(xiàn)期的極值風(fēng)速。只是在構(gòu)造虛擬臺(tái)風(fēng)的具體過(guò)程上兩種方法有所不同。本文分別采用Monte-Carlo模擬方法以及經(jīng)驗(yàn)路徑模擬方法預(yù)測(cè)了中國(guó)東南沿海區(qū)域臺(tái)風(fēng)的極值風(fēng)速，并對(duì)兩種方法的預(yù)測(cè)結(jié)果以及產(chǎn)生差異的原因進(jìn)行了分析。

2 數(shù)據(jù)與方法

Monte-Carlo模擬與經(jīng)驗(yàn)路徑模型都需要以歷史臺(tái)風(fēng)數(shù)據(jù)為統(tǒng)計(jì)基礎(chǔ)，本文采用的歷史臺(tái)風(fēng)數(shù)據(jù)來(lái)源于中國(guó)氣象局-上海臺(tái)風(fēng)研究所（ China Meteorological Administration and Shanghai Typhoon Institute, CMASTI）西北太平洋熱帶氣旋最佳路徑數(shù)據(jù)集（1949-2017年，www.typhoon.gov.cn）。該數(shù)據(jù)集由中國(guó)氣象局負(fù)責(zé)維護(hù)并進(jìn)行更新。數(shù)據(jù)集包括臺(tái)風(fēng)每6 h間隔的位置與強(qiáng)度信息，具體包括臺(tái)風(fēng)名稱、編號(hào)、臺(tái)風(fēng)中心位置（經(jīng)度與緯度）、中心壓強(qiáng)、近臺(tái)風(fēng)中心的2 min平均最大風(fēng)速。

臺(tái)風(fēng)模型通常用一組參數(shù)進(jìn)行參數(shù)化，稱之為臺(tái)風(fēng)關(guān)鍵參數(shù)，主要包括臺(tái)風(fēng)年發(fā)生率（λ）、中心氣壓差（Δp）、最大風(fēng)速半徑（Rmax）、移動(dòng)速度（VT）、移動(dòng)方向（θ）、臺(tái)風(fēng)移動(dòng)路徑與模擬點(diǎn)間的最小距離（Dmin），這些參數(shù)表征了臺(tái)風(fēng)的結(jié)構(gòu)和強(qiáng)度。

圖 1 虛擬臺(tái)風(fēng)過(guò)程示意圖Fig. 1 A schematic graph for simulation of the virtual typhoon

2.1 Monte-Carlo模擬方法

采用Monte-Carlo模擬方法分析臺(tái)風(fēng)危險(xiǎn)性可分為以下5個(gè)步驟：（1）確定研究地點(diǎn)，采用模擬圓方法從歷史臺(tái)風(fēng)事件集中提取對(duì)該研究點(diǎn)有影響的臺(tái)風(fēng)事件。模擬圓方法是以研究點(diǎn)為圓心，劃定一定的半徑，當(dāng)一個(gè)臺(tái)風(fēng)路徑與模擬圓相交或是位于模擬圓里面的時(shí)候，則該臺(tái)風(fēng)被提取出來(lái)作為對(duì)該研究點(diǎn)有影響的臺(tái)風(fēng)[8,13]。Li和Hong[22]以及Vickery等[15]通過(guò)參數(shù)化實(shí)驗(yàn)，建議模擬圓半徑采用250 km。在本文中我們也設(shè)定模擬圓的半徑為250 km，在該范圍內(nèi)認(rèn)為臺(tái)風(fēng)氣候是一致的。（2）基于提取的臺(tái)風(fēng)事件，擬合臺(tái)風(fēng)各個(gè)關(guān)鍵參數(shù)的概率模型。基于前人研究[11,13,35]，臺(tái)風(fēng)年發(fā)生率可以用泊松分布、二項(xiàng)分布或負(fù)二項(xiàng)分布來(lái)描述；中心氣壓差可以用對(duì)數(shù)正態(tài)分布、伽馬分布或是威布爾分布來(lái)擬合；臺(tái)風(fēng)中心移動(dòng)速度可以用正態(tài)分布、對(duì)數(shù)正態(tài)分布或伽馬分布來(lái)描述；臺(tái)風(fēng)中心移動(dòng)方向可以用正態(tài)分布、雙峰分布或是Von Mises分布來(lái)擬合；最小距離可以用均勻分布或梯形分布來(lái)描述。（3）采用Monte-Carlo隨機(jī)抽樣從關(guān)鍵參數(shù)概率模型中抽取一組臺(tái)風(fēng)關(guān)鍵參數(shù)，組成一個(gè)虛擬臺(tái)風(fēng)事件，以此類推可以產(chǎn)生大量的虛擬臺(tái)風(fēng)事件。（4）采用臺(tái)風(fēng)風(fēng)場(chǎng)模型計(jì)算虛擬臺(tái)風(fēng)的風(fēng)速與風(fēng)向，從中提取極值風(fēng)速序列。（5）利用極值分布模型對(duì)極值風(fēng)速序列進(jìn)行擬合進(jìn)而估計(jì)臺(tái)風(fēng)不同重現(xiàn)期的極值風(fēng)速。

其中通過(guò)Monte-Carlo隨機(jī)抽樣產(chǎn)生虛擬臺(tái)風(fēng)的步驟如下：（1）判斷從概率分布隨機(jī)抽取的關(guān)鍵參數(shù)是否在合理的范圍內(nèi)：Δp為[0, 135] hPa，VT為[2, 65]km/h，θ為 [-180°, 180°]，|Dmin|≤250 km。（2）用Dmin和θ確定臺(tái)風(fēng)的直線路徑，該直線路徑與模擬圓的交點(diǎn)即為虛擬臺(tái)風(fēng)的起始點(diǎn)，如圖1所示。（3）在臺(tái)風(fēng)移動(dòng)過(guò)程中移動(dòng)速度與方向保持不變，以此來(lái)確定臺(tái)風(fēng)其他時(shí)刻的中心位置。（4）在臺(tái)風(fēng)登陸前，中心壓強(qiáng)保持不變，登陸后采用衰減模型來(lái)描述中心壓差Δp的變化。（5）當(dāng)臺(tái)風(fēng)中心離開(kāi)模擬圓范圍則認(rèn)為臺(tái)風(fēng)消失。該過(guò)程示意圖如圖1所示。

表 1 不同地形地貌的粗糙度Table 1 Roughness length of different geomorphology

2.2 衰減模型

臺(tái)風(fēng)登陸后，由于海洋能量來(lái)源的損失和地面摩擦的增加，臺(tái)風(fēng)強(qiáng)度會(huì)減弱。臺(tái)風(fēng)衰減的速率與臺(tái)風(fēng)本身、登陸的位置以及登陸時(shí)的中心壓差有關(guān)。采用合適的衰減模型描述臺(tái)風(fēng)的衰減規(guī)律，對(duì)于預(yù)測(cè)臺(tái)風(fēng)登陸后的風(fēng)速有重要影響。本文采用Vickery和Twisdale[9]提出的衰減模型，即

式中，Δp（t）是臺(tái)風(fēng)登陸后t時(shí)刻的中心壓差（單位：hPa）；Δp0是登陸時(shí)刻的中心壓差（單位：hPa）；a是特定點(diǎn)衰減系數(shù)；ε是正態(tài)分布的誤差項(xiàng)，均值為0；標(biāo)準(zhǔn)偏差為σε。通過(guò)擬合登陸的臺(tái)風(fēng)數(shù)據(jù)，可以得到研究點(diǎn)或區(qū)域的臺(tái)風(fēng)衰減模型的系數(shù)。

2.3 風(fēng)場(chǎng)模型

臺(tái)風(fēng)風(fēng)場(chǎng)模型是進(jìn)行臺(tái)風(fēng)數(shù)值求解的關(guān)鍵。本文采用了Yan Meng (YM)風(fēng)場(chǎng)模型，該模型是Meng等[10]在1995年提出的解析模型，其數(shù)值形式如下：

式中，V代表臺(tái)風(fēng)風(fēng)速；ρ代表空氣密度；k代表垂向單位向量；f代表科氏參數(shù)；F代表邊界層摩擦力。該模型以擾動(dòng)平衡方程和經(jīng)過(guò)邊界層摩擦力修正的壓力梯度平衡方程為基礎(chǔ)，引入了可以考慮地貌等復(fù)雜因素的“等效粗糙長(zhǎng)度”，能夠很好地表示出臺(tái)風(fēng)風(fēng)場(chǎng)中各種力之間的相互關(guān)系，并且能夠比較直觀地求解梯度風(fēng)方程，得到完整的解析解。基于之前的研究[3，36-37]不同地形地貌的粗糙度如表1所示。

YM模型具有足夠的精度并被Matsui等[38]、謝汝強(qiáng)等[39]所采用。關(guān)于風(fēng)場(chǎng)模型的詳細(xì)內(nèi)容請(qǐng)參考文獻(xiàn)[10]。

氣壓是臺(tái)風(fēng)風(fēng)場(chǎng)模型的重要輸入?yún)?shù)，YM風(fēng)場(chǎng)模型采用Holland[40]氣壓模型，其形式為

式中，p(r)是距臺(tái)風(fēng)中心徑向距離為r處的海平面壓強(qiáng)；p0是臺(tái)風(fēng)中心的氣壓值。式中有兩個(gè)未知的參數(shù)，分別是最大風(fēng)速半徑（Rmax）和氣壓剖面參數(shù)（B）。這兩個(gè)參數(shù)是臺(tái)風(fēng)風(fēng)場(chǎng)模型中的重要參數(shù)，它們的計(jì)算公式直接決定風(fēng)場(chǎng)模型計(jì)算風(fēng)速的大小，進(jìn)而影響極值風(fēng)速的預(yù)測(cè)。本文采用的計(jì)算公式取自Vickery和Wadhera[23]，該公式同時(shí)被Li和Hong[22,33,41]以及Hong等[42]所使用，形式如下：

2.4 經(jīng)驗(yàn)路徑模型

Vickery等[14]發(fā)展了經(jīng)驗(yàn)路徑模型，對(duì)臺(tái)風(fēng)的移動(dòng)速度、移動(dòng)方向以及相對(duì)強(qiáng)度進(jìn)行了建模。模型表達(dá)式如下：

式中，c表示臺(tái)風(fēng)移動(dòng)速度；θ表示臺(tái)風(fēng)移動(dòng)方向；I表示臺(tái)風(fēng)相對(duì)強(qiáng)度 (臺(tái)風(fēng)壓強(qiáng))；ai(i=1, 2,···,5)是移動(dòng)速度 系 數(shù) ；bi(i=1,2,···,6)是 移 動(dòng) 方 向 系 數(shù) ；di(i=1,2,···,6)是相對(duì)強(qiáng)度系數(shù)；ψ和λ是臺(tái)風(fēng)的緯度與經(jīng)度；ci是i時(shí)刻的移動(dòng)速度；θi、θi-1是i與i-1 時(shí)刻的移動(dòng)方向；Ii+1、Ii,、Ii-1、Ii-2分別是時(shí)刻i+1、i、i-1 以及i-2的相對(duì)強(qiáng)度；Tsi+1、Tsi是i+1與i時(shí)刻海表面溫度（單位：K）；ε是誤差項(xiàng)。

臺(tái)風(fēng)相對(duì)強(qiáng)度的概念是由Darling[43]建立的，

式中，pda和pdc分別表示環(huán)境和最小可持續(xù)干空氣分子壓強(qiáng)。臺(tái)風(fēng)中心壓強(qiáng)可以用相對(duì)強(qiáng)度I來(lái)表示，反之也是一樣。關(guān)于相對(duì)強(qiáng)度的具體計(jì)算方法可以查看Darling[43]附錄。

Vickery等[14]首先將大西洋海盆劃分成 5°×5°的網(wǎng)格，然后通過(guò)對(duì)各個(gè)網(wǎng)格單元中歷史臺(tái)風(fēng)的移動(dòng)速度、移動(dòng)方向和相對(duì)強(qiáng)度數(shù)據(jù)回歸分析得到東向與西向運(yùn)動(dòng)臺(tái)風(fēng)的ai、bi、di的兩套數(shù)值。然后采用直接抽樣方法從HURDAT[44]臺(tái)風(fēng)數(shù)據(jù)集中抽取臺(tái)風(fēng)起始時(shí)刻的位置、日期、移動(dòng)方向和移動(dòng)速度，以此來(lái)初始化模擬過(guò)程。接著根據(jù)公式（7）至公式（9）計(jì)算得到的臺(tái)風(fēng)移動(dòng)速度與方向，可以計(jì)算下一個(gè)6 h間隔的臺(tái)風(fēng)位置與強(qiáng)度，以此類推可以得到完整的臺(tái)風(fēng)路徑。對(duì)于那些沒(méi)有臺(tái)風(fēng)數(shù)據(jù)或是臺(tái)風(fēng)數(shù)據(jù)比較少的網(wǎng)格點(diǎn)，可以用最鄰近網(wǎng)格中的模型來(lái)代替。當(dāng)臺(tái)風(fēng)登陸后，采用Vickery和Twisdale[9]提出的衰減模型描述中心壓強(qiáng)的衰減。如果臺(tái)風(fēng)再次回到海里，則采用公式（10）繼續(xù)計(jì)算臺(tái)風(fēng)在新的網(wǎng)格位置的中心壓強(qiáng)。

在Vickery等[14]模型中，每個(gè)網(wǎng)格單元需要確定的系數(shù)比較多。Li和Hong[22]在模型回歸分析中消除了一些次要的解釋變量，對(duì)經(jīng)驗(yàn)路徑模型進(jìn)行了簡(jiǎn)化，克服了變量多重相關(guān)性問(wèn)題，形式如下：

Li和Hong[22]對(duì)簡(jiǎn)化的路徑模型進(jìn)行了充分的驗(yàn)證，并將其用于分析中國(guó)東南沿海區(qū)域臺(tái)風(fēng)的危險(xiǎn)性[33]。

在本文中，我們將簡(jiǎn)化的經(jīng)驗(yàn)路徑方法應(yīng)用到整個(gè)西北太平洋區(qū)域，構(gòu)造了該海域1 000年熱帶氣旋事件集。臺(tái)風(fēng)年頻數(shù)服從泊松分布，均值為27.27 個(gè)/a。臺(tái)風(fēng)發(fā)生的起始位置和日期是從CMA數(shù)據(jù)集中直接抽樣而得到。圖2給出了CMA數(shù)據(jù)集中所有臺(tái)風(fēng)的起始點(diǎn)分布以及西北太平洋區(qū)域5°×5°的網(wǎng)格劃分。基于歷史臺(tái)風(fēng)數(shù)據(jù)，臺(tái)風(fēng)的初始移動(dòng)速度符合伽馬分布，形狀參數(shù)為3.74，尺度參數(shù)為5.07；臺(tái)風(fēng)初始移動(dòng)方向符合雙峰分布，均值為-69.28 和-17.71，標(biāo)準(zhǔn)偏差為18.80和62.78，權(quán)重系數(shù)為0.514，擬合的結(jié)果分別如圖3a和圖3b所示。泊松分布、伽馬分布以及雙峰分布都在5%的顯著性水平上通過(guò)了Kolmogorov-Smirnov（KS）擬合優(yōu)度檢驗(yàn)。臺(tái)風(fēng)的初始移動(dòng)速度與移動(dòng)方向分別從上述伽馬分布與雙峰分布中隨機(jī)抽樣而得到，臺(tái)風(fēng)初始相對(duì)強(qiáng)度是從經(jīng)驗(yàn)分布抽取得到。

給定臺(tái)風(fēng)初始移動(dòng)方向、移動(dòng)速度以及強(qiáng)度后，路徑模型可以基于當(dāng)前時(shí)刻的移動(dòng)速度與移動(dòng)方向估計(jì)臺(tái)風(fēng)新的位置。基于CMA數(shù)據(jù)集我們擬合的公式（11）中東向與西向運(yùn)動(dòng)臺(tái)風(fēng)的系數(shù)ai和εc如圖4所示。同樣我們也擬合了公式（12）以及公式（13）中的系數(shù)b和d，限于文章篇幅，不再給出這些系數(shù)的數(shù)值分布。

為了擬合衰減模型，將西北太平洋臺(tái)風(fēng)登陸區(qū)域劃分成5部分：30°N以北溫帶氣旋區(qū)域（區(qū)域1），25°～30°N 臺(tái)灣島以北區(qū)域（區(qū)域 2），20°～25°N 包含臺(tái)灣島區(qū)域（區(qū)域3），菲律賓群島區(qū)域（區(qū)域5），剩下的區(qū)域（區(qū)域4）。對(duì)衰減模型中擬合的系數(shù)如表2所示。

基于臺(tái)風(fēng)的初始化過(guò)程以及擬合的經(jīng)驗(yàn)路徑模型與衰減模型構(gòu)造了西北太平洋1 000年的虛擬臺(tái)風(fēng)數(shù)據(jù)集。為了驗(yàn)證該數(shù)據(jù)集的可靠性，對(duì)比了中國(guó)沿岸46個(gè)站點(diǎn)（圖5）觀測(cè)與虛擬的臺(tái)風(fēng)關(guān)鍵參數(shù)。用于驗(yàn)證的觀測(cè)臺(tái)風(fēng)來(lái)源于CMA數(shù)據(jù)集。首先將每個(gè)站點(diǎn)250 km范圍內(nèi)的觀測(cè)與虛擬臺(tái)風(fēng)數(shù)據(jù)分別提取出來(lái)，然后統(tǒng)計(jì)這些臺(tái)風(fēng)事件在距離各個(gè)站點(diǎn)最近時(shí)（除了中心壓差）的關(guān)鍵參數(shù)，包括臺(tái)風(fēng)年發(fā)生率的平均值、移動(dòng)速度、移動(dòng)方向與最小距離（當(dāng)臺(tái)風(fēng)經(jīng)過(guò)研究點(diǎn)的右側(cè)時(shí)，最小距離取正值，反之取負(fù)值）的平均值與標(biāo)準(zhǔn)差。中心壓差取臺(tái)風(fēng)事件在距各站點(diǎn)250 km范圍內(nèi)的壓差最小值。驗(yàn)證結(jié)果如圖6所示，可以看出模擬的我國(guó)沿海臺(tái)風(fēng)的特征與實(shí)測(cè)數(shù)據(jù)吻合較好，說(shuō)明經(jīng)驗(yàn)路徑模型能夠較好地再現(xiàn)西北太平洋臺(tái)風(fēng)的統(tǒng)計(jì)特征。

圖 2 CMA數(shù)據(jù)中臺(tái)風(fēng)的起始點(diǎn)以及西北太平洋區(qū)域5°×5°的網(wǎng)格Fig. 2 Spatial distribution of genesis based on the CMA best-track dataset and the grid cell of 5°×5° in the Northwest Pacific Ocean

圖 3 觀測(cè)的臺(tái)風(fēng)初始移動(dòng)速度（a）與移動(dòng)方向（b）的概率分布擬合Fig. 3 Modeled and observed statistical distributions of translation speed (a) and storm heading (b)

3 特定點(diǎn)虛擬臺(tái)風(fēng)的模擬與計(jì)算

我們將中國(guó)東南沿海臺(tái)風(fēng)多發(fā)區(qū)域作為本文的研究區(qū)域。采用網(wǎng)格法來(lái)確定研究點(diǎn)，網(wǎng)格的分辨率是 0.25°×0.25°（參考 Hong 等[42]）。一共有 579 個(gè)網(wǎng)格點(diǎn)（研究點(diǎn)）覆蓋整個(gè)研究區(qū)域（圖5）。

圖 4 對(duì)于東向與西向運(yùn)動(dòng)的臺(tái)風(fēng)擬合的公式（11）中的系數(shù)ai和εc的分布Fig. 4 Illustration of regression coefficients ai and εc in Eq.(11) for easterly and westerly headed storms

表 2 公式(2)的衰減常數(shù)aTable 2 Decay constant a in Eq. (2)

首先采用Monte-Carlo方法模擬每個(gè)研究點(diǎn)的虛擬臺(tái)風(fēng)事件，模擬的時(shí)間長(zhǎng)度為1 000年。采用模擬圓方法從CMA歷史臺(tái)風(fēng)數(shù)據(jù)集中提取對(duì)各個(gè)點(diǎn)有影響的臺(tái)風(fēng)事件，通過(guò)統(tǒng)計(jì)，每個(gè)研究站點(diǎn)有25～267個(gè)不等的歷史臺(tái)風(fēng)樣本。提取這些臺(tái)風(fēng)的關(guān)鍵參數(shù)，用相應(yīng)的概率分布進(jìn)行擬合，采用卡方檢驗(yàn)與KS檢驗(yàn)兩個(gè)非參數(shù)檢驗(yàn)方法對(duì)臺(tái)風(fēng)關(guān)鍵參數(shù)的擬合優(yōu)度進(jìn)行檢驗(yàn)，檢驗(yàn)的顯著性水平為5%。如果理論分布沒(méi)有通過(guò)檢驗(yàn)，則采用經(jīng)驗(yàn)分布。由此我們確定了每個(gè)研究點(diǎn)各個(gè)臺(tái)風(fēng)關(guān)鍵參數(shù)的最優(yōu)概率模型。利用第2節(jié)介紹的Monte-Carlo方法以及衰減模型，為每個(gè)研究點(diǎn)構(gòu)造了1 000年間的虛擬臺(tái)風(fēng)事件。再利用YM風(fēng)場(chǎng)模型可以計(jì)算每一個(gè)虛擬臺(tái)風(fēng)的風(fēng)速與風(fēng)向。在利用YM模型計(jì)算風(fēng)場(chǎng)時(shí)，根據(jù)實(shí)際地形標(biāo)記了每個(gè)研究點(diǎn)的粗糙度，如圖7所示。

為了比較Monte-Carlo方法與經(jīng)驗(yàn)路徑方法預(yù)測(cè)的極值風(fēng)速的異同，我們從第2.4節(jié)采用經(jīng)驗(yàn)路徑方法構(gòu)造的西北太平洋1 000年的臺(tái)風(fēng)數(shù)據(jù)集中提取了對(duì)各個(gè)研究點(diǎn)有影響的虛擬臺(tái)風(fēng)事件。同樣采用YM風(fēng)場(chǎng)模型計(jì)算了每一個(gè)虛擬臺(tái)風(fēng)在研究點(diǎn)處產(chǎn)生的風(fēng)速與風(fēng)向。

圖 5 研究站點(diǎn)（黑色星號(hào)）、用于驗(yàn)證經(jīng)驗(yàn)路徑模型的站點(diǎn)（藍(lán)色方框）以及8個(gè)沿岸城市（紅點(diǎn)）的位置分布Fig. 5 Locations of research points (black asterisks), stations for validation of empirical track model (blue squares) and 8 coastal cities (red dots)

4 極值風(fēng)速的預(yù)測(cè)與對(duì)比

基于第3節(jié)計(jì)算的虛擬臺(tái)風(fēng)風(fēng)速，可以從中提取每個(gè)臺(tái)風(fēng)在研究點(diǎn)產(chǎn)生的最大風(fēng)速，從而構(gòu)成極值風(fēng)速序列。極值風(fēng)速序列通常會(huì)用威布爾（Weibull）分布、耿貝爾（Gumbel）分布或經(jīng)驗(yàn)分布來(lái)進(jìn)行描述，采用這些極值分布可以預(yù)測(cè)不同重現(xiàn)期的極值風(fēng)速。

對(duì)于經(jīng)驗(yàn)分布，推導(dǎo)出的重現(xiàn)期為T的風(fēng)速表達(dá)式為

式中，m是極值風(fēng)速序列的樣本容量，大括號(hào)中的公式是下標(biāo)i的計(jì)算公式，當(dāng)風(fēng)速樣本從小到大排列時(shí)，i是重現(xiàn)期為T時(shí)極值風(fēng)速對(duì)應(yīng)的秩序。

對(duì)于Weibull分布，推導(dǎo)出的重現(xiàn)期為T的風(fēng)速表達(dá)式為

式中，β為形狀參數(shù)；η為尺度參數(shù)；γ為位置參數(shù)；λ為臺(tái)風(fēng)發(fā)生率。

對(duì)于Gumbel分布，推導(dǎo)出的重現(xiàn)期為T的風(fēng)速表達(dá)式為

式中，μ為位置參數(shù)；α為尺度參數(shù)；λ為 臺(tái)風(fēng)發(fā)生率。

為了驗(yàn)證本文采用Monte-Carlo方法與經(jīng)驗(yàn)路徑方法預(yù)測(cè)的極值風(fēng)速的準(zhǔn)確性，我們將本文預(yù)測(cè)的8個(gè)沿岸重點(diǎn)城市（上海、寧波、溫州、福州、廈門、廣州、深圳、湛江，如圖5所示）的極值風(fēng)速與Li和Hong[33]的結(jié)果以及中國(guó)《建筑結(jié)構(gòu)荷載規(guī)范》[45]中的數(shù)值進(jìn)行了對(duì)比，結(jié)果如圖8所示。本文的預(yù)測(cè)結(jié)果與Li和Hong[33]以及設(shè)計(jì)規(guī)范的值吻合較好。

本文對(duì)Monte-Carlo方法與經(jīng)驗(yàn)路徑方法預(yù)測(cè)的不同重現(xiàn)期極值風(fēng)速進(jìn)行了對(duì)比。圖9a為采用Monte-Carlo方法以及經(jīng)驗(yàn)分布預(yù)測(cè)的50年重現(xiàn)期極值風(fēng)速等值線圖，圖9b是相同條件下采用經(jīng)驗(yàn)路徑方法預(yù)測(cè)的極值風(fēng)速與Monte-Carlo方法預(yù)測(cè)的風(fēng)速的差值。正值表示經(jīng)驗(yàn)路徑方法預(yù)測(cè)結(jié)果高于Monte-Carlo方法預(yù)測(cè)的結(jié)果。可以看出在研究區(qū)域的內(nèi)陸側(cè)，經(jīng)驗(yàn)路徑方法預(yù)測(cè)的結(jié)果略高于Monte-Carlo方法預(yù)測(cè)的結(jié)果；在余下的大部分區(qū)域尤其是海岸沿線一帶，經(jīng)驗(yàn)路徑方法預(yù)測(cè)的結(jié)果比Monte-Carlo方法預(yù)測(cè)的結(jié)果低0～5 m/s，最大差值為-13.3 m/s，最大相對(duì)差值為31.5%。圖10為采用Monte-Carlo方法以及經(jīng)驗(yàn)分布預(yù)測(cè)的100年重現(xiàn)期極值風(fēng)速等值線圖以及經(jīng)驗(yàn)路徑方法與Monte-Carlo方法預(yù)測(cè)結(jié)果的差值。可以看出100年重現(xiàn)期極值風(fēng)速普遍比50年重現(xiàn)期極值風(fēng)速高0～5 m/s，同樣在內(nèi)陸側(cè)少數(shù)地區(qū)經(jīng)驗(yàn)路徑方法預(yù)測(cè)的結(jié)果略高，在余下的大部分區(qū)域尤其是海岸沿線一帶經(jīng)驗(yàn)路徑方法預(yù)測(cè)的結(jié)果偏低，最大差值為-16.8 m/s，最大相對(duì)差值為33.5%。

圖 6 沿岸46個(gè)站點(diǎn)虛擬與觀測(cè)的臺(tái)風(fēng)關(guān)鍵參數(shù)的對(duì)比Fig. 6 Comparison of key parameters of simulated and observed typhoons at 46 coastal stations

圖11、圖12與圖9類似，只是極值分布分別采用的是威布爾分布和耿貝爾分布。對(duì)比圖9a、圖11a以及圖12a可以看出威布爾分布預(yù)測(cè)的極值風(fēng)速略低于經(jīng)驗(yàn)分布的預(yù)測(cè)結(jié)果，最大差值為-5.1 m/s，最大相對(duì)差值為11.4%；而耿貝爾分布預(yù)測(cè)的極值風(fēng)速略高于經(jīng)驗(yàn)分布的結(jié)果，最大差值為5.3 m/s，最大相對(duì)差值為11.6%。圖11b與圖12b對(duì)比的經(jīng)驗(yàn)路徑方法與Monte-Carlo方法預(yù)測(cè)的極值風(fēng)速的差異與圖9b是類似的，最大差值分別為-11.0 m/s和-10.2 m/s，最大相對(duì)差值分別為27.4%和25.3%。

為了解釋經(jīng)驗(yàn)路徑方法與Monte-Carlo方法預(yù)測(cè)的極值風(fēng)速的差異，計(jì)算了對(duì)各個(gè)站點(diǎn)產(chǎn)生影響的臺(tái)風(fēng)在距離各站點(diǎn)最近時(shí)的中心壓差的平均值。因?yàn)橹行膲翰畋碚髁伺_(tái)風(fēng)的強(qiáng)度，直接影響各個(gè)站點(diǎn)極值風(fēng)速的大小。其中對(duì)各個(gè)站點(diǎn)產(chǎn)生影響的臺(tái)風(fēng)數(shù)據(jù)分別來(lái)自Monte-Carlo方法與經(jīng)驗(yàn)路徑方法構(gòu)造的虛擬臺(tái)風(fēng)事件集。計(jì)算了基于這兩種方法得到的579個(gè)研究站點(diǎn)的中心壓差的平均值的差值，結(jié)果如圖13a所示，其中正值表示經(jīng)驗(yàn)路徑方法得到的中心壓差的平均值大于Monte-Carlo方法的結(jié)果。對(duì)比圖13a與圖9b可以看出，兩者有著類似的正負(fù)分布，圖13b給出了兩者正負(fù)號(hào)一致的站點(diǎn)（紅點(diǎn)）與不一致的站點(diǎn)（藍(lán)點(diǎn)）分布。可以看出對(duì)于絕大多數(shù)研究站點(diǎn)，經(jīng)驗(yàn)路徑方法得到的中心壓差的平均值偏大的地方，其預(yù)測(cè)的極值風(fēng)速也偏大。因此我們認(rèn)為經(jīng)驗(yàn)路徑方法與Monte-Carlo方法預(yù)測(cè)的極值風(fēng)速的差異，主要是由兩者構(gòu)造的虛擬臺(tái)風(fēng)的中心壓強(qiáng)存在差異造成的。這可能是由于Monte-Carlo方法假定臺(tái)風(fēng)登陸前壓強(qiáng)保持不變，因此在臺(tái)風(fēng)登陸時(shí)中心壓強(qiáng)差偏大，導(dǎo)致其預(yù)測(cè)的靠近海邊的站點(diǎn)極值風(fēng)速也偏大。對(duì)于還有少數(shù)站點(diǎn)不能用壓強(qiáng)差解釋的風(fēng)速差異可能是模型本身的不確定性造成的，比如臺(tái)風(fēng)路徑的不確定性、衰減模型的差異等。

圖 7 不同研究點(diǎn)的粗糙度Fig. 7 The roughness length for each site of interest

圖 8 本文預(yù)測(cè)的8個(gè)沿岸城市極值風(fēng)速與Li和Hong[33]以及規(guī)范值[45]的對(duì)比Fig. 8 Comparison of extreme wind speeds of 8 coastal cities between this study, Li and Hong[33] and the design code[45]

圖 9 采用Monte-Carlo方法以及經(jīng)驗(yàn)分布預(yù)測(cè)的50年重現(xiàn)期的極值風(fēng)速（a）以及經(jīng)驗(yàn)路徑方法與Monte-Carlo方法預(yù)測(cè)結(jié)果的差值(b)Fig. 9 Extreme wind speeds for 50-year return period by Monte-Carlo simulation and empirical distribution (a) and the wind speed difference between empirical track model and Monte-Carlo simulation (b)

圖 10 采用Monte-Carlo方法以及經(jīng)驗(yàn)分布預(yù)測(cè)的100年重現(xiàn)期的極值風(fēng)速（a）以及經(jīng)驗(yàn)路徑方法與Monte-Carlo方法預(yù)測(cè)結(jié)果的差值（b）Fig. 10 Extreme wind speeds for 100-year return period by Monte-Carlo simulation and empirical distribution (a) and the wind speed difference between empirical track model and Monte-Carlo simulation (b)

5 結(jié)論

本文分別采用Monte-Carlo模擬方法以及經(jīng)驗(yàn)路徑模擬方法預(yù)測(cè)了中國(guó)東南沿海區(qū)域臺(tái)風(fēng)的極值風(fēng)速。將整個(gè)東南沿海區(qū)域劃分為0.25°×0.25°的網(wǎng)格，將每個(gè)網(wǎng)格點(diǎn)作為研究點(diǎn)。基于歷史臺(tái)風(fēng)數(shù)據(jù)，采用Monte-Carlo方法產(chǎn)生每個(gè)研究點(diǎn)1 000年間的虛擬臺(tái)風(fēng)事件。采用經(jīng)驗(yàn)路徑模型方法構(gòu)建了整個(gè)西北太平洋1 000年的熱帶氣旋事件集，并對(duì)構(gòu)造的虛擬臺(tái)風(fēng)的特征進(jìn)行了驗(yàn)證。從構(gòu)造的熱帶氣旋事件集中可以提取對(duì)各個(gè)研究站點(diǎn)有影響的臺(tái)風(fēng)事件。采用YM風(fēng)場(chǎng)模型計(jì)算虛擬臺(tái)風(fēng)的風(fēng)場(chǎng)，采用極值分布模型預(yù)測(cè)每個(gè)研究點(diǎn)不同重現(xiàn)期的極值風(fēng)速。對(duì)比兩種方法預(yù)測(cè)的極值風(fēng)速，研究發(fā)現(xiàn)在研究區(qū)域內(nèi)陸側(cè)，經(jīng)驗(yàn)路徑方法預(yù)測(cè)的結(jié)果略高于Monte-Carlo方法預(yù)測(cè)的結(jié)果；在大部分海岸沿線一帶，經(jīng)驗(yàn)路徑方法預(yù)測(cè)的結(jié)果比Monte-Carlo方法預(yù)測(cè)的結(jié)果略低，這主要是由兩種方法構(gòu)造的虛擬臺(tái)風(fēng)的中心壓強(qiáng)存在差異以及模型本身的不確定性造成的。本文研究結(jié)果可以為相關(guān)部門防災(zāi)減災(zāi)提供一定的參考。

圖 11 采用Monte-Carlo方法以及威布爾分布預(yù)測(cè)的50年重現(xiàn)期的極值風(fēng)速（a）以及經(jīng)驗(yàn)路徑方法與Monte-Carlo方法預(yù)測(cè)結(jié)果的差值（b）Fig. 11 Extreme wind speeds for 50-year return period by Monte-Carlo simulation and Weibull distribution (a) and the wind speed difference between empirical track model and Monte-Carlo simulation (b)

圖 12 采用Monte-Carlo方法以及耿貝爾分布預(yù)測(cè)的50年重現(xiàn)期的極值風(fēng)速（a）以及經(jīng)驗(yàn)路徑方法與Monte-Carlo方法預(yù)測(cè)結(jié)果的差值（b）Fig. 12 Extreme wind speeds for 50-year return period by Monte-Carlo simulation and Gumbel distribution (a) and the wind speed difference between empirical track model and Monte-Carlo simulation (b)

圖 13 采用Monte-Carlo方法以及經(jīng)驗(yàn)路徑方法得到的579個(gè)研究站點(diǎn)的中心壓差的平均值的差值（a）以及圖13a與圖9b正負(fù)號(hào)一致的站點(diǎn)（紅點(diǎn)）與不一致的站點(diǎn)（藍(lán)點(diǎn)）（b）Fig. 13 The difference of the mean value of the central pressure difference for 579 research stations between Monte-Carlo simulation and empirical track method (a), and stations with the same (red dots) and different (blue dots) signs in Fig.13a and Fig.9b (b)