林 青 潘 寧 劉 銘

?

光流法以及結(jié)合機(jī)器學(xué)習(xí)在臨近預(yù)報(bào)中的應(yīng)用

林 青 潘 寧 劉 銘

福建省氣象臺(tái)

通過(guò)采用Gunnar Farneback的稠密光流跟蹤法對(duì)雷達(dá)反射率進(jìn)行2h外推，并結(jié)合半拉格朗日外推方法，保持天氣系統(tǒng)回波的旋轉(zhuǎn)效應(yīng)，可有效提高回波短臨預(yù)測(cè)的精準(zhǔn)度。結(jié)合該方法外推的回波結(jié)果及歷史降水觀測(cè)，采取機(jī)器學(xué)習(xí)中ConvLSTM(空間深度學(xué)習(xí))的模型結(jié)合Z-R關(guān)系進(jìn)行1h降水預(yù)測(cè)，結(jié)果表明其降水預(yù)測(cè)的空漏報(bào)率低于SWAN業(yè)務(wù)系統(tǒng)QPF預(yù)測(cè)結(jié)果。

光流法 機(jī)器學(xué)習(xí) 臨近預(yù)報(bào)

1 概述

強(qiáng)對(duì)流天氣是一種多發(fā)生在中小尺度范圍的、伴有劇烈氣流上升運(yùn)動(dòng)的災(zāi)害性天氣。當(dāng)空氣發(fā)生強(qiáng)烈的垂直運(yùn)動(dòng)時(shí)，就可能發(fā)生強(qiáng)對(duì)流天氣，并且常伴隨有雷雨、大風(fēng)、冰雹、龍卷風(fēng)等危害性很強(qiáng)的現(xiàn)象。然而，由于強(qiáng)對(duì)流天氣發(fā)生的速度快、范圍小、時(shí)間短，如何準(zhǔn)確預(yù)報(bào)出強(qiáng)對(duì)流天氣發(fā)生的時(shí)間和地點(diǎn)一直是氣象學(xué)研究的重點(diǎn)和難點(diǎn)。

為了對(duì)強(qiáng)對(duì)流天氣進(jìn)行預(yù)測(cè)，在氣象研究和業(yè)務(wù)中，臨近預(yù)報(bào)會(huì)根據(jù)短時(shí)間內(nèi)雷達(dá)回波、云圖演變等，進(jìn)行推演和判斷，預(yù)測(cè)未來(lái)強(qiáng)對(duì)流天氣可能發(fā)生的區(qū)域和時(shí)間。一方面，外推判斷依賴人的經(jīng)驗(yàn)；另一方面，研究如何實(shí)現(xiàn)氣象數(shù)據(jù)外推自動(dòng)化。

早期有拉格朗日持續(xù)性預(yù)報(bào)、最優(yōu)化相關(guān)系數(shù)的交叉相關(guān)法等，這些方法主要是基于臨近短時(shí)間的氣象數(shù)據(jù)做簡(jiǎn)單線性的外推，沒(méi)有考慮到氣象變化的一些非線性關(guān)系。為了改進(jìn)回波外推技術(shù)，引入了在計(jì)算機(jī)視覺(jué)領(lǐng)域應(yīng)用較廣的光流法來(lái)追蹤回波的運(yùn)動(dòng)。該方法在1950 年由Gibson[1]首先提出，當(dāng)被觀測(cè)的目標(biāo)和傳感器之間有相對(duì)運(yùn)動(dòng)時(shí)，所觀察到的亮度模式運(yùn)動(dòng)稱之為光流。對(duì)于變化較快的強(qiáng)對(duì)流降水天氣過(guò)程，光流法具有明顯的預(yù)報(bào)優(yōu)勢(shì)，曹春燕等[2]研究認(rèn)為光流法預(yù)報(bào)效果優(yōu)于SWAN（Severe Weather Automatically Nowcast System）業(yè)務(wù)系統(tǒng)中使用的交叉相關(guān)法。光流法可以彌補(bǔ)傳統(tǒng)交叉相關(guān)法的缺陷，提升對(duì)流臨近預(yù)報(bào)系統(tǒng)的性能。因?yàn)樵跒槠骄彽膶訝钤平邓到y(tǒng)中，交叉相關(guān)法的效果相對(duì)理想。但在強(qiáng)對(duì)流天氣中，跟蹤失敗的情況會(huì)顯著增加，從而影響到最終的預(yù)報(bào)結(jié)果，而光流法采用全局平滑的假設(shè)使得風(fēng)暴單體的追蹤計(jì)算在整個(gè)區(qū)域內(nèi)都得到了準(zhǔn)確的運(yùn)動(dòng)矢量。張蕾等[3]利用金字塔分層思想和RPM-SL方案改進(jìn)HS光流法后發(fā)現(xiàn)，對(duì)強(qiáng)對(duì)流回波預(yù)報(bào)具有優(yōu)勢(shì)。

