高云，郭艷萍，周建慧，張葉娥，楊澤民

山西大同大學(xué)計(jì)算機(jī)與網(wǎng)絡(luò)工程學(xué)院，山西大同，037009

0 引言

隨著全球氣候變暖，工業(yè)污染加重，全球生態(tài)遭到嚴(yán)重破壞，引發(fā)的氣候異常問(wèn)題也越來(lái)越嚴(yán)重，中國(guó)的自然災(zāi)害是世界上最嚴(yán)重的國(guó)家之一。在各種重大自然災(zāi)害中，冰雹、暴雨、雷暴等突發(fā)性強(qiáng)的氣象災(zāi)害，發(fā)生頻率最高；又因其程度強(qiáng)而持續(xù)時(shí)間短而難以提前預(yù)知和防范，導(dǎo)致危害程度較大。如2005年07月，印度西部連日特大暴雨引發(fā)多起洪水和泥石流等災(zāi)害，造成至少900人死亡。2005年06月，黑龍江省寧安市洪災(zāi)死亡人數(shù)達(dá)106人，給人民的生命財(cái)產(chǎn)帶來(lái)很大危害。

強(qiáng)對(duì)流天氣是具有重大破壞性的災(zāi)害性天氣之一。目前，通常使用雷達(dá)回波外推的方法對(duì)強(qiáng)對(duì)流天氣進(jìn)行監(jiān)測(cè)和預(yù)警。雷達(dá)回波外推是指對(duì)回波運(yùn)行的速度、方向、強(qiáng)度和形態(tài)變化等特征進(jìn)行跟蹤和預(yù)測(cè)[1-2]。傳統(tǒng)雷達(dá)回波外推方法存在預(yù)測(cè)回波變化較快的天氣過(guò)程時(shí)準(zhǔn)確率大幅下降，預(yù)測(cè)時(shí)效較短，隨著預(yù)測(cè)時(shí)長(zhǎng)的增加準(zhǔn)確性快速下降等不足。

早期的外推方法最常使用TREC技術(shù)，即交叉相關(guān)法。Hilst等先后在文獻(xiàn)[3-6]中應(yīng)用TREC技術(shù)對(duì)整體移動(dòng)的風(fēng)暴簇進(jìn)行了跟蹤，結(jié)果表明該技術(shù)能夠?qū)︼L(fēng)暴簇的內(nèi)部運(yùn)動(dòng)進(jìn)行較好的反演。后經(jīng)過(guò)多位學(xué)者的改進(jìn)，徐亞欽等在文獻(xiàn)[7]中提出了多重動(dòng)態(tài)區(qū)域TREC算法，對(duì)TREC方法進(jìn)行改進(jìn)。單體質(zhì)心法是較常用的一種回波外推法。喬春貴等在文獻(xiàn)[8]中對(duì)整幅雷達(dá)回波圖通過(guò)單體質(zhì)心法進(jìn)行線(xiàn)性外推，得出對(duì)于穩(wěn)定的層狀云該方法雖然有較好的外推能力，但對(duì)其他類(lèi)型的云則外推結(jié)果較差。1950年Gibson在文獻(xiàn)[9]中首先提出了光流法，張蕾等在文獻(xiàn)[10]中對(duì)傳統(tǒng)光流法進(jìn)行了改進(jìn)。陳家慧等在文獻(xiàn)[11]中使用BP神經(jīng)網(wǎng)絡(luò)進(jìn)行混合型降水回波的外推，表現(xiàn)出了較好的外推能力。

雷達(dá)回波外推方法近年來(lái)雖然取得了較大的進(jìn)步，外推預(yù)報(bào)的準(zhǔn)確性逐漸提高，但目前的雷達(dá)回波外推算法，對(duì)于1小時(shí)以上的外推預(yù)報(bào)，準(zhǔn)確性下降很快。雷達(dá)回波外推方法還需要進(jìn)一步地探索和改進(jìn)。

1 ConvLSTM模型

ConvLSTM是LSTM結(jié)構(gòu)的變體，通過(guò)將W的權(quán)值計(jì)算變成了卷積運(yùn)算來(lái)提取出圖像的特征。圖1所示即為其單元結(jié)構(gòu)圖。

