基于時空深度學習模型的數值降水預報后處理

鄭超昊，尹志偉，曾鋼鋒，許月萍，周鵬，劉莉

(1.浙江大學 建筑工程學院，浙江 杭州 310058；2.臺州市水利局 浙江 臺州 318000；3.臺州市水利水電勘測設計院有限公司 浙江 臺州 318000)

降水是水循環的重要組成部分，也是地表水文過程的基本驅動元素[1].降水預報是進行流域水文預報、水資源規劃以及防洪防汛等研究的基礎資料.降水時空差異性大，因此如何獲取高精度、高分辨率的降水預報是水文領域一大難題[2].確定性降水預報無法完全模擬大氣過程和初始場資料觀測誤差的影響，降雨預報的模式逐漸從確定性預測轉向多成員集合預報，以集合的形式來描述未來降雨發生的概率情況[3].由于初值擾動、模式設計、集合成員數等方面的不完善因素，導致集合預報產品存在預報誤差、集合成員離散度低的問題[4]而無法直接使用，必須運用后處理方法提升集合預報的可靠性和準確性，才能保證高質量的降水預報信息供給[5].

在對數值降水預報的統計后處理方法研究上，代刊等[6]系統介紹了定量降水集合預報的統計學后處理技術；Zhao等[7]評估了分位數映射方法對季節性降雨預測的校正；Andrew等[8-9]發現，參數化后處理方法貝葉斯聯合概率建模的校正效果較全面.深度學習方法的學習和泛化能力出色[10]，具有非線性關系數據處理能力和數據稀疏地區較好表現的獨特優勢，在偏差糾正后處理研究中被廣泛應用[11]：Qiu等[12]以多任務卷積神經網絡考慮時間序列和站點之間的相關關系來進行降水預報，Krasnopolsky等[13]使用人工神經網絡(artificial neural network, ANN)對降水集合預報進行后處理，皇甫江等[14]利用深度學習開展偏振雷達定量降水估測研究.在降水預報后處理的應用研究中，CNN和LSTM模型鮮有提及.

本研究提出結合CNN和LSTM的時空深度后處理模型(CNN-LSTM).該模型具有對時間序列數據進行分析的能力[15]，能夠有效提取數據的空間信息特征特點，還能夠同時考慮預報降水數據空間特征和時間相關性.本研究以浙江省椒江流域為對象，將不同預報時效的歐洲中期天氣預報中心（European centre for medium-range weather forecasts,ECMWF）數值預報網格數據作為輸入數據，以CMA-CMORPH實測雨量數據作為參照，對后處理前后降水資料進行全面評估.

1 研究區域與數據

1.1 研究區域概況

椒江流域地處浙江中東部，東臨東海，介于東經120°21′～121°47′，北緯28°30′～29°31′之間，主要位于臺州市境內，屬于亞熱帶季風氣候，多年平均氣溫為17.1 ℃，多年平均降雨量為1 634.9 mm.該區域的降水具有明顯的季節性，全年降水主要集中在4—6月梅雨季和7—10月臺風期；徑流量年內分配主要集中在4—10月.永安溪和始豐溪是位于椒江上游的2條最主要河流[16].永安溪為椒江的源頭，長144 km，流域面積為2 704 km2.始豐溪是椒江流域二級支流，為椒江水系最大支流，主流長132.7 km，流域面積為1 610 km2.如圖1為椒江流域所處浙江省位置分布圖，以及永安溪、始豐溪流域位置以及水系分布圖.

圖1 椒江流域的流域分布、水系、降水格點數據示意圖Fig.1 Schematic diagram of basin distribution, water system and precipitation grid data in Jiaojiang river basin

1.2 研究數據

1.2.1 ECMWF降水預報數據 使用的預報數據為ECMWF數值降水預報產品，從ECMWF的官網（www.ecmwf.int）下載得到.該數據由51個預報成員組成，空間分辨率為0.5°×0.5°，選取協調世界時零點（北京時間08:00）起預報的未來0～72 h逐6 h的降水預報進行后處理.考慮到降水預報可能存在的落區誤差，選取的經緯度范圍為東經120°～122°，北緯27°30′～29°30′，該范圍完全覆蓋椒江流域及其周邊區域.以2008—2018年汛期4—10月逐6 h降水資料作為輸入數據，椒江流域在2012、2015和2016年遭受典型臺風事件，并且這些年份的降水數據完整性高，無需插補，故選定2012、2015和2016年作為測試期，其余時間的數據作為率定期.