2 雷達(dá)回波外推方法

2.1 光流法簡(jiǎn)介

光流是指空間運(yùn)動(dòng)物體在觀察成像平面上的像素運(yùn)動(dòng)的瞬時(shí)速度，是利用圖像序列中像素在時(shí)間域上的變化以及相鄰幀之間的相關(guān)性來(lái)找到上一幀跟當(dāng)前幀之間存在的對(duì)應(yīng)關(guān)系，從而計(jì)算出相鄰幀之間物體的運(yùn)動(dòng)信息的一種方法。

光流場(chǎng)就是類剛體物體的速度矢量場(chǎng)，光流場(chǎng)與交叉相關(guān)法所求得的運(yùn)動(dòng)矢量場(chǎng)在物理意義上是等價(jià)的：如果把傳感器從攝像機(jī)換成多普勒雷達(dá)，被探測(cè)的目標(biāo)從一般的運(yùn)動(dòng)物體換成雷達(dá)回波，則可以用光流法替代交叉相關(guān)法得到運(yùn)動(dòng)矢量場(chǎng)。

2. 2 半拉格朗日外推方法

由于臺(tái)風(fēng)、強(qiáng)天氣中的渦度大，單體運(yùn)動(dòng)軌跡都有一定旋轉(zhuǎn)效應(yīng)，為了解決線性外推不考慮旋轉(zhuǎn)的缺點(diǎn)，采用了半拉格朗日平流方案對(duì)雷達(dá)回波進(jìn)行外推。半拉格朗日方法具有良好的穩(wěn)定性和精度，目前已廣泛應(yīng)用到數(shù)值預(yù)報(bào)模式和氣候模式，采用RPM-SL ( remapped particle-mesh semi- Lagrangian) 方案外推回波，該方案可保持回波的旋轉(zhuǎn)性，提高回波預(yù)報(bào)效果。

2.3 應(yīng)用分析

本文選取兩次對(duì)流過(guò)程資料，利用光流法對(duì)雷達(dá)組合反射率因子進(jìn)行0～2h的外推預(yù)報(bào)，并與SWAN產(chǎn)品（即交叉相關(guān)法）對(duì)比分析。

2.3.1 2017年4月11日福建強(qiáng)對(duì)流天氣過(guò)程

2017年4月11日福建強(qiáng)對(duì)流天氣過(guò)程是一般性對(duì)流天氣過(guò)程，11日，湖南南部到福建中部出現(xiàn)帶狀的強(qiáng)對(duì)流天氣分布區(qū)；南面廣東上空出現(xiàn)分散性強(qiáng)對(duì)流天氣，閩粵交界地帶也出現(xiàn)明顯的強(qiáng)對(duì)流天氣，15時(shí)（北京時(shí)間）前后回波從江西進(jìn)入福建省，最強(qiáng)的回波強(qiáng)度在40～45dBZ。30dBZ以上的強(qiáng)回波區(qū)呈現(xiàn)弓形分布（圖略），回波后側(cè)有偏北風(fēng)侵入，南側(cè)為西南風(fēng)，總體為弓形回波的結(jié)構(gòu)。

