基于ECMWF集合預(yù)報的海河流域降水概率預(yù)報應(yīng)用和檢驗

徐姝，熊明明，陳法敬

(1.天津市氣象臺，天津 300074；2.天津市氣候中心，天津 300074；3.國家氣象中心，北京 100081)

引言

大氣是一個混沌系統(tǒng)，天氣預(yù)報的不確定性是其根本屬性之一(葉篤正等，2006)，預(yù)報不確定性通常用概率方式來表達(王東海等，2011)。概率天氣預(yù)報的任務(wù)，是要為預(yù)報產(chǎn)品使用者提供預(yù)報制作者對天氣狀況未來發(fā)展趨勢認(rèn)知狀態(tài)的較完整表達；其中，很重要的一點就是對預(yù)報不確定性的定量化表達。同確定性預(yù)報相比，概率預(yù)報可提供更加豐富的預(yù)報信息(俞小鼎等，2004；李俊，2015)。集合預(yù)報是由單一確定論向多值隨機論轉(zhuǎn)換的關(guān)鍵技術(shù)(王晨稀和姚建群，2008；杜鈞和鄧國，2010；杜鈞和陳靜，2010；杜鈞和李俊，2014)。集合預(yù)報方法產(chǎn)生之前，概率預(yù)報通常是通過統(tǒng)計方法得到的(Wilks，2005)，但統(tǒng)計概率預(yù)報只反映過去出現(xiàn)過的天氣狀態(tài)。常規(guī)統(tǒng)計方法可減小系統(tǒng)性誤差，但對隨機性誤差無可奈何，集合預(yù)報雖能消除隨天氣系統(tǒng)變化帶來的隨機性誤差，其各成員分別表征一種可能的天氣狀況，預(yù)報的不確定性可采用各成員的概率密度來體現(xiàn)，實現(xiàn)更具有動力學(xué)意義的概率預(yù)報(Stephenson et al.，2005)，但歷史預(yù)報資料中的集合預(yù)報性能信息沒有得到有效利用。貝葉斯概率預(yù)報方法可融合歷史資料和集合預(yù)報系統(tǒng)生成的數(shù)據(jù)，給出連續(xù)的概率預(yù)報曲線，受到集合預(yù)報應(yīng)用方面的重點關(guān)注。在氣象領(lǐng)域，馬培迎(1999)較早地應(yīng)用貝葉斯原理，將天氣分為有、無降水，以其期望概率作為先驗概率，對降水概率預(yù)報進行修正，提高了降水預(yù)報精度。陳朝平等(2010)利用貝葉斯公式求得暴雨的后驗分布，并建立了四川境內(nèi)147 個站點的暴雨預(yù)警模型。陳法敬等(2011)將地面氣溫作為預(yù)報量，利用BPO(Bayesian Processor of Output)方法，將單一集合成員預(yù)報轉(zhuǎn)化為概率預(yù)報，得到預(yù)報效果優(yōu)于單一成員的集成貝葉斯概率預(yù)報。劉建國等(2013)利用TIGGE 構(gòu)成的多中心超級集合系統(tǒng)，結(jié)合淮河流域地面觀測建立貝葉斯模型平均(Bayesian Model Average，BMA)概率預(yù)報模型，獲得的地面日均氣溫概率預(yù)報比原始集合預(yù)報效果好。韓焱紅等(2013)通過分別對廣州、武漢、南京、成都和北京集合預(yù)報成員建立貝葉斯產(chǎn)品處理技術(shù)降水概率預(yù)報模型，發(fā)現(xiàn)其可靠性高于基于集合預(yù)報生成的直接概率預(yù)報結(jié)果。張宇彤等(2013)基于貝葉斯原理和貝葉斯降水概率預(yù)報模型，使用中國T213全球集合預(yù)報歷史資料作為先驗信息，對中國不同氣候區(qū)代表站建立貝葉斯降水概率預(yù)報模型，對預(yù)報結(jié)果有一定改進。石嵐等(2017)利用ECMWF將BMA概率預(yù)報模型與確定性數(shù)值模式預(yù)報進行融合，對風(fēng)能預(yù)測產(chǎn)品具有明顯改進。