1.2.2 CMA-CMORPH融合數據 CMA-CMOPRH融合降水數據作為驗證數據，從國家氣象信息中心下載得到，為中國30 000個自動站與CMORPH降水產品融合的逐時降水量網格數據集，擁有2008—2018年分辨率為0.1°的逐小時降水數據.融合數據能夠有效減少系統偏差和隨機誤差，CMA-CMOPRH融合降水數據總體誤差水平不超過10%，優于同類型產品在研究區域的精度[17].椒江流域47個實測雨量基本站點逐小時降水數據不連續且不充分，不足以支持深度學習模型的建模，而CMA-CMOPRH融合降水數據精度較高，數據較長，因此利用該降水產品替代實測數據對ECMWF數據進行后處理.

2 研究方法

2.1 數據處理方法

在CMA-CMORPH融合數據中，選取東經120°30′～122°，北緯28°～29°30′，2008—2018年汛期4—10月逐小時降水格點數據.將數據進行6 h累計計算（北京時間8:00點作為起始），篩出數據中的異常值和缺失值，利用min-max標準化方法進行數據歸一化，將6 h累計數據按照日期排序.將ECMWF數值預報降水數據轉化成0～72 h預報時效的6 h累計降水預報.篩出各預報時效內數據的異常值和缺失值，利用min-max標準化對51個成員的各數據進行歸一化處理，再將處理后的數據按照日期排序，以保證數據的一致性.單個成員無法很好地代表ECMWF數據的整體性，而集合平均（等權重算術平均）數據可以涵蓋所有成員的預報能力，因此采用51個成員集合平均數據作為指標評判.將集合平均數據進行降尺度，使其分辨率為0.1°.利用實測降水資料的空間分布進行降尺度可以更好地捕捉降水的空間特征，但實際預報應用中不存在實測降水資料，因此本研究選用反距離權重插值法（IDW）進行數據降尺度，得到與CMA-CMORPH數據具有相同空間尺度的ECMWF數據用于后續評估.

2.2 結合CNN和LSTM的降水時空后處理模型

將CMA-CMORPH柵格數據和ECMWF降水預報數據作為輸入數據，CNN提取分析空間特性，LSTM獲取時間相關性[18]，CNN-LSTM進行預報降水數據的時空誤差后處理.

2.2.1 模型基本原理 CNN基礎構成包含3種類型的層級：卷積層、池化層以及全連接層.卷積層通過卷積運算來提取輸入數據中的特征[18].在池化層中，濾波器將卷積層提取的特征信息進行壓縮，并且保留提取到的最顯著特征.CNN可以通過多層網絡學習不同感知域對目標網格和周圍區域的空間依賴性[19].LSTM是特殊的循環神經網絡（RNN），它時序前饋連接的特性，使過去時刻的信息可以對當前時刻的輸出產生影響.相較于傳統ANN，LSTM具有4種結構：更新門、輸入門、遺忘門和輸出門.更新門用于儲存過去信息的累計；遺忘門可以選擇將過去信息中影響較小的因素“遺忘”，使模型累計的系統誤差減少，得到比傳統RNN更優的結果[19].

如圖2所示為CNN-LSTM的模型結構，該模型是時間序列數據集從輸入到結果輸出的降水后處理過程.

圖2 CNN-LSTM的模型結構Fig.2 Model structure of CNN-LSTM

2.2.2 模型結構 CNN-LSTM包括CNN空間特征提取、LSTM時間相關性分析2個模塊.每個空間位置的降水不孤立存在，而是與周圍空間相關聯，利用CNN模塊提取目標周圍區域的降水信息，為當前目標點的修正提供參考.過去時間的降水信息會影響當前時刻降水，為此，利用LSTM模塊時序前饋的特點提取降水前后方向的時間相關性.由此可知，CNN-LSTM能夠兼顧時間、空間，可以運用于降水后處理中.