利用光流法結(jié)合拉格朗日跟蹤外推算法，外推預(yù)報(bào)（圖1c）得到的強(qiáng)對(duì)流的形狀、位置和雷達(dá)觀測(cè)到的實(shí)況（圖1a）比較接近。系統(tǒng)較好地預(yù)報(bào)出了福建中部、南部的回波及福建西部的強(qiáng)回波。當(dāng)然，比較圖1a、1c也可以看出，回波南段的移動(dòng)速度偏慢。對(duì)比預(yù)報(bào)和實(shí)況可知，光流法基本預(yù)報(bào)出了最大反射率因子特征的位置。實(shí)況回波上的強(qiáng)中心主要分布在西段、中段和南段，西段主要為大面積回波，最強(qiáng)回波約40dBZ左右，中段和南段為類弓形回波，最強(qiáng)回波約為50dBZ。光流法預(yù)報(bào)后的強(qiáng)回波中心也主要在西段、中段和南段，南段的最強(qiáng)回波比實(shí)況偏南，中段相當(dāng)，但最大值即50dBZ范圍大于實(shí)況回波。從本次過(guò)程１h預(yù)報(bào)以及１h后的回波實(shí)況對(duì)比可知，預(yù)報(bào)1h以后強(qiáng)回波位置與回波實(shí)況比較一致，形狀也與實(shí)況比較接近，而SWAN產(chǎn)品（圖1b）盡管回波位置相近，中段強(qiáng)度明顯偏弱偏北。因此利用光流法對(duì)大范圍的對(duì)流系統(tǒng)1h外推預(yù)報(bào)也有較大的參考價(jià)值。

圖1 2017年4月11日19時(shí)(北京時(shí)間)福建強(qiáng)對(duì)流天氣過(guò)程

a——實(shí)況回波；b——SWAN 10時(shí)預(yù)報(bào)1h外推回波；c——光流法10時(shí)預(yù)報(bào)1h后的外推回波

2.3.2 2017年4月24日福建對(duì)流天氣過(guò)程

2017年4月24日，福建省上空出現(xiàn)分散性強(qiáng)對(duì)流天氣，強(qiáng)度相對(duì)一般，全省多地有局地的對(duì)流過(guò)程，較11日強(qiáng)度小，風(fēng)暴主體在浙江境內(nèi)，福建以40dBZ以下回波為主，光流法外推預(yù)報(bào)（圖2c）得到的強(qiáng)對(duì)流的形狀、位置與雷達(dá)觀測(cè)到的實(shí)況（圖2a）比較接近。SWAN產(chǎn)品（圖2b）明顯強(qiáng)度較弱，35dbz以上區(qū)域預(yù)報(bào)偏少，由此可見(jiàn)對(duì)于分散性對(duì)流，光流法效果也較SWAN方法理想。

圖2 2017年4月24日22時(shí)(北京時(shí)間)福建對(duì)流天氣過(guò)程

a——實(shí)況回波；b——SWAN21時(shí)預(yù)報(bào)1h外推回波；c——光流法21時(shí)預(yù)報(bào)1h外推回波

3 臨近降水估測(cè)

3.1 算法介紹

近些年來(lái)，隨著數(shù)據(jù)的不斷豐富，深度學(xué)習(xí)應(yīng)運(yùn)而生。深度學(xué)習(xí)可以從大量數(shù)據(jù)中自動(dòng)學(xué)習(xí)數(shù)據(jù)的特征，并映射到輸出空間中，在圖像識(shí)別、文本預(yù)測(cè)等領(lǐng)域已經(jīng)取得了非常好的效果。

對(duì)于雷達(dá)、氣象這樣的時(shí)空數(shù)據(jù)，近年來(lái)也有一些研究工作。例如，Srivastava等人提出了基于整個(gè)圖像的推演模型FC-LSTM[4]，通過(guò)將全圖展開(kāi)為一維向量，并輸入到LSTM模型中，再通過(guò)卷積神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)模型。

雷達(dá)回波圖中，可以得到一定范圍內(nèi)每個(gè)格點(diǎn)的雷達(dá)回波率。根據(jù)Wilson等人[5]的研究，雷達(dá)回波率和降水率的關(guān)系可以表示為：