ConvLSTM的思路就是使用卷積來(lái)代替矩陣乘法。它在輸入到狀態(tài)以及狀態(tài)到狀態(tài)轉(zhuǎn)換中都采用了卷積結(jié)構(gòu)，通過(guò)堆疊多個(gè)ConvLSTM層并形成編碼預(yù)測(cè)結(jié)構(gòu)來(lái)建立更一般的時(shí)空序列預(yù)測(cè)模型。

該設(shè)計(jì)的一個(gè)顯著特點(diǎn)是所有輸入X1,...,Xt，細(xì)胞輸出C1,...,Ct，隱藏狀態(tài)H1,...,Ht，和ConvLSTM的3個(gè)門(mén)it，ft，ot都是3維張量，2維網(wǎng)格向3維張量的轉(zhuǎn)換如圖2所示。

ConvLSTM內(nèi)部結(jié)構(gòu)如圖3所示。

其對(duì)應(yīng)公式如下：

其中星號(hào)*表示卷積，小圓圈○表示哈達(dá)瑪乘積?？梢钥闯?，這里在輸入門(mén)、遺忘門(mén)和輸出門(mén)這3個(gè)門(mén)的輸入上，都考慮了細(xì)胞狀態(tài)Ct-1。

為了確保細(xì)胞狀態(tài)Ct-1具有與輸入相同的行數(shù)和列數(shù)，在應(yīng)用卷積運(yùn)算之前需要先進(jìn)行padding。將邊界點(diǎn)上隱藏狀態(tài)的填充視為使用外部世界的狀態(tài)進(jìn)行計(jì)算。通常，在第一個(gè)輸入到來(lái)之前，將LSTM的所有狀態(tài)初始化為零，這對(duì)應(yīng)于對(duì)于未來(lái)的“完全無(wú)知”。

ConvLSTM模型用于預(yù)測(cè)網(wǎng)絡(luò)時(shí)其編解碼結(jié)構(gòu)如圖4所示。

ConvLSTM也可以作為更復(fù)雜結(jié)構(gòu)的構(gòu)建塊。對(duì)于時(shí)空序列預(yù)測(cè)問(wèn)題，使用圖 4所示的結(jié)構(gòu)。它包括編碼網(wǎng)絡(luò)和預(yù)測(cè)網(wǎng)絡(luò)2個(gè)網(wǎng)絡(luò)，將編碼網(wǎng)絡(luò)的最后狀態(tài)復(fù)制到預(yù)測(cè)網(wǎng)絡(luò)的初始狀態(tài)和單元輸出，2個(gè)網(wǎng)絡(luò)都是通過(guò)堆疊多個(gè)ConvLSTM層形成的。預(yù)測(cè)目標(biāo)與輸入具有相同的維度，將預(yù)測(cè)網(wǎng)絡(luò)中的所有狀態(tài)連接起來(lái)并發(fā)送到1×1卷積層以生成最終預(yù)測(cè)。

2 實(shí)驗(yàn)過(guò)程

2.1 數(shù)據(jù)準(zhǔn)備

本算法進(jìn)行分析的原始數(shù)據(jù)序列為山西省大同市2019-6-1～2019-6-6共1191幅雷達(dá)站基本反射率圖。雷達(dá)站名：大同；雷達(dá)掃描間隔時(shí)間為5′47"；數(shù)據(jù)范圍200km；顯示仰角：0.5/1.5/2.4。以原始基數(shù)據(jù)的命名方式進(jìn)行命名，方便后期進(jìn)行按順序索引。其中2019-6-5 16:14:55BJT基本反射率圖如圖5所示。

（1）切割

從圖 5可以看出，雷達(dá)站大同掃描范圍為200km，除覆蓋大同地區(qū)以外，還覆蓋了山西省其他地區(qū)以及內(nèi)蒙古自治區(qū)和河北省等部分地區(qū)。本實(shí)驗(yàn)只對(duì)大同地區(qū)的強(qiáng)對(duì)流氣候進(jìn)行預(yù)測(cè)，排除其他區(qū)域的干擾，先設(shè)定圖片輸出分辨率和初始化最終數(shù)組，處理的第一步是對(duì)原始雷達(dá)回波圖像進(jìn)行切割，只切割原始圖像（200,150, 310, 350）的區(qū)域作為本次實(shí)驗(yàn)的數(shù)據(jù)。切割后的圖像如圖 6所示，只保留了大同市周邊地區(qū)。同時(shí)減少比色卡、坐標(biāo)、圖片標(biāo)題等干擾項(xiàng)，增加模型訓(xùn)練的精度。