上述研究表明，應(yīng)用貝葉斯概率預(yù)報方法提取模式產(chǎn)品有效信息對預(yù)報量的氣候(先驗)概率進行修正，得到模式預(yù)報與氣候概率最佳融合的預(yù)報信息可提高預(yù)報準(zhǔn)確率，但先前貝葉斯產(chǎn)品處理技術(shù)的研究主要在做方法的探索，試驗大都只是針對單站，并且試驗時段較短,難以客觀全面地反映貝葉斯預(yù)報的效果，特別是缺乏空間分布上面的表現(xiàn)。隨著智能網(wǎng)格預(yù)報業(yè)務(wù)在氣象領(lǐng)域的深入推廣，急需開展基于網(wǎng)格降水概率預(yù)報的應(yīng)用研究，相比于站點試驗僅能獲得固定站點位置的信息，格點化研究，可以為覆蓋的范圍提供所有位置上概率預(yù)報信息，豐富擴展了該方法可以提供改進信息的空間范圍。本文利用海河流域范圍內(nèi)降水實況和ECMWF 集合預(yù)報數(shù)據(jù)，采用BPO的方法(Krzysztofowicz，1999；Krzysztofowicz et al.，2001)對流域內(nèi)289 個格點建立降水概率預(yù)報模型，根據(jù)預(yù)報能力對模式51 個集合成員進行貝葉斯降水概率預(yù)報信息融合，從而生成集成的海河流域貝葉斯降水概率預(yù)報產(chǎn)品，采用RPS 和BS 評分方法檢驗其在2018年汛期的預(yù)報效果，為概率預(yù)報在業(yè)務(wù)中的應(yīng)用提供有益的指導(dǎo)。

1 方法與資料

1.1 方法

貝葉斯方法的核心是貝葉斯公式，即

上式表達了根據(jù)已有的經(jīng)驗推斷一個先驗概率P(Y)，然后在新信息或者說似然模型P(X|Y)的更新下不斷地重新調(diào)整先驗概率，進而得到事件發(fā)生的后驗概率P(Y|X)的過程。

BPO 方法(韓焱紅等，2013)是通過建立數(shù)值模式預(yù)報的概率化模型，把確定性預(yù)報量轉(zhuǎn)變?yōu)楦怕暑A(yù)報，通過概率來定量化地表達模式預(yù)報值的不確定性。對降水而言，總體可分為有無兩種情況，而在降水發(fā)生后，降水量則變成了一個連續(xù)型的氣象要素。BPO方法的降水概率預(yù)報模型如式(2)或(3)所示，公式表明概率預(yù)報由兩部分構(gòu)成：無降水的概率和條件降水量的概率分布，其中x 為降水預(yù)報因子，y為預(yù)報量，π 表示有降水的后驗概率，Φ()y和φ(y)分別表示了條件降水量的后驗概率分布和后驗密度函數(shù)，δ為狄拉克函數(shù)。其概率分布形式為

由式(4)得到有降水的后驗概率π，其中g(shù)為有降水的先驗概率(降水閾值為0.1mm)，可由歷史降水觀測數(shù)據(jù)統(tǒng)計得到。f0、f1分別為無和有降水時的條件密度函數(shù)，代表預(yù)報因子和預(yù)報量之間的隨機相關(guān)性。

后驗概率分布Φ(y)和后驗概率密度φ(y)由式(5)、(6)計算得到。其概率分布函數(shù)形式為

其中，Q為標(biāo)準(zhǔn)正態(tài)分布函數(shù)；G(y)、K(x)是降水歷史觀測值和預(yù)報值估計得到的條件降水量y的先驗概率分布、預(yù)報因子x的邊緣分布，本文采用威布爾函數(shù)估計G(y)和K(x)，即

