閆 軍, 王新舒, 韓旭日
(內(nèi)蒙古自治區(qū)氣象科學(xué)研究所, 內(nèi)蒙古 呼和浩特 010051)
雷暴等強對流天氣是世界公認的嚴重威脅航空器飛行安全的天氣系統(tǒng),該類系統(tǒng)所產(chǎn)生的風切變、強降水、冰雹等天氣現(xiàn)象是造成飛行事故的直接原因。因此,一直以來強對流都是航空氣象的研究重點,同時也是難點。一方面,針對高頻次發(fā)生的強對流天氣,國內(nèi)外學(xué)者雖然已經(jīng)進行了長期的研究,但是對其形成演變的機理認知十分有限,系統(tǒng)的復(fù)雜性使得諸多問題依然沒有被解決;另一方面,社會及行業(yè)重點領(lǐng)域?qū)υ擁棙I(yè)務(wù)服務(wù)的要求越來越高,預(yù)報預(yù)警的準確率、時空分辨率、時效性需求極高,航空安全在其中尤為突出。因此,如何利用新的技術(shù)方法快速、準確地獲取強對流天氣影響范圍、精度及程度是臨近預(yù)報業(yè)務(wù)的研究熱點,得到的研究成果在航空氣象勤務(wù)中能發(fā)揮有效作用[1]。
層狀云與對流云自動識別算法的研究一直是天氣雷達應(yīng)用研究領(lǐng)域的熱點,并且已經(jīng)產(chǎn)生了一批切實可行的成果。識別算法主要分為背景場閾值技術(shù)(BET)和概率匹配方法兩大類。D.D.Churchill 等人在雷達反射率數(shù)據(jù)基礎(chǔ)上使用BET 方法,該方法集中在核的識別上,為每一個已識別核施加一個固定的影響半徑,以此確定對流云區(qū)域[2]。M.Steiner 等人認為文獻[2]所使用的固定影響半徑依據(jù)是不充分的,他們將影響半徑作為變量,它是背景場平均反射率的函數(shù),設(shè)定的門限閾值取決于背景場平均反射率,該算法已被TRMM(the Tropical Rainfall Measuring Mission)項目采納為降水分類的業(yè)務(wù)運行算法[3]。DeMot 等人將上述方法擴展到三維,通過使用采自雙極化雷達的垂直速度,證明了三維分析改進了降水分類的精度[4]。Houze 等人也將他們的研究建立在三維分析的基礎(chǔ)上,但不同于其他的BET 方法,他們將重點放在通過窗區(qū)概率匹配方法來改進降水估測,將層狀云和對流云不同降水增長機制的影響與雷達數(shù)據(jù)固有的限制相結(jié)合,由3 個參數(shù)代表降水分割結(jié)果:徑向反射率梯度、亮帶小數(shù)、有效性(測量云的深度)[5]。國內(nèi)從事這一領(lǐng)域研究的人員相對較少。王改利在Lakshmanan 的對雷暴識別聚類分析方法基礎(chǔ)上提出一種K-means 方法,該方法基于K-means 的分層聚類能夠?qū)崿F(xiàn)暴雨回波的多尺度識別。肖艷姣選取了組合反射率及其水平梯度、反射率因子等于35 dBZ 的回波頂高及其水平梯度、垂直累積液態(tài)水含量及其密度等6 個候選參數(shù),并統(tǒng)計這6 個候選識別參數(shù)分布的概率密度特征,利用模糊邏輯法進行層狀云和對流云的識別[6-7]。
本研究是在充分汲取前人成果的基礎(chǔ)上,嘗試使用深度學(xué)習方法來試驗對流云區(qū)和層狀云區(qū)的識別[8-10]。深度學(xué)習的技術(shù)理論已經(jīng)在多個領(lǐng)域展示出較高的應(yīng)用價值,并且不斷在計算機視覺、自然語言處理、自動駕駛等多個方面取得重大突破。其算法模型的建立依賴大量的樣本數(shù)據(jù),天氣雷達具備數(shù)量龐大的高時空分辨率觀測數(shù)據(jù),這就為將深度學(xué)習理論引入天氣系統(tǒng)識別和預(yù)測提供了可能。在眾多的機器學(xué)習研究方向中,表征學(xué)習的重點是如何自動找出表示數(shù)據(jù)的合適方式,以便能更好地將輸入變換為輸出;而本研究所采用的深度學(xué)習理論是具有多級表示的表征學(xué)方法,它可以逐級表示越來越抽象的概念或模型。在該模型的框架下,構(gòu)建數(shù)據(jù)集合和特征標簽,以此來訓(xùn)練端到端的多層神經(jīng)網(wǎng)絡(luò),雷達數(shù)據(jù)中的對流云區(qū)可以逐級表征為特定的位置和角度的邊緣,即通過具備一定模型深度的多個隱藏層的神經(jīng)網(wǎng)絡(luò),逐層組合特征非線性變換,將原始特征映射到高維度的特征空間,從而語義分割出對流云區(qū)。
訓(xùn)練、驗證和測試數(shù)據(jù)集的建立采用呼和浩特新一代天氣雷達2015—2022 年(呼和浩特C 波段HHHTRC)基數(shù)據(jù),該雷達饋源海拔高度為2 062 m,預(yù)警業(yè)務(wù)觀測使用VCP21 模式(9 層立體掃描),數(shù)據(jù)流間隔6 min,數(shù)據(jù)庫長為250 m。雷達站所處位置凈空條件良好,無明顯遮擋及干擾源,數(shù)據(jù)質(zhì)量較高。
1) 依據(jù)雷達觀測記錄挑選出較為顯著的對流性天氣過程的基數(shù)據(jù);
2) 由于雷達基數(shù)據(jù)獲取過程中使用的是徑向極坐標采集和存儲的方式,加之地球曲面的變化,需要對數(shù)據(jù)做等距保角投影轉(zhuǎn)換,使之投影到與地球相平行的平面上[6]。因此,需要將極坐標系數(shù)據(jù)通過坐標變換的方式轉(zhuǎn)換成地理經(jīng)緯坐標系三維數(shù)據(jù)場,坐標轉(zhuǎn)換公式為:
式中:λ0、φ0為雷達天線所在位置的經(jīng)緯度;λ、φ為目標物的經(jīng)緯度;L為目標物到天線之間的距離;R為等效地球半徑;h為天線到達地面的高度;H為回波所在高度。
3) 原始基數(shù)據(jù)是沿雷達掃描線徑向分布的,為保持探測數(shù)據(jù)的連續(xù)性,相鄰掃描線之間和相鄰掃描層之間存在數(shù)據(jù)點,需要填充以免造成信息的缺失。本文采用四線插值算法,該算法在確定掃描盲區(qū)的一個空間點物理量數(shù)值時需4 條雷達掃描線上的8 個極坐標點。其算法為:將坐標系中的空間某一數(shù)據(jù)點(xn,yn,zn),通過公式轉(zhuǎn)換到極坐標系中去,設(shè)它在極坐標系中的相應(yīng)位置為(Ri,θj,φk),對其取整數(shù),從而獲得與它靠近的8 個極坐標數(shù)據(jù)點(Ri1,θj1,φk1)、(Ri2,θj1,φk1)、(Ri1,θj2,φk1)、(Ri2,θj2,φk1)、(Ri1,θj1,φk2)、(Ri2,θj1,φk2)、(Ri1,θj2,φk2)、(Ri2,θj2,φk2)。所對應(yīng)的回波數(shù)據(jù)如圖1 所示。
圖1 本研究邏輯結(jié)構(gòu)
根據(jù)(Ri,θj,φk)與相鄰8 個數(shù)據(jù)點的距離遠近,對回波數(shù)據(jù)進行四線性插值,得出(Ri,θj,φk)點上的數(shù)據(jù),這也就是(xn,yn,zn)上的數(shù)據(jù);最后,綜合雷達數(shù)據(jù)和深度學(xué)習算法的可行性,將經(jīng)上述處理后雷達數(shù)據(jù)場數(shù)據(jù)轉(zhuǎn)化為垂直最大回波強度顯示分析產(chǎn)品(CR)建立模型數(shù)據(jù)集,CR 產(chǎn)品利用數(shù)據(jù)場強度數(shù)據(jù),在以1 km×1 km(或2 km×2 km)為底面積的數(shù)據(jù)點垂直向上到回波頂?shù)拇怪敝w中,對所有位于該柱體中的回波強度數(shù)據(jù)點進行比較,挑選出最大回波強度。再用測高公式計算最大回波強度的所有高度,從而得到最大回波強度分布圖像。這種產(chǎn)品有助于用戶快速查看最大回波強度的分布,對有些突發(fā)性的強對流天氣,其初始回波往往出現(xiàn)在中空,用戶使用這種產(chǎn)品可以有效監(jiān)測這一類的強對流天氣,在穩(wěn)定性降水條件下,還有助于用戶識別0 ℃層亮帶。由于在降雹區(qū)域,相應(yīng)的中空可能存在水分累積區(qū),所以CR 產(chǎn)品還可以作為監(jiān)測冰雹的發(fā)生發(fā)展的工具之一。
主要選取呼和浩特雷達站汛期(5 月—10 月)預(yù)警區(qū)域內(nèi)所發(fā)生的強對流天氣過程的探測基數(shù)據(jù),經(jīng)預(yù)處理后組成算法數(shù)據(jù)集,其中231 個數(shù)據(jù)樣本作為訓(xùn)練集,69 個數(shù)據(jù)樣本作為驗證集,85 個數(shù)據(jù)樣本作為測試集,比例設(shè)置為3∶1∶1。數(shù)據(jù)集的時間跨度為2015—2022 年,其中2022 年的數(shù)據(jù)作為測試集樣本,在實驗過程中2022 年的樣本不參與模型訓(xùn)練和調(diào)參,這有利于客觀地驗證模型的泛化和遷移能力以及識別效果,如表1 所示。
表1 數(shù)據(jù)集詳細信息
基于雷達強度數(shù)據(jù)對流云區(qū)識別算法模型的建立可以充分借鑒計算機視覺中圖像語義分割的理論,后者著重像素級的圖像目標識別,即通過算法標注出每個像素所屬的對象類別,而對流云區(qū)的識別實質(zhì)上是通過強度數(shù)據(jù)點的識別劃分出不同目標區(qū)域。
鑒于雷達數(shù)據(jù)量的原因,模型設(shè)計之初應(yīng)該考慮到在此種分割場景下算法要足夠輕量。一般的分割模型網(wǎng)絡(luò)結(jié)構(gòu)通常:輸出是一個N×H×W的概率圖,其中N是一個向量,向量的維數(shù)代表類別的數(shù)量,即向量的每個元素對應(yīng)一個類別,而這個元素的數(shù)值是概率值,表示對應(yīng)識別類別的概率,給出每一個數(shù)值點類別,所以H×W大小要與原數(shù)據(jù)場保持一致。分割網(wǎng)絡(luò)實質(zhì)上是一個沒有全連接層的全卷積網(wǎng)絡(luò),該類型網(wǎng)絡(luò)結(jié)構(gòu)可以分為Encoder 和Decoder 兩個部分。Encoder 部分特征圖的長寬會逐步縮小,需要下采樣;Decoder 部分特征圖的長寬會逐步增大,需要上采樣。
一個切實可行的分割網(wǎng)絡(luò)需要解決如下問題:
1) 通過盡可能少的Encoder 層數(shù)獲得較大的感受野。所謂感受野是指輸出的二維數(shù)據(jù)場上每個點的數(shù)值,是由輸入二維數(shù)據(jù)場上大小為kh×kw的區(qū)域的元素與卷積核每個元素相乘再相加得到的,所以輸入數(shù)據(jù)場區(qū)域內(nèi)每個元素數(shù)值的改變都會影響輸出點的像素值,將這個區(qū)域叫作輸出數(shù)據(jù)場上對應(yīng)點的感受野。感受野的大小是評價分割網(wǎng)絡(luò)的特征表示能力的重要指標之一。
2) 下采樣過程中雖然得到了高階特征,但其長和寬在不斷縮小,如果直接作為上采樣輸入,則會丟失信息。解決這一問題的方法就是通過跳躍連接,使得上采樣過程不僅使用高階特征,也可以使用同維度的低階特征。
3) 多尺度問題。在雷達探測資料中,對流性云區(qū)目標有大有小,如何保證不同尺度的目標均能被很好的處理呢?解決的方式是不同尺度的目標交給不同感受野的卷積層來處理。
本研究借鑒分割理論在生物醫(yī)學(xué)領(lǐng)域的應(yīng)用成果,其中U-Net 網(wǎng)絡(luò)是采用編碼器(下采樣)-解碼器(上采樣)結(jié)構(gòu)和跳躍連接的一種設(shè)計方法,也是比較基礎(chǔ)的方法,后續(xù)的許多卷積神經(jīng)網(wǎng)絡(luò)仍延續(xù)了該算法的核心思想。U-Net 網(wǎng)絡(luò)如圖2 所示,其基本特點是使用全卷積網(wǎng)絡(luò)(FCN)。它與卷積神經(jīng)網(wǎng)絡(luò)(CNN)的區(qū)別在于將CNN 最后的全連接層(FC)替換為卷積層(Conv),因此該網(wǎng)絡(luò)是一個端到端的網(wǎng)絡(luò);網(wǎng)絡(luò)的Contracting path(如 圖2 左 側(cè))使 用Conv 和MaxPooling;網(wǎng) 絡(luò) 的Expansive path(如圖2 右側(cè))使用上采樣(upsampling)與左側(cè)的Contracting path、poolingde featuremap 相結(jié)合,然后逐層上采樣到heatmap;最后經(jīng)過兩次Conv,達到最終的heatmap,再用1×1 的Conv 做分類。
圖2 U-Net 網(wǎng)絡(luò)結(jié)構(gòu)圖
對流云和層狀云無論從外部形態(tài)還是內(nèi)部機理都存在很大的區(qū)別:對流云是由大氣層結(jié)構(gòu)不穩(wěn)定引起的氣流垂直運動所致,在雷達探測分析產(chǎn)品中對流性降水表現(xiàn)為較高的反射率數(shù)值,云體水平反射率梯度較大,云區(qū)形態(tài)分布面積較小,較大的垂直厚度,云頂部不平整,云頂高度梯度大;層狀云在雷達探測分析產(chǎn)品中所反映的反射率強度相對較弱,分布較為均勻,水平反射率梯度較小,其外部形態(tài)表現(xiàn)為云區(qū)分布面積較大,比較薄的垂直厚度, 頂部比較平整, 云頂高度梯度不大,常出現(xiàn)零度層亮帶等特征。根據(jù)上述云體特征與模型技術(shù)要求進行算法實驗。
整個模型網(wǎng)絡(luò)有23 層卷積層,其中左側(cè)為收縮網(wǎng)絡(luò),右側(cè)是擴張網(wǎng)絡(luò)。收縮網(wǎng)絡(luò)由多組卷積神經(jīng)網(wǎng)絡(luò)組成,每組由兩層 3×3 unpadded 卷積層組成,每個卷積層后使用ReLU。每組之間使用size=2、stride=2 的最大池化進行下采樣。每組第一次卷積時,將特征圖通道數(shù)增加1 倍。擴張網(wǎng)絡(luò)也是由多組卷積層組成,每組先使用2×2 的up-convolution(轉(zhuǎn)置卷積)將分辨率翻倍且通道數(shù)減半,然后在與收縮網(wǎng)絡(luò)中對應(yīng)的層,輸出特征結(jié)果concatenation。但是因為U-Net 全部使用unpadded 卷積,對應(yīng)的收縮網(wǎng)絡(luò)的特征圖的分辨率更高。由于算法實驗中采用了Overlap-tile 策略,可以將收縮網(wǎng)絡(luò)對應(yīng)的特征結(jié)果四周裁剪后,再拼接到擴張網(wǎng)絡(luò)特征圖中,拼接后接兩個unpadde 的3×3 卷積層,每個卷積后接ReLU,擴張網(wǎng)絡(luò)的最后一層是1×1 卷積層,將通道數(shù)轉(zhuǎn)換為類別數(shù)。
為了減小開銷,U-Net 網(wǎng)絡(luò)在訓(xùn)練參數(shù)上將batch size 設(shè)置為1,所使用的SGD 的動量系數(shù)為0.99。加權(quán)損失部分最后一層使用逐強度點的softmax 和交叉熵損失,公式如下:
式中:ak(X)是二維數(shù)據(jù)場中數(shù)據(jù)點的位置X處通道k的激活值;K是通道總數(shù)(也就是類別數(shù));?∈{1 ,2,…,K},表示每一個數(shù)據(jù)點的類別;ω(X)是X數(shù)據(jù)點的損失權(quán)重。
實驗過程中還涉及數(shù)據(jù)增強,這是由于樣本數(shù)量相對較少,采用的方法是使用3×3 組網(wǎng)格上的隨機唯一向量生成平滑形變。
算法模型把對流云區(qū)識別分割后,需要將模型預(yù)測的結(jié)果與人工標注的驗證標簽進行比對,來衡量分割算法的優(yōu)劣。本研究采用該領(lǐng)域通用的評價指標:準確率(Accuracy)、召回率(Recall)、Dice 系數(shù)。
準確率(Accuracy)為在雷達分析產(chǎn)品強度數(shù)據(jù)總數(shù)中被正確分類的數(shù)據(jù)點數(shù),該指標通常被用來衡量圖像分類網(wǎng)絡(luò)模型對像素分類的性能;召回率(Recall)是真實類別為對流云區(qū)的數(shù)據(jù)點被網(wǎng)絡(luò)模型正確分類的概率;Dice 系數(shù)通常在醫(yī)學(xué)圖像分割領(lǐng)域用作評價指標,該指標用于對網(wǎng)絡(luò)模型預(yù)測值與真實標注的相似度進行評估。Dice 系數(shù)取值范圍為[0,1],其評估的分值越高,說明該模型的分割效果越好,用于本研究則表示算法識別結(jié)果與真實標注越接近。準確率(Accuracy)、召回率(Recall)、Dice 系數(shù)計算公式分別為:
式中:TP 為測試集中對流云區(qū)真實的數(shù)據(jù)點總數(shù)(被標記為對流云區(qū)的數(shù)據(jù)點,并被正確預(yù)測為對流云區(qū)的數(shù)據(jù)點);TN 為測試集中非對流云區(qū)的數(shù)據(jù)點總數(shù)(被標記為非對流云區(qū)的數(shù)據(jù)點,并被正確預(yù)測為非對流云區(qū)的數(shù)據(jù)點);FP 為測試集中對流云區(qū)真實的數(shù)據(jù)點總數(shù)(被標記為非對流云區(qū)的數(shù)據(jù)點,但被預(yù)測為對流云區(qū)的數(shù)據(jù)點);FN 為測試集中對流云區(qū)真實的數(shù)據(jù)點總數(shù)(被標記為對流云區(qū)的數(shù)據(jù)點,但被預(yù)測為非對流云區(qū)的數(shù)據(jù)點)[10-11]。
對實驗的測試結(jié)果采用上述評價指標進行計算,最終統(tǒng)計評價指標,如表2 所示。由表2 可知:準確率平均值為76.63%,最高值為82.02%;召回率與準確率相近,分別為76.34%和81.89%。
表2 實驗結(jié)果統(tǒng)計 %
Dice 系數(shù)則反映出一個分割效果較好的數(shù)值,從圖中總體識別效果看:
1) 算法能夠準確地識別并分割出對流云區(qū)域,特別是相對面積較大的區(qū)域(如圖3a)所示);
圖3 識別結(jié)果
2) 對于正在發(fā)展,強度不是很強的對流性云區(qū)域識別率不高,未能識別出整個大的區(qū)域,將整片區(qū)域分割為若干區(qū)域(如圖3b)所示);
3) 盡管大的區(qū)域捕捉較為準確,但區(qū)域的邊緣識別不夠精準,邊緣識別結(jié)果較為粗糙;
4) 面積小的對流性云區(qū)域識別度較低,這與雷達數(shù)據(jù)本身存在噪聲點、雜波及地物干擾項有關(guān)(如圖3c)、圖3d)方框區(qū)域所示)。圖中,白色區(qū)域為識別出新一代天氣雷達對流云區(qū)大致位置。
本文嘗試將卷積神經(jīng)網(wǎng)絡(luò)模型用于天氣雷達探測特征天氣系統(tǒng)領(lǐng)域,解決對流云區(qū)域的識別及分割問題。在數(shù)據(jù)集的建立過程中用到了坐標轉(zhuǎn)換、四線插值和組合反射率等算法,并建立了U-Net 網(wǎng)絡(luò)架構(gòu)的算法模型。經(jīng)過訓(xùn)練和調(diào)參的實驗過程,從而找到識別效果最佳的模型網(wǎng)絡(luò)組成。最后利用準確率、召回率和Dice系數(shù)等評價指標建立識別分割效果的統(tǒng)計數(shù)值支撐。從研究成果看,總體達到了設(shè)計之初的實驗效果;但從識別結(jié)果中也發(fā)現(xiàn)了諸多問題,在后續(xù)的研究中仍有較大的提升空間。