（2）灰度轉(zhuǎn)換

第二步是將三通道的彩色圖片轉(zhuǎn)換為單通道灰度圖片。由于計(jì)算機(jī)硬件限制，原始繪制的300dpi高分辨率圖像不適合輸入到神經(jīng)網(wǎng)絡(luò)中直接進(jìn)行訓(xùn)練，首先要將裁切后的圖像降低分辨率（壓縮）至100*100。原始雷達(dá)回波圖像為RGB三通道圖像，為減少訓(xùn)練時(shí)間可將原始圖片轉(zhuǎn)換為灰度圖像，每個(gè)像素用8個(gè)bit表示，RGB轉(zhuǎn)換為灰度圖像的公式為：

將圖像抗鋸齒，取消字體平滑，抗混疊。每個(gè)雷達(dá)站每個(gè)時(shí)刻的體掃數(shù)據(jù)均可得到一個(gè)100*100分辨率大小的灰度圖像，轉(zhuǎn)換為灰度的圖像如圖7所示。

（3）得到訓(xùn)練集

增加一個(gè)維度備用，每幅圖像載入內(nèi)存后即為一個(gè)100*100*1大小的三維數(shù)組。處理完單個(gè)體掃數(shù)據(jù)之后再依次遍歷所有圖片，將每個(gè)圖片得到的數(shù)組連接成一個(gè)四維的數(shù)組序列。為了排除雜點(diǎn)的干擾，再次將數(shù)組轉(zhuǎn)化為2維（圖像數(shù)量，100*100），依次遍歷數(shù)組中每個(gè)點(diǎn)的值，如果該點(diǎn)的值＜50，則置為0。得到所需的圖像數(shù)組序列。最后將所得數(shù)組保存為.npz格式的numpy數(shù)組方便后續(xù)取用，這就是最終可以輸入到網(wǎng)絡(luò)中進(jìn)行訓(xùn)練的樣本集。

2.2 ConvLSTM模型建立與分析

使用ConvLSTM模型進(jìn)行雷達(dá)回波圖像外推步驟如下：

（1）訓(xùn)練集和測(cè)試集劃分

ConvLSTM模型樣本數(shù)據(jù)的輸入尺寸為如圖8所示的 5D 張量（samples，time, rows，cols，channels），要提前將訓(xùn)練集和測(cè)試集reshape成如上形式的tensor張量，即樣本總數(shù)、各樣本幀數(shù)、圖像寬度、圖像高度和顏色通道數(shù)。輸出尺寸為，如果返回全部序列，則返回5D張量，即（samples, timesteps, output_row, output_col, filters）；否則，返回4D張量，即（samples,output_row, output_col, filters）。o_row和o_col取決于filter和padding的尺寸。

假如上一層是ConvLSTM2D layer，那么其輸出為以上形式的4D張量或5D張量，當(dāng)后面再接另外一個(gè)layer時(shí)，就要考慮該layer是否能接受4D張量或5D張量，即要考慮ConvLSTM2D的輸出能否作為該layer的輸入。

數(shù)據(jù)準(zhǔn)備生成的.npz格式的numpy數(shù)組是2維的，首先將其讀入并reshape為4維數(shù)據(jù)，即（NUMBER, WIDTH, HEIGHT, 1），樣本集總數(shù)為NUMBER，數(shù)據(jù)集為灰度圖像，將顏色通道數(shù)設(shè)為1。將數(shù)據(jù)集設(shè)置為2個(gè)5D張量的數(shù)組，BASIC_SEQUENCE和NEXT_SEQUENCE，分別存放當(dāng)前時(shí)刻圖像和下一時(shí)刻圖像作為訓(xùn)練集和結(jié)果集，數(shù)據(jù)為5維(NUMBER - FRAMES,FRAMES, WIDTH, HEIGHT, 1)，樣本總量=原樣本總數(shù)-各樣本幀數(shù)。

訓(xùn)練時(shí)將整個(gè)樣本集切割成長(zhǎng)度統(tǒng)一的樣本，16（FRAMES）幀圖像為一個(gè)單元。共分為（NUMBER –FRAMES）個(gè)單元，即可直接訓(xùn)練。

（2）ConvLSTM模型建立