由于降水這一要素不符合高斯分布隨機依賴結(jié)構(gòu)，但符合亞高斯結(jié)構(gòu)，因此這里采用Kelly 等(1995)和Krzysztofowicz(1997；2008)提出的亞高斯似然模型，即在正態(tài)-線性似然模型的基礎(chǔ)上，引入正態(tài)分位數(shù)轉(zhuǎn)換(Normal Quantile Transform，NQT)，得到亞高斯似然模型，正態(tài)分位數(shù)轉(zhuǎn)換表示為

Q-1代表標(biāo)準(zhǔn)正態(tài)分布函數(shù)的逆，轉(zhuǎn)換后的變量符合標(biāo)準(zhǔn)正態(tài)分布。c1，c0，T 表示BPO 預(yù)報模型中的后驗參數(shù)，；其中a、b、σ是似然參數(shù)，由正態(tài)分布變量u、z 的均值μ0、μ1，方差以及協(xié)方差σ10通過公式得到。

1.2 資料

本文選取海河流域范圍內(nèi)(111°—120°E，34.5°—43°N)289個格點進行試驗(圖1)。將集合預(yù)報中單一集合成員24 h 降水預(yù)報值作為預(yù)報因子x，格點24 h降水量作為預(yù)報量y。

圖1 海河流域國家氣象站、ECMWF集合預(yù)報格點、代表格點1和2分布圖Fig.1 Distribution of precipitation observation stations,ECMWF ensemble forecast grid points and the representative grid 1 and 2 in Haihe River Basin.

由海河流域1961—2018 年6—8 月國家氣象觀測站逐日20—20 時24 h 降水觀測資料，采用最鄰近法，插值到鄰近格點上，并計算得到各格點的降水先驗概率。為研究BPO方法在集合預(yù)報中的應(yīng)用，本文預(yù)報因子x 選取ECMWF 集合預(yù)報逐日預(yù)報的海河流域格點預(yù)報資料，試驗時段為2012—2018年6—8月期間，其中2012—2017 年為模型建立樣本，采用2018 年樣本進行驗證。

2 建立基于單一預(yù)報成員的貝葉斯概率預(yù)報模型

2.1 正態(tài)分位數(shù)轉(zhuǎn)換

以格點1、格點2為例，圖2a、c表示降水量的控制預(yù)報與觀測值的隨機依賴結(jié)構(gòu)，散點分布圖為聯(lián)合樣本x與y的分布形式，可見降水分布特征不服從標(biāo)準(zhǔn)正態(tài)分布。圖2b、d中散點是預(yù)報量y和預(yù)報因子x經(jīng)過NQT轉(zhuǎn)換后的u和v的分布形式，紅線為y經(jīng)NQT轉(zhuǎn)換后的變量u的概率密度曲線，藍線為x經(jīng)NQT轉(zhuǎn)換后變量z的概率密度曲線，可見，轉(zhuǎn)換后的變量均符合正態(tài)分布。

圖2 降水觀測值y與預(yù)報值x(a.格點1；c.格點2)散點分布圖和NQT轉(zhuǎn)換量u、z(b.格點1；d.格點2)的標(biāo)準(zhǔn)正態(tài)分布Q-1(G(y))及Q-1(K(x))Fig.2 The scatter plots of predict and y and predictor x of(a)grid 1 and(c)grid 2,and the standard normal distribution of NQT transformed u and z of(b)grid 1 and(d)grid 2.

2.2 單一集合成員降水概率預(yù)報結(jié)果分析

基于ECMWF 集合預(yù)報控制預(yù)報，將轉(zhuǎn)換后的先驗分布和預(yù)報值的邊緣分布函數(shù)帶入公式(2)、(3)，建立流域內(nèi)格點2018年6—8月逐24 h降水BPO概率預(yù)報模型，進而生成流域內(nèi)各格點逐24 h的降水概率分布(或概率密度)預(yù)報。