從CNN模塊來看，由于預報數據的不確定性會導致單獨成員的預報結果有較大的偏差，將51個預報成員的數據融合，生成數據集{51, 5, 5}，并將率定期所有數據組連接成連續時間序列數據集作為耦合模型的輸入數據，即該數據集包括不同成員的預報信息要素，還能包含空間和時間信息.輸入數據先通過二維卷積層Conv2D，再通過最大值池化層MaxPool提取最大值特征，利用Dropout方法防止模型過擬合.先將數據根據設定參數進行卷積、最大值池化，再向量平鋪成一維向量.將CNN模塊輸出得到的向量與神經節點為256個的全連接層相連，再將全連接層轉化成（16×16）的中間結果輸入LSTM模塊.從LSTM模塊來看，16×16的二維中間向量，分別經過神經元結點為64、32、16的雙向LSTM模型，每一步都通過Dropout方法，再由激活函數Relu輸出16×16的二維向量與實測數據進行訓練.

本研究模型的優化器采用Adam，訓練迭代次數為1 000，每組訓練的數據量為64，耦合模型的學習率最初為0.001，每迭代400次其學習率變為原來的0.1，模型損失函數為均方誤差MSE、兩幅圖像相似度的指標SSIM兩者加權的結果，權重設置分別為0.7、0.3.

式中：x、y均為對比圖像，μx、μy分別為x、y像素的平均值，σx、σy分別為x、y像素方差，σxy為x、y的像素協方差，c1、c2均為穩定常數.

2.3 模型評估方法

以下4個性能指標用于定量評估CNN-LSTM的后處理能力：均方根誤差RMSE、相對偏差RB、相關系數CC、平均絕對誤差MAE.指標的計算式[20]分別為

式中：n為各個網格點觀測值的數量，為網格預測降水，為網格實測觀察降水，為網格預測降水的平均值，為網格實測觀察降水的平均值.

為了進行分級檢驗評估，在ECMWF降水預報資料中提取24 h累計降水，計算永安溪、始豐溪、平原區域面雨量.根據《中短期天氣預報質量檢驗辦法》降水等級劃分如表1所示，將降水預報和降水實況按照24 h累計降水量劃分成對應等級.對各降水分級進行檢驗，各等級TS（threat score）評分為

表1 降水等級劃分表Tab.1 Classification of precipitation levelsmm

式中：k為1～4級，分別代表小雨、中雨、大雨、暴雨及以上降水預報；NAk為k級降水預報正確次數；NBk為k級降水空報次數；NCk為k級降水漏報次數.ECMWF預報數據采用集合平均的數據進行TS評分計算.

3 結 果

評估驗證期內不同預報時效的預報后處理數據，對比分析不同預報時效內永安溪、始豐溪、平原區域面雨量；選取驗證期內的多場典型暴雨進行分析，著重分析永安溪、始豐溪、平原區域面雨量.

3.1 不同預報時效后處理效果評估

計算12個預報時效內每個網格點的數據，對比后處理前、后預報數據的指標來評估CNNLSTM.流域全部格點不同預報時效的各指標計算結果如圖3(小提琴圖)所示，其中虛線為25%、50%、75%分位線，t為預報時效.原始預報RMSE普遍超過3 mm，最大值超過6 mm；經CNNLSTM后處理，各個預報時效的RMSE均有所降低，基本不超過2 mm，最大值為3 mm.經CNNLSTM后處理，RB絕大部分落在-0.25附近區域，而原始預報的RB跨度大.可以看出，原始預報并不是預報時效越短預報技能越高，且預報技能也不隨著預報時效增長呈基本的單調增減，這是氣象預報領域較為常見的“預報不一致性”現象，即在連續多次預報中，前后兩次預報的預報結果差異較大.該現象被認為主要與初始誤差和模式誤差有關，相關分析可參考文獻[21].

圖3 模型評估指標在所有格點不同預報時效后處理前、后的對比圖Fig.3 Comparison of model assessment indexes before and after post-processing for all grids with different lead times