卷積長(zhǎng)短時(shí)記憶神經(jīng)網(wǎng)絡(luò)自提出以來(lái)逐步完善已經(jīng)成為了較為成熟的圖像序列預(yù)測(cè)模型，目前Python中也有不少對(duì)其模型的封裝，方便研究人員直接構(gòu)建網(wǎng)絡(luò)對(duì)自己的數(shù)據(jù)進(jìn)行訓(xùn)練和預(yù)測(cè)。

本文所建模型為如圖 9所示的Keras中的序貫?zāi)Ｐ停粗饘忧短滓来芜B接，數(shù)據(jù)流由輸入端輸入，逐層流動(dòng)，反向傳播時(shí)更新參數(shù)，逐步降低損失函數(shù)。該網(wǎng)絡(luò)由四層ConvLSTM2D網(wǎng)絡(luò)層堆疊而成，最后為一個(gè)三維卷積層用以格式化輸出數(shù)組以便求取損失函數(shù)或獲得預(yù)測(cè)結(jié)果。

本文所建模型由4層ConvLSTM2D網(wǎng)絡(luò)層+1層ConvLSTM3D網(wǎng)絡(luò)層堆疊而成，每層卷積核數(shù)目不同。4層ConvLSTM2D網(wǎng)絡(luò)層卷積核數(shù)目分別為40、40、60、40；卷積核大小均為3*3；補(bǔ)0策略為′same′，即輸出shape與輸入shape相同；使用線(xiàn)性激活函數(shù)；返回輸出序列中的全部序列。ConvLSTM3D網(wǎng)絡(luò)層卷積核數(shù)目分別為1；卷積核大小為3*3*3；補(bǔ)0策略為′same′，即輸出shape與輸入shape相同；使用′sigmoid′激活函數(shù)；輸出中維度的順序?yàn)橥ǖ涝诤竺妗?/p>

為了保證在線(xiàn)性向非線(xiàn)性轉(zhuǎn)變時(shí)，權(quán)重的尺度不變，所以使用BatchNormalization在激活函數(shù)前對(duì)輸入進(jìn)行了批標(biāo)準(zhǔn)化操作，優(yōu)化神經(jīng)網(wǎng)絡(luò)，將分散的數(shù)據(jù)統(tǒng)一。

由于sigmod函數(shù)具有梯度飽和現(xiàn)象的缺點(diǎn)，所以本實(shí)驗(yàn)使用的是學(xué)習(xí)率為0.01的隨機(jī)梯度下降優(yōu)化器optimizers.SGD，將所有參數(shù)梯度裁剪到數(shù)值范圍為[-0.5,0.5]的區(qū)間內(nèi)。

2.3 ConvLSTM模型訓(xùn)練與預(yù)測(cè)

（1） 模型訓(xùn)練

在訓(xùn)練模型之前，需要通過(guò)compile來(lái)對(duì)學(xué)習(xí)過(guò)程進(jìn)行配置。compile接收三個(gè)參數(shù)，分別為優(yōu)化器、損失函數(shù)和指標(biāo)列表。優(yōu)化器optimizer為一個(gè)已預(yù)定義的優(yōu)化器名，或者一個(gè)類(lèi)型為Optimizer的對(duì)象；損失函數(shù)loss為目標(biāo)函數(shù)，目的模型優(yōu)化，可以使用預(yù)定義的損失函數(shù)名或一個(gè)自定義的損失函數(shù)；指標(biāo)列表metrics的設(shè)置是針對(duì)分類(lèi)問(wèn)題，一般設(shè)metrics值為′accuracy′。指標(biāo)可以是一個(gè)預(yù)定義指標(biāo)的名字，也可以是一個(gè)用戶(hù)定制的函數(shù)。指標(biāo)函數(shù)應(yīng)該返回單個(gè)張量，或一個(gè)完成metric_name - ＞ metric_value映射的字典。

本模型的損失函數(shù)為對(duì)數(shù)損失函數(shù)′binary_crossentropy′，以adadelta為優(yōu)化方法進(jìn)行編譯。為了降低訓(xùn)練時(shí)間，在Keras中創(chuàng)建一個(gè)多GPU模型，設(shè)定4GPU并行處理。輸入的數(shù)據(jù)即為2.1中處理的BASIC_SEQUENCE和NEXT_SEQUENCE作為訓(xùn)練集和結(jié)果集序列；batch_size值為10，即計(jì)算梯度所需的樣本數(shù)量為10；epochs值為10，即代表樣本集內(nèi)所有的數(shù)據(jù)經(jīng)過(guò)了10次訓(xùn)練；validation_split=0.05，即樣本集中5%的數(shù)據(jù)為驗(yàn)證集。訓(xùn)練樣本集，取其中一回合的部分訓(xùn)練結(jié)果如圖10所示。