由BPO 模型計算獲得289 個格點的6—8 月逐日有降水的后驗概率πk。這里選取格點1和格點2作為代表，將通過BPO建模修訂得到的后驗概率πk與直接由歷史觀測資料統(tǒng)計得到的先驗概率gk進行比較(圖3)，從中看出，后驗概率πk較先驗概率gk逐日變化起伏更加明顯，說明將模式預(yù)報作為預(yù)報因子提供的有效信息發(fā)揮了一定修訂作用。

圖3 格點1(a)、格點2(b)的6—8月逐日有降水先驗概率(gk)與后驗概率(πk)Fig.3 The prior probability of precipitation occurrence gk and posterior probability πk of each day from June to August at(a)grid 1 and(b)grid 2.

以ECMWF 對格點1 和格點2 的2018 年7 月24 日降水預(yù)報為例，得到基于ECMWF 單一集合成員的貝葉斯概率化降水概率預(yù)報顯示，格點1 的觀測值為117 mm，ECMWF 控制預(yù)報為49.6 mm，修訂后概率密度最大值為60 mm，與確定性預(yù)報相比，預(yù)報準(zhǔn)確性有所提高(圖4b)；格點2 的觀測值為0.4 mm，ECMWF 控制預(yù)報為7.3 mm，修訂后概率密度最大值為1 mm，與觀測值更接近(圖4d)。可見，相比確定性降水預(yù)報，貝葉斯概率預(yù)報給出更多預(yù)報信息，且以概率的形式連續(xù)、定量化給出降水預(yù)報不確定性(韓焱紅等，2013)。

圖4 貝葉斯與氣候降水概率預(yù)報結(jié)果對比:(a)基于ECMWF控制預(yù)報對格點1的2018年7月24日降水預(yù)報得到的概率分布預(yù)報；(b)同(a)，但為概率密度預(yù)報；(c)基于ECMWF控制預(yù)報對格點2的2018年7月24日降水預(yù)報得到的概率分布預(yù)報；(d)同(c)，但為概率密度預(yù)報Fig.4 Comparisons of Bayesian forecast and climate forecast based on deterministic prediction.(a)Cumulative probability prediction taking ECMWF precipitation forecast on 24 July 2018 at grid 1 for example,and(b)the same as(a)but for probability density.(c)Cumulative probability prediction taking ECMWF precipitation forecast on 24 July 2018 at grid 2 for example,and(d)the same as(c)but for probability density.

3 集合預(yù)報各成員的集成降水概率預(yù)報模型

由上文分析可知，通過貝葉斯方法可以將EC?MWF集合預(yù)報的一個成員轉(zhuǎn)化為概率預(yù)報形式，EC?MWF 集合預(yù)報有51 個成員，對應(yīng)有51 個貝葉斯概率預(yù)報，各成員的概率分布各不相同，預(yù)報效果存在差異，需要有一種合理的方法將所有成員的預(yù)報能力差異性或者說不確定信息定量地表現(xiàn)出來。

Krzysztofowicz(1992)提出了表征預(yù)報因子有效信息的評分方法(Informativeness Score,IS)，即

IS值越大，表示預(yù)報因子中有效信息越多，該成員的預(yù)報技巧越高，其值在0～1 之間。2018 年7 月24 日這次過程為例，基于ECMWF 集合預(yù)報成員對格點1的24,48,…,216,240 h 的降水預(yù)報值建立亞高斯似然模型，得到各集合預(yù)報成員的IS 值，結(jié)果如圖5a 所示。不同時效內(nèi)，24 h降水預(yù)報得到的IS 值明顯高于240 h；對于24 h 預(yù)報，各成員中集合成員20 的預(yù)報能力最優(yōu)，集合成員5 最弱，其對應(yīng)的IS 值為0.71 和0.56?？梢姡琁S 值能較好地表征不同集合成員所含有效預(yù)報信息的差異。進而，采用正比于IS3的值作為權(quán)重系數(shù)來體現(xiàn)各集合成員在所有成員中預(yù)報價值的差異。

其中，min函數(shù)用于計算各集合成員IS3的最小值，n為所有集合成員的數(shù)量。通過式(10)融合ECMWF的51個成員的IS值，得到各集合成員的24 h權(quán)重系數(shù)(圖5b)。