=2001.6

雷達(dá)回波圖可以很好地反映地區(qū)的降水量，對(duì)預(yù)測(cè)強(qiáng)對(duì)流具有非常重要的作用。但在實(shí)際過(guò)程中，雷達(dá)回波和降水的關(guān)系可能受地理位置、降水類型、發(fā)展過(guò)程的影響，系數(shù)值200、指數(shù)值1.6往往需要依靠人工的經(jīng)驗(yàn)。本方法利用深度學(xué)習(xí)神經(jīng)網(wǎng)絡(luò)，通過(guò)歷史數(shù)據(jù)學(xué)習(xí)來(lái)找到雷達(dá)回波和降水強(qiáng)度之間的關(guān)系，實(shí)驗(yàn)結(jié)果證明，該方法可以取得非常好的效果。

卷積神經(jīng)網(wǎng)絡(luò)（CNN）是一種前饋神經(jīng)網(wǎng)絡(luò)，通過(guò)疊加卷積層、池化層、RELU層等操作來(lái)對(duì)圖像進(jìn)行處理，并生成圖像的特征。

卷積層通過(guò)不同的卷積對(duì)圖像進(jìn)行處理，可以生成圖像的局部特征。一幅圖像可能存在各種各樣不同的局部特征，傳統(tǒng)的圖像處理算法SIFT等通過(guò)手工提取特征來(lái)對(duì)圖像進(jìn)行分類，然而，這一方法不僅需要大量人工的參與，也很難取得很好的效果。近年來(lái)，隨著數(shù)據(jù)的不斷增多和GPU處理能力的增速，通過(guò)誤差反向傳播的BP算法可以自動(dòng)學(xué)習(xí)神經(jīng)網(wǎng)絡(luò)中的卷積核，從而提取更好的圖像特征。圖3展示了一個(gè)3×3卷積的示例。

圖3 卷積過(guò)程示例

3.2 模型建立

收集2015年12月～2017年6月共19個(gè)月的雷達(dá)回波圖，建立數(shù)據(jù)集。雷達(dá)回波圖為從0dBZ到70dBZ以上的不同等級(jí)的基本回波率值?；夭实膁BZ值越高，一般表示該區(qū)域的降雨量越大。在雷達(dá)圖中，用前17個(gè)月的數(shù)據(jù)作為訓(xùn)練集，倒數(shù)第二個(gè)月的數(shù)據(jù)作為驗(yàn)證集，最后一個(gè)月作為測(cè)試集。在對(duì)未來(lái)天氣進(jìn)行預(yù)測(cè)中，我們用過(guò)去5個(gè)小時(shí)的歷史雷達(dá)回波圖作為輸入，并預(yù)測(cè)未來(lái)3個(gè)小時(shí)的雷達(dá)回波圖。

實(shí)際模型計(jì)算時(shí)，設(shè)置了以下一些神經(jīng)網(wǎng)絡(luò)參數(shù)：將卷積神經(jīng)網(wǎng)絡(luò)設(shè)置為Alexnet的結(jié)構(gòu)，并去掉了最后的兩層全連接層；模型中一共有6層LSTM（長(zhǎng)短期記憶網(wǎng)絡(luò)）層。使用前3層LSTM對(duì)輸入進(jìn)行編碼，并用后3層來(lái)對(duì)編碼后的隱藏單元進(jìn)行解碼。對(duì)于編碼網(wǎng)絡(luò)層，每層的隱藏單元數(shù)分別為128、64、64。對(duì)于解碼網(wǎng)絡(luò)層，每層的隱藏單元數(shù)也為128、64、64。輸入數(shù)據(jù)為訓(xùn)練過(guò)去三個(gè)小時(shí)每6min（12min）雷達(dá)回波，過(guò)去三個(gè)小時(shí)自動(dòng)站站點(diǎn)累積降水量。通過(guò)模型計(jì)算未來(lái)三個(gè)小時(shí)每6min雷達(dá)回波矩陣，進(jìn)而換算未來(lái)三小時(shí)各個(gè)站點(diǎn)累積降水量預(yù)測(cè)。具體步驟為：在該模型中，首先對(duì)輸入回波進(jìn)行兩次卷積，這里假設(shè)在每個(gè)卷積核大小的范圍內(nèi)，降水和回波具有相近的關(guān)系。然后對(duì)卷積后的結(jié)果進(jìn)行Z-R積分，積分的參數(shù)是動(dòng)態(tài)可以優(yōu)化的。最后在用兩層反卷積得到降水輸出結(jié)果。模型的結(jié)果如圖4所示。