（2）預(yù)測(cè)

模型訓(xùn)練完成后存入“nice_model.h5”文件中，取樣本集中第600組圖片中的12幀進(jìn)行預(yù)測(cè)，預(yù)測(cè)16次，訓(xùn)練結(jié)果連接形成預(yù)測(cè)集。

（3）可視化結(jié)果

將預(yù)測(cè)結(jié)果與訓(xùn)練所用樣本集及其對(duì)應(yīng)的結(jié)果集對(duì)比形成可視化結(jié)果。結(jié)果分為16組圖片，前8組圖片顯示初始軌跡與標(biāo)記數(shù)據(jù)的對(duì)比結(jié)果，如圖11所示為第4幅圖片。

后8組圖片顯示預(yù)測(cè)結(jié)果與標(biāo)記數(shù)據(jù)對(duì)應(yīng)結(jié)果的對(duì)比結(jié)果，如圖 12所示為第12幅圖片。

2.4 算法評(píng)價(jià)

從圖10所示的模型訓(xùn)練結(jié)果來(lái)看，損失函數(shù)趨于減小的趨勢(shì)，最終穩(wěn)定在一個(gè)較小的值，說(shuō)明模型的擬合程度較好。從圖11的初始軌跡對(duì)比來(lái)看，軌跡運(yùn)動(dòng)與標(biāo)記數(shù)據(jù)吻合度很好，從圖12的預(yù)測(cè)結(jié)果來(lái)看，預(yù)測(cè)圖像與標(biāo)記數(shù)據(jù)的結(jié)果集吻合度基本一致，說(shuō)明預(yù)測(cè)結(jié)果正確性好。

值得注意的是，模型圖像的處理對(duì)后期結(jié)果的影響很大，由于不同時(shí)段雷達(dá)回波反射率因子分布區(qū)間不同，在進(jìn)行灰度轉(zhuǎn)換時(shí)，若使用原始數(shù)據(jù)繪制圖像可能使比色卡（color map）的范圍不同，這將導(dǎo)致不同圖像中的相同強(qiáng)度的反射率因子值對(duì)應(yīng)的RGB色彩不同，若不進(jìn)行處理將極大影響訓(xùn)練效果。

另外，為了得到更為準(zhǔn)確的外推結(jié)果，可在模型訓(xùn)練完成后，進(jìn)行3次以上的外推，可以提高預(yù)測(cè)精度。

3 結(jié)語(yǔ)

本算法使用大同市歷史氣象雷達(dá)數(shù)據(jù)，重點(diǎn)介紹了ConvLSTM模型在雷達(dá)回波外推短期降水預(yù)測(cè)的應(yīng)用方法，并詳細(xì)介紹了ConvLSTM模型在雷達(dá)回波外推即時(shí)空序列分析建模的整個(gè)過(guò)程，同時(shí)對(duì)相應(yīng)的算法和整個(gè)建模流程給出了實(shí)現(xiàn)及結(jié)果說(shuō)明。本算法可以從以下方面進(jìn)行拓展思考：

（1）本實(shí)驗(yàn)在建模時(shí)使用的損失函數(shù)是對(duì)數(shù)損失函數(shù)，損失函數(shù)有很多種，實(shí)際應(yīng)用中均方誤差函數(shù)也是常用的一種損失函數(shù)。后續(xù)實(shí)驗(yàn)中應(yīng)考慮如何優(yōu)化損失函數(shù)，選擇使模型擬合性更好的損失函數(shù)。

（2）ConvLSTM是Xingjian Shi等提出的，后面他針對(duì)雷達(dá)圖的天氣預(yù)測(cè)又提出了TrajGRU，基于運(yùn)行軌跡對(duì)圖像做更精準(zhǔn)的捕捉。本實(shí)驗(yàn)如何使用TrajGRU模型進(jìn)行更好的預(yù)測(cè)，這些問(wèn)題都需要進(jìn)一步地進(jìn)行探討。