權(quán)重系數(shù)同IS值的變化一致，集合成員20的24 h預(yù)報效果最好，故其權(quán)重系數(shù)也最大(0.039)；同理，集合成員5 的24 h 預(yù)報效果最差，其權(quán)重系數(shù)為0 (圖5b)，權(quán)重系數(shù)可較好地表征各集合成員在融合預(yù)報結(jié)果中的地位。因此，可通過權(quán)重系數(shù)ri(ISi)參考式(11)的計算方法對集合成員的概率預(yù)報信息進行歸一化加權(quán)平均，從而得到融合后概率密度預(yù)報結(jié)果。概率分布預(yù)報結(jié)果，進一步由式(12)對概率密度函數(shù)積分得到。

圖5 基于ECMWF集合預(yù)報對格點1降水預(yù)報得到的IS值(a)和權(quán)重系數(shù)r(ISi)(b)Fig.5 The(a)IS and(b) r(ISi) obtained based on ECMWF ensemble forecast for grid 1.

4 集成貝葉斯降水概率預(yù)報結(jié)果分析和檢驗

4.1 應(yīng)用檢驗

分別采用RPS和BS評分，以模式區(qū)域內(nèi)289個格點的24 h降水概率預(yù)報結(jié)果作為統(tǒng)計對象，對試驗結(jié)果進行評價。

RPS評分是一個對所有降水量級的概率預(yù)報結(jié)果綜合評價與衡量的量度，RPS 值在0～1 范圍內(nèi)，且其越小越好。從圖6中看到，對于2018年6—8月每天的降水情況，格點的集成貝葉斯降水概率預(yù)報綜合效果比集合直接概率預(yù)報要好，尤其是“七下八上”華北主汛期；另外，對于48 h和72 h預(yù)報時效，集成貝葉斯降水概率預(yù)報的可靠性較集合直接概率預(yù)報效果更加明顯，隨著預(yù)報時效增加，各個集合成員之間的離散度也隨之增加，導(dǎo)致集合概率預(yù)報的不確定加大，而貝葉斯概率預(yù)報方法對歷史數(shù)值預(yù)報信息的有效融合，所提供的預(yù)報因子有效信息對最終預(yù)報結(jié)果具有改進效果。

BS評分是對降水的某一個等級的概率預(yù)報結(jié)果評價的一個指標(biāo)，在0～1的范圍內(nèi)，其值越小越好。圖7給出了2018年6—8月各降水量級24 h、48 h和72 h降水概率預(yù)報的平均BS 評分情況，3個預(yù)報時效、5 個閾值中集成貝葉斯降水概率預(yù)報的BS 評分均低于集合直接概率預(yù)報，尤其是0.1 mm 和10 mm 兩個閾值，分別減小了19%和31%，且隨著時效的增加，集成貝葉斯降水概率預(yù)報較集合直接概率預(yù)報的優(yōu)勢更加明顯，48 h的0.1 mm和10 mm的BS評分分別減小了35%和40%，72 h分別減小了35%和35%，這與RPS評分效果一致。從圖7 結(jié)果還可看到，對于25 mm 以上降水量級的BS評分明顯降低，這與統(tǒng)計方法有關(guān)，并非降水量級越大概率預(yù)報效果越好。由于參評格點數(shù)不變，而降水范圍是變化的，導(dǎo)致25 mm 以下降水范圍總比25 mm以上降水范圍大，而降水區(qū)域大的要比區(qū)域小的BS 評分大，不同量級降水的BS 評分不具可比性(王晨稀，2005)。

圖7 2018年6—8月海河流域所有格點在不同降水閾值下24 h(a)、48 h(b)、72 h(c)集成貝葉斯降水概率預(yù)報與集合直接概率預(yù)報的BS評分Fig.7 The 24h(a),48h(b)and 72h(c)BS value of integration Bayesian precipitation probability forecast(Integrated Forecast)and ensemble direct probability forecast(Ensemble Forecast)on the grids based on ECMWF ensemble forecast under different precipitation thresholds from July to August in 2018.