圖4 降水量卷積模型

3.3 應(yīng)用分析

同樣以4月11日強(qiáng)對(duì)流天氣過(guò)程為例。

圖5為4月11日過(guò)程小時(shí)降水預(yù)測(cè)，可見(jiàn)模型預(yù)測(cè)與實(shí)況相對(duì)接近，隨著起報(bào)時(shí)間臨近，強(qiáng)降水范圍與實(shí)況更接近，集中在福建南側(cè)，圖5b為預(yù)報(bào)2h后結(jié)果，強(qiáng)降水范圍接近，但強(qiáng)度偏大，但隨著起報(bào)時(shí)間接近，明顯強(qiáng)降水范圍減小到與實(shí)況一致。

選取2017年4月～6月資料，計(jì)算了逐月的TS和漏報(bào)率評(píng)分（表1），與SWAN產(chǎn)品對(duì)比，模型預(yù)測(cè)TS接近，但漏報(bào)率遠(yuǎn)小于SWAN產(chǎn)品，因此對(duì)模型的業(yè)務(wù)化應(yīng)用有一定價(jià)值。

圖5 2017年4月11日20時(shí)1小時(shí)降水

a——實(shí)況19～20時(shí)1h累計(jì)降水；b——18時(shí)起報(bào)預(yù)報(bào)19～20時(shí)1h累計(jì)降水；c——19時(shí)起報(bào)預(yù)報(bào)19～20時(shí)1h累計(jì)降水

表1 2017年4～6月逐月的TS和漏報(bào)率評(píng)分

4 結(jié)論

本文提出了聯(lián)用光流法和拉格朗日外推模型對(duì)雷達(dá)回波圖進(jìn)行推演預(yù)測(cè)。通過(guò)兩個(gè)對(duì)比，發(fā)現(xiàn)方案優(yōu)于SWAN產(chǎn)品，在強(qiáng)度預(yù)測(cè)上有明顯優(yōu)勢(shì)；同時(shí)，利用CNN+LSTM的深度學(xué)習(xí)模型對(duì)降水進(jìn)行預(yù)測(cè)，結(jié)果也較有參考價(jià)值，在后續(xù)的研究中需要更多個(gè)例驗(yàn)證。同時(shí)光流法外推和交叉相關(guān)均存在局限，即短時(shí)間局地生消預(yù)測(cè)還有一定難度，需要引入模式動(dòng)力機(jī)制影響改進(jìn)回波生消的預(yù)測(cè)方案。

[1] Gibson JJ. The ecological approach to visual perception[M]. Boston: Houghton Mifflin, 1950: 332.

[2] 曹春燕，陳元昭，劉東華，等．光流法及其在臨近預(yù)報(bào)中的應(yīng)用[J]．氣象學(xué)報(bào)，2015, 73(3)：471-480

[3] 張蕾, 魏鳴, 李南, 等. 改進(jìn)的光流法在回波外推預(yù)報(bào)中的應(yīng)用[J]. 科學(xué)技術(shù)與工程, 2014,14(32): 133-137.

[4] Srivastava N, Mansimov E, Salakhutdinov R. Unsupervised learning of video representations using LSTMs[C]//International Conference on Machine Learning, 2015: 843-852.

[5] Wilson J W, Brandes E A. Radar measurement of rainfall—A summary[J]. Bulletin of the American Meteorological Society, 1979, 60:1048-1058.

[6] Krizhevsky A, Sutskever I, Hinton G E. ImageNet classification with deep convolutional neural networks[J]. Advances in Neural Information Processing Systems, 2012(25): 1106-1114.