如圖4所示為不同預報時效后處理前、后各指標的空間分布圖.可以看出，不同預報時效后處理后的RMSE、RB、MAE、CC在全流域均有很大的提升.其中RMSE從2.6～6.0 mm下降為1.0～3.0 mm，后處理前、后RMSE下降的格點比例為100%；RB從-0.6～1.0穩定至-0.35～-0.15；MAE從0.7～1.8 mm下降為0.3～0.7 mm，CC從0.2～0.5上升至0.70～0.90；后處理前、后MAE、CC指標改善的格點比例為100%.綜合所述，CNN-LSTM對于數值預報降水產品不同預報時效的預報精度均有顯著提升，并且對于整個椒江流域空間精度提高具有普遍性.

圖4 不同預報時效后處理前、后模型評估指標的空間分布圖Fig.4 Spatial distribution of model assessment indexes before and after post-processing for different lead times

如圖5所示，永安溪、始豐溪流域以及平原區域不同預見期面雨量結果相似.圖中，下標數字0代表后處理前，下標數字1代表處理后.各個預報時效的RMSE均落在基準線下側，說明后處理后數據的RMSE均小于后處理前數據的RMSE，且后處理前RMSE為2.2～4.0 mm，處理后RMSE為1.0～2.0 mm.各個預報時效的CC均落在基準線上方，后處理前的CC為0.0～0.4，處理后的CC＞0.8.各預報時效指標的MAE均落在基準線下方，后處理前的MAE為0.75～1.75 mm，處理后的MAE為0.25～0.75 mm.各個預報時效的RB分布在基準線兩側，其中后處理后的RB穩定在-0.25附近，后處理前的RB為-0.75～0.25，可見處理后的RB穩定在更合理的范圍.綜上所述，在永安溪、始豐溪流域，4個模型評估指標經過CNN-LSTM模型后處理后的預報面雨量數據較于后處理之前在數據精度方面有很大程度的提升.

圖5 模型評估指標在永安溪、始豐溪流域面雨量不同預報時效后處理前、后的對比圖Fig.5 Comparison of model assessment indexes of areal rainfall before and after post-processing over Yonganxi and Shifengxi watershed for different lead times

3.2 分等級降雨評估

如表2所示為計算所得永安溪、始豐溪、平原區域的各量級降水出現的累計頻次N.如表3所示為ECMWF降水預報后處理前、后永安溪、始豐溪、平原區域的TS評分.由表可知，經過后處理3個區域的TS分別為0.900、0.906、0.894，均高于后處理之前的0.754、0.755、0.752，分別提升了19.4%、20%和18.9%.結果證實經CNN-LSTM后處理的雨季預報正確率有顯著提升.

表2 永安溪、始豐溪、平原區域降水分級累計出現頻次Tab.2 Cumulative occurrence frequency of precipitation classification in Yongan, Shifeng and plain regions

表3 ECMWF預報產品后處理前、后雨季3個區域的TS評分項Tab.3 TS scores items of ECMWF forecast products before and after post-processing at three regions during rainy season

如圖6所示為各區域不同降水級別的TS評分圖.后處理后的各等級累計降水TS評分較之后處理前的均有很大程度提升，尤其是小雨級別的TS，數值由0.77提升至0.91，中雨級別的TS從0.48提升至0.6，大雨級別的TS從0.3提升至0.4，暴雨級別的TS評分也有所提升.在后處理前、后的ECMWF預報產品在降水等級評估中，小雨的預報正確率最高，中雨次之，暴雨的預報正確率較低，這可能是暴雨事件較少導致的.

圖6 ECMWF預報產品后處理前、后的等級累計降水TS評分對比圖Fig.6 Comparison of TS scores of ECMWF forecast products before and after post-processing

3.3 典型暴雨后處理效果分析