4.2 強降水過程檢驗

為檢驗BPO 降水概率預(yù)報模型對智能網(wǎng)格預(yù)報業(yè)務(wù)支撐能力，評估該模型格點預(yù)報場的預(yù)報效果，選取2018 年7 月23 日一次強降水過程做試報檢驗。本次過程受臺風(fēng)“安比”影響，流域自南向北先后出現(xiàn)強降水，流域東部出現(xiàn)大到暴雨，灤河下游、北三河下游、大清河下游以及海河干流出現(xiàn)大暴雨。以2018年7月23日未來24h降水預(yù)報為例，首先對ECMWF各集合降水預(yù)報成員建立BPO概率模型，獲得一組集合成員貝葉斯降水概率預(yù)報，然后根據(jù)式(6)和(7)對51個成員的概率預(yù)報進行融合，最終得到集成的貝葉斯降水概率預(yù)報(圖8)，不同降水量級的高概率區(qū)(圖8a—e)與實況(圖8f)的分布基本一致。從表1也可看出，此次過程24 h的RPS和BS評分，集成貝葉斯降水概率預(yù)報均好于集合直接概率預(yù)報。

表1 24 h降水概率預(yù)報評分Table 1 Value of precipitation probability forecast for 24 h.

圖8 海河流域2018年7月23日20時—24日20時小雨(a)、中雨(b)、大雨(c)、暴雨(d)、大暴雨(e)概率(單位:%)分布及實況降水(單位:mm)分布(f)Fig.8 Probability(unit:%)distributions of(a)light rain,(b)moderate rain,(c)heavy rain,(d)torrential rain,(e)downpour and(f)precipitation(unit:mm)distribution in Haihe River Basin from 20:00 BT 23 July to 20:00 BT 24 July 2018.

綜上所述，基于ECMWF 集合預(yù)報對流域內(nèi)格點2018 年逐日24 h 降水建立降水概率預(yù)報模型并采用RPS和BS評分檢驗。所得結(jié)果表明，格點的集成貝葉斯降水概率預(yù)報的預(yù)報效果優(yōu)于集合直接概率預(yù)報，這種基于貝葉斯理論的集合降水概率預(yù)報方法適用于海河流域。

5 結(jié)論

將BPO 方法應(yīng)用于海河流域進行降水量概率預(yù)報試驗，動態(tài)建立了流域內(nèi)基于ECMWF集合預(yù)報289個格點的BPO概率預(yù)報模型，并對其在研究區(qū)域的預(yù)報技巧進行評估。ECMWF集合預(yù)報每個集合成員性能存在差異，利用正比于IS3的值作為權(quán)重系數(shù)來體現(xiàn)不同集合成員的預(yù)報能力，進行各成員貝葉斯概率預(yù)報信息的融合，得到連續(xù)概率密度函數(shù)形式的集成貝葉斯概率預(yù)報。2018 年6—8 月及“安比”個例的RPS和BS評分檢驗結(jié)果，均表明集成貝葉斯降水概率預(yù)報效果優(yōu)于集合直接概率預(yù)報，因此，貝葉斯降水概率預(yù)報方法可在海河流域網(wǎng)格概率降水預(yù)報業(yè)務(wù)中進行應(yīng)用。

BPO概率預(yù)報給出的概率密度函數(shù)，將預(yù)報的不確定性定量化表達，通過概率密度函數(shù)的解釋應(yīng)用，可在一定程度上提高預(yù)報精度，使得對極端事件具有一定的預(yù)報能力，但仍有可能存在部分極端降水的漏報和空報。本文只選取數(shù)值預(yù)報模式的降水量預(yù)報作為預(yù)報因子，下一步工作還可以嘗試采用多個預(yù)報因子，如：加入比濕、相對渦度等作為預(yù)報因子，利用它們構(gòu)成的一個預(yù)報因子向量進行定量降水概率預(yù)報，從而得到更好的預(yù)報效果。