選取2012年8月“?？迸_風、2015年8月“蘇迪羅”臺風作為典型暴雨，分析永安溪流域、始豐溪流域、平原區域未來3 d累計6 h面雨量.選取2012年8月6日8:00起未來72 h數據作為“?？迸_風的預報數據，并且利用分位線代表集合預報中95%置信區間.由圖7可知(圖中，p為降水量)，本場臺風實際暴雨中心位于流域以北區域，流域內最大實測降水量為362.8 mm，原始預報降水中心偏離，最大雨量僅為226.6 mm，小于實測降水量；后處理后預報暴雨中心與實測降水一致，流域內累計降水量為332.6 mm.選取2015年8月7日8:00起未來72 h數據作為“蘇迪羅”臺風的預報數據.本場臺風的暴雨中心落在椒江流域，實測累計降水的降水值為338.1 mm，原始預報降水最大雨量為308.6 mm，后處理后的預報降水最大值約為310.4 mm.由對比可知，在2場臺風期間，后處理后的預報數據不管是降雨累計，還是暴雨中心都比處理前的更為精確.后處理后的預報數據對比原數據有相當程度的提升.

圖7 2場臺風72 h累計降水分布圖Fig.7 Distribution of 72h cumulative precipitation for two typhoons

如圖8所示為臺風期間永安溪、始豐溪流域和平原區域逐6 h面雨量情況，為了進一步評估后處理前、后預報的差異，圖中增加原始集合預報的95%的置信區間結果.在“海葵”臺風期間，3個區域最大累計實測降水量分別為35、40、35 mm，預報后處理后累計降水分別為32、36、33 mm；三者預報后處理后的累計降水誤差均不超過10%.在“蘇迪羅”臺風期間，3個區域最大累計實測降水量分別為50、60、70 mm，后處理后累計降水分別為45、54、64 mm.三者預報后處理后的累計降水誤差均不超過10%，后處理前的累計降水誤差較大.

圖8 2場臺風在3個區域的逐6 h面雨量對比圖Fig.8 Comparison of 6 h areal rainfall of two typhoons at three regions

對2場臺風預報進行評估，結果如表4所示.在“海葵”臺風期間，3個區域后處理后的RMSE下降均超過90%，MAE在后處理后也有大幅度下降， CC在后處理后穩定在0.99，RB在后處理后提升效果均超過89%.在“蘇迪羅”臺風期間，3個區塊后處理后的RMSE分別下降了91.3%、90.7%、91.6%；MAE分別下降了91.3%、90.2%、91.7%；后處理前的CC不足0.55，后處理后的均穩定在0.99；RB分別提升了91.3%、90.1%、91.5%.結果表明CNN-LSTM后處理的數據模擬典型臺風的效果良好，面雨量模擬精度大幅度超過未經后處理的數據.結果也表明原始數值預報降水產品對臺風降水預報能力較弱，通過耦合模型后處理后的降水產品對臺風降水預報能力較強.

4 結 論

（1）相較于原始數值天氣預報數據，后處理后的RMSE的最大值從6 mm降至3 mm；CC的最大值從0.6升至0.9，RMSE、CC的改善格點比例為100%.說明CNN-LSTM能夠有效提高降水預報精度.

（2）后處理后永安溪、始豐溪流域和平原區域不同預報時效面雨量的指標有顯著提升：RMSE全面提升，最大值從4 mm降至2 mm；后處理后的CC均超過0.8；MAE從0.75～1.75 mm降至0.25～0.75 mm；RB后處理后的值優于后處理前的，并穩定在-0.25.

（3）后處理后永安溪、始豐溪以及平原區域的雨季降水預報的TS評分均大于0.89，高于后處理前的TS評分；后處理后的各降水等級TS評分均好于后處理前的，其中小雨TS評分從不足0.8提升至0.91，中雨的TS評分從不足0.5提升至0.6，大雨TS評分從不足0.4提升至0.4.

（4）2場典型臺風降水后處理后的累計降水量和暴雨中心評估準確性高.對于永安溪、始豐溪流域、平原區域逐6 h面雨量降水峰值的預報數據結果誤差均不超過10%，各類指標評提升均超過85%，說明經過CNN-LSTM后處理的ECMWF降水產品可以用于臺風降水預報.

（5）CNN-LSTM適用于數值降水預報后處理，能夠提高定量降水預報的準確性，對防汛減災都具有十分重要的意義.本研究僅針對南方濕潤地區得出結論，后續將針對北方地區開展相關研究.