基于修正小波變換插值-TAN的雷達降水粒子分類

李 海, 白 錦, 孫 研, 任嘉偉

(中國民航大學天津市智能信號與圖像處理重點實驗室, 天津 300300)

0 引 言

降水粒子分類的研究在很多領域都有關鍵的指導意義。識別不同類型的降水粒子對研究粒子的微物理特性極其重要,降水粒子的精確分類對天氣預報、危險預警等也有極其重要的參考價值。雙偏振天氣雷達作為新型天氣雷達,常被用來探測降水目標,其與常規天氣雷達相比,可以獲得除降水粒子回波強度外的其他各種偏振參量。不同的偏振參量可以反映不同的降水粒子信息,利用這些信息可以更好地實現降水粒子的檢測。

高分辨率的雙偏振氣象雷達獲取的目標回波對天氣過程的描述更加精細,可以獲得更準確的氣象數據,故可以提高降水粒子的探測精度。要提高雙偏振氣象雷達距離分辨率和方位分辨率,可以采用擴大天線直徑或者是擴大發射信號帶寬的方法,但此類做法需要對雷達天線、發射機等硬件進行改造,由此會帶來改造成本高、改造時間長等問題。目前，地基氣象雷達在遠距離工作時會因波束展寬導致獲得分辨率低的雷達數據;星載氣象雷達雖可提供較大范圍的降水觀測,且不易受波束阻擋的影響,但其時空分辨率較低;因此探究如何在數據分辨率比較低的情況下提高降水粒子探測精度依然是值得研究的課題。

針對低分辨率數據,常用的插值算法包括線性插值算法與非線性插值算法。線性插值方法通常會對分布不平滑的數據有均衡化影響,會對數據高頻信息造成損失,而非線性插值方法在一定程度上會改善這種缺陷,這當中比較多見的非線性插值方法是基于小波變換的插值方法。目前插值方法中基于小波變換的插值方法在圖像插值方面應用頗多,但很少涉及雷達數據處理中的應用。

基于統計決策的降水粒子分類方法、基于模糊邏輯的降水粒子分類方法、基于支持向量機等的降水粒子分類方法是目前常見的應用于各個領域的降水粒子分類方法。其中模糊邏輯方法是當前被廣泛使用的分類方法,而基于機器學習的降水粒子分類方法是近幾年的研究熱點。利用隸屬函數來綜合評定偏振參量是模糊邏輯法的特點,這種方法可以避免分類閾值的設定,故比傳統統計模式識別方法的分類精度高。但此方法通常將權威研究者在各種不同情形下總結的經驗值作為模糊邏輯隸屬函數的相關參數,所以當環境變化時分類結果會有較大的不穩定性。在現有的眾多分類算法中,貝葉斯算法因為能綜合各種先驗信息和有效利用數據樣本信息成為當下機器學習領域的研究熱點之一。何偉等應用樸素貝葉斯算法于預測降雨量,模型預測精度與目前的短期氣候預測精度比有明顯提高。郭雅芬等運用貝葉斯分類法預測中尺度對流系統的移動路徑。劉虹利等利用貝葉斯方法對降水空間的插值模型進行改進,以提升插值精度。

針對低分辨率雷達偏振參量數據的降水粒子分類問題,本文提出一種雙偏振氣象雷達體制下基于修正小波變換插值-樹擴展樸素貝葉斯(tree augmented naive Bayesian, TAN)的降水粒子分類方法,這種方法可以提高數據分辨率,也可以做到對各種不同類型的降水粒子進行分類。首先需要通過修正小波變換插值算法對原始低分辨率雷達偏振參量數據進行插值處理以獲得插值后的高分辨率偏振參量數據;其次針對得到的高分辨率偏振參量數據再做進一步離散化操作,并利用得到的離散化數據進行TAN網絡的結構和參數訓練以獲得最佳的TAN降水粒子分類網絡;最后,使用TAN分類網絡對修正小波變換插值后的高分辨率雷達偏振數據進行降水粒子分類。此方法在數據插值時可以保留更多高頻信息,且在構造TAN網絡時是通過對數據樣本進行訓練的方式得到網絡參數的,所以即使數據采集條件發生變化，也可以重新對數據進行相應的操作,故比模糊邏輯方法的分類效果更穩定。

1 修正小波變換插值-TAN算法

本文所提修正小波變換插值-TAN算法的整體思路為:首先將修正小波變換插值應用到低分辨率雷達數據以獲得高分辨率的雷達數據,然后利用獲得的高分辨率雷達數據進行TAN網絡結構和參數訓練,最后將訓練好的TAN網絡用于高分辨率雷達數據的降水粒子分類。下面對所提方法進行具體論述。

1.1 修正小波變換插值方法

本文中修正小波變換插值方法首先對原低分辨率偏振參量數據進行雙線性插值得到雙線性插值后高分辨率偏振參量數據。然后對進行小波變換,得到其對應的低頻分量數據與高頻分量數據。這時, 低頻分量數據近似等于原低分辨率偏振參量數據,則由除以可求得修正參數矩陣。考慮到各個分解部分之間具有相似性,利用求得的修正參數矩陣對高頻分量數據做同樣的修正,可得到修正高頻分量數據。最后由原低分辨率偏振參量數據替代修正低頻分量數據,與修正高頻分量數據一起經過小波逆變換得到所需的修正小波變換插值后的高分辨率偏振參量數據。該方法的示意圖如圖1所示。

圖1 修正小波變換插值方法示意圖Fig.1 Schematic diagram of modified wavelet transform interpolation method

圖1中,符號÷為除法符號,表示(,)除以(,),即用數據矩陣的第行第列元素除以數據矩陣的第行第列元素,相應得到數據矩陣的第行第列元素;符號×為乘法符號,表示(,)乘以(,),即用數據矩陣的第行第列元素乘以數據矩陣的第行第列元素,最后得到數據矩陣的第行第列元素。

下面以雷達偏振參量中的反射率因子為例(其他偏振參量同理),對修正小波變換插值方法具體過程進行介紹。

(1) 計算低頻分量數據矩陣_high_(×)與高頻分量數據矩陣_high_(×)

假設反射率因子的低分辨數據矩陣為行列的矩陣_low(×),對其進行雙線性插值可得到插值后高分辨率數據矩陣_high(×2)。

根據小波變換的mallat算法進行信號的分解與重構。對數據矩陣_high(×2)按行進行小波變換,()表示數據矩陣_high(×2)的第一行數據, 其信號分解公式如下：

(1)

(2)

式中：+1()和+1()分別表示離散信號()經分解后得到的低頻成分與高頻成分；()和()分別表示低通濾波器和高通濾波器對應的濾波系數。 按照式(1)和式(2)對_high(×2)的每一行數據進行信號分解,每一行的低頻成分組合即可得到低頻分量數據矩陣_high_(×),每一行的高頻成分組合即可得到高頻分量數據矩陣_high_(×)。

(2) 計算修正參數矩陣_(×)

修正參數矩陣_(×)實際上是由雙線性插值后高分辨率數據的低頻分量數據矩陣_high_(×)計算原低分辨率偏振參量數據矩陣_low(×)時所需的系數矩陣。修正參數矩陣_(×)可用下式進行計算:

(3)

(3) 計算修正高頻分量數據矩陣_high_(×)

因為小波分解過程中各個分解部分間具有相似性,根據此性質,利用修正參數矩陣_(×)與雙線性插值后高分辨率偏振參量數據的高頻分量數據矩陣_high_(×)可以計算得到修正高頻分量數據矩陣_high_(×),計算公式如下:

_high_(,)=_(,)×_high_(,)，=1,2,…,;=1,2,…,

(4)

(4) 計算修正小波變換插值后的高分辨率反射率因子數據矩陣_new(×2)

根據小波變換mallat算法的信號重構公式,如下所示：

(5)

對原低分辨率偏振參量數據矩陣_low(×)和修正高頻分量數據矩陣_high_(×)按行進行小波逆變換,即可得到最終修正小波變換插值后的高分辨率反射率因子數據矩陣_new(×2)。其中，+1()和+1()分別表示數據矩陣_low(×)和數據矩陣_high_(×)的第1行的數據;()表示小波逆變換得到的數據矩陣_new(×2)的第1行數據,對其他行數據同樣按照式(5)進行小波逆變換即可得到數據矩陣_new(×2)。

1.2 TAN降水粒子分類方法

針對修正小波變換插值后的高分辨率降水粒子雷達回波數據,本文采用TAN網絡進行降水粒子分類。為了得到離散的雷達偏振參量數據集,首先對高分辨率的雷達偏振參量數據進行離散化操作,然后以降水粒子類型變量作為類別屬性,偏振參量,,作為條件屬性搭建TAN網絡框架,接著用離散化的高分辨率雷達偏振參量數據集進行TAN網絡結構和參數訓練以完成TAN降水粒子分類網絡的構建,最后用TAN網絡實現降水粒子的分類。

1.2.1 數據離散化處理

TAN降水粒子分類方法需要數據的屬性是離散化的,否則學習準確率比較低,因此本文在構建離散屬性TAN網絡前首先對雷達偏振參量數據做相應的離散化處理。關于數據離散化處理的各類方法中,等寬算法較為常見,本文即采用等寬區間離散化方法對雷達偏振參量數據進行離散化處理以得到所需的離散數據集。

以條件屬性為例,對第21節修正小波變換插值后得到的高分辨率數據集_new進行離散化處理以得到離散數據集。

首先，生成離散化標準。假設該數據集離散化前的數據總量為,利用史特吉斯公式可以計算離散區間數:

=1+332lg

(6)

假設該數據集中最大的屬性值為_max,最小的屬性值為_min,根據式(7)可計算得到斷點間隔:

=(_max-_min)

(7)

接著根據()=_min+(=1,2，…,)求斷點,得到的斷點集為=[(1),(2),…,(),…,()]，此處得到的斷點集稱為離散化標準。

其次,生成離散數據集。離散數據集是按最小歐式距離準則,根據上述離散化標準生成的。方法如下:依次取出高分辨率數據集_new中的各個數據,計算當前取出的數據_new()(=1,2,…,)與離散化標準(=[(1),(2),…,(),…,()])中每一個元素(∈)的歐式距離=[,1,,2,…,,],找到集合中最小歐式距離min,,(,∈),確定與min,對應的斷點()(()∈),該斷點即原數據經離散化后得到的數據值。將該數據放入離散數據集中,令()=(),作為離散數據集中的第個元素,直至遍歷原數據集_new中的所有元素即可得到最終的離散數據集。

按照上述方法對偏振參量和做同樣的處理,得到和的離散數據集為與,離散化標準為(1×)與(1×)。

122 TAN網絡結構訓練

對雷達偏振參量數據做離散化處理后,接下來進行TAN網絡結構訓練。首先需要確定TAN網絡的類節點和屬性節點,本文以降水粒子類型變量作為TAN網絡的類節點,屬性節點由偏振參量屬性變量，,來充當,以此構造初始TAN網絡結構。接著根據互信息理論繼續確定各屬性節點間的依賴關系并通過添加擴展弧的方式對初始TAN網絡結構進行優化。最后得到與樣本數據集匹配度最好的TAN網絡結構。初始TAN網絡結構圖如圖2所示。

圖2 初始TAN網絡結構圖Fig.2 Initial TAN network structure diagram

圖2中,圓形稱為網絡的節點,網絡的邊是圖中的有向線段,邊所指向的節點稱為子節點,邊的出發端節點是父節點,父節點的數據值會對子節點的數據取值有一定的影響。

接下來根據互信息理論對TAN網絡結構繼續進行優化。在圖2的基礎上通過計算屬性節點,,兩兩之間的互信息并與給定互信息門限比較,從而繼續確定各偏振參量屬性之間的依賴關系。

以反射率因子和差分反射率為例,首先計算與之間的互信息(,)。該互信息可以通過下式得到:

(,)=

(8)

式中：,為反射率因子數據和差分反射率數據的個數;(_,_)表示偏振參量屬性(,)的狀態為(_,dr_)時的聯合概率;(_)和(dr_)分別表示偏振參量=_時的邊緣概率以及=dr_時的邊緣概率,這些概率均可通過對數據樣本進行統計得到。

接下來將互信息(,)與互信息門限(一般取001≤002)比較,若(,)≤,則認為偏振參量與之間不存在關聯性;若(,)>,則表明偏振參量與之間存在關聯性,此時要在屬性節點與間新增一條邊,暫不定向。接著，根據式(8)計算屬性節點和與類節點之間的互信息(,)與(,),比較(,)與(,),若(,)>(,),則新增邊的箭頭由屬性節點指向屬性節點;若(,)<(,),則新增邊的箭頭由屬性節點指向屬性節點;若(,)=(,),則隨機定向。最后更換屬性節點,用上述同樣的方法確定其他屬性節點兩兩之間的依賴關系,并通過加帶箭頭邊的方式體現在TAN網絡結構中。TAN網絡結構的具體訓練過程如下。

構建初始TAN網絡結構,初始TAN網絡結構的屬性節點為雷達偏振參量,,,類節點為降水粒子類型變量。

給定互信息門限(一般取001≤≤002)。

計算屬性節點與之間的互信息并與門限比較,若(,)>則執行步驟4;若(,)≤則執行步驟5。

在屬性節點對和之間建立一條無向邊,然后根據屬性節點和與類節點之間互信息值的大小來確定無向邊的方向,由互信息值大的屬性節點指向互信息值小的屬性節點,若兩個互信息值一樣大小則屬性節點間無向邊的方向可隨機指向任一節點。

更換屬性節點并重新執行步驟3,直到兩兩屬性節點對均遍歷完為止。

以上TAN網絡結構訓練的操作流程圖如圖3所示。

圖3 TAN網絡結構訓練操作流程圖Fig.3 Flow chart of TAN network structure training operation

為了得到最終的TAN分類網絡結構圖,將各離散化高分辨率雷達偏振參量數據集按照上述流程圖進行TAN網絡結構訓練,最終得到的結構圖如圖4所示。

圖4 TAN降水粒子分類網絡結構圖Fig.4 Structure diagram of TAN hydrometeor classification network

1.2.3 TAN網絡參數訓練

離散屬性TAN網絡使用條件概率表描述雷達偏振參量對不同降水粒子的分布,條件概率表為一種分布列,為了得到各個節點的條件概率表,故進行相應的TAN網絡參數訓練。

圖4的TAN網絡結構圖中,屬性節點的條件概率與節點和的數值都有關系,因為屬性節點的父節點有兩個,即和。根據第221節得到的偏振參量的離散化標準的維數,以及偏振參量的離散化標準的維數,同時結合類標簽=(=1,2,…,9)(代表9種粒子類型)的維數可以確定屬性節點的條件概率表為××9的三維分布列。

將第221節經過離散化處理后得到的各偏振參量的離散數據集作為訓練數據集,用()表示訓練數據集的樣本個數,將訓練數據集中滿足=()(=1,2,…,),=()(=1,2,…,),=(=1,2,…,9)的樣本個數記為(()∩()∩),那么條件概率表中條件概率計算方式如下:

(9)

在訓練數據集中可能存在(()∩()∩)=0這種情形,因為后續分類運算涉及到乘法運算,所以這種情況會對分類結果產生消極影響,故通過拉普拉斯平滑對式(9)做相應的處理,具體操作如下所示:

(10)

以屬性節點為例,說明參數訓練過程,如圖5為的條件概率表示意圖,其中每個單元格表示對應各雷達偏振參量取值的概率。

圖5 ZDR條件概率表示意圖Fig.5 ZDR conditional probability indication intention

假設訓練數據集中樣本個數為,滿足圖5中標出單元格=-0410 5,=8489 7,=的樣本個數為,則該單元內概率計算方式如下所示：

(11)

根據這種計算方式,計算出每個雷達偏振參量節點的條件概率表,可完成參數訓練。

124 TAN網絡的降水粒子分類原理

通過上述TAN網絡結構訓練和參數訓練構建好TAN網絡后,下面對TAN網絡的具體分類原理進行論述。

用(=1,2,…,9)表示類節點的各種可能取值,即各類降水粒子類型,其中表示地雜波(將地雜波也看作一種降水粒子)、表示冰晶、表示干雪、表示濕雪、表示雨、表示暴雨、表示大雨滴、表示霰、表示雨夾雹。根據貝葉斯分類算法的基本原理,若以上TAN網絡中3種偏振參量,,的取值分別為,,,則分類結果為的概率(|,,)可表示為

(12)

式中：()代表的是第類降水粒子的先驗概率,在未知的情況下一般采用等概率,即()=1(=1,2,…,9);(,,)表示偏振參量取值為,,的數據出現的概率,在樣本確定的情況下該概率值是固定的;(,,|)稱為類條件概率,根據訓練得到的TAN網絡結構,該概率值可計算如下:

(,,|)=(|)(|,)(|)，=1,2,…,9

(13)

式中：(|)表示在約束條件為時取值為的數據出現的概率;(|,)與(|)同理。此處的概率值通過對TAN網絡進行參數訓練記錄在條件概率表中。

最終的分類結果可以計算如下:

(14)

式(14)所表達的含義是:計算偏振參量值,,對應每種分類結果(=1,2,…,9)的概率(|,,),取其中概率值最大的作為分類結果輸出。

2 算法流程與步驟

本文中雙偏振天氣雷達體制下的修正小波變換插值-TAN的降水分類方法詳細操作步驟如下。

基于修正小波變換插值方法對低分辨率雷達偏振參量數據進行插值處理以獲得高分辨率雷達偏振參量數據。

在步驟1得到的對高分辨率數據的基礎上實行等距離散化操作,最終得到所需的離散數據集。

構造初始TAN網絡結構,初始TAN網絡結構的屬性節點由各偏振參量變量充當,類節點由降水粒子類別變量充當。

引入互信息理論,使用離散化雷達偏振參量數據對初始TAN網絡結構重新進行結構訓練,接著進行參數訓練,獲得最終的TAN網絡。

使用最終的TAN網絡對離散的高分辨率雷達偏振參量數據進行降水粒子分類。

綜上所述,雙偏振天氣雷達體制下修正小波變換插值-TAN的降水粒子分類方法整體流程圖如圖6所示。

圖6 基于修正小波變換插值-TAN的降水粒子分類方法流程圖Fig.6 Flow chart of hydrometeor classification method based on modified wavelet transform interpolation-TAN

3 實驗結果及分析

為了驗證所提方法的性能,本文選取美國國家海洋和大氣管理局(National Oceanic and Atmospheric Administration, NOAA)數據庫中的KTLX雷達提供的偏振參量數據進行實驗。KTLX雷達為美國WSR-88D雷達網中位于俄克拉荷馬市(35.195 8 N°, 97.164 0 W°)的一部雙偏振氣象雷達。這部雷達的波束寬度是1.25°,天線增益是45 dB,脈沖重復頻率是250～1 200 Hz,帶寬是0.3 MHz,該雷達的回波分辨率為360×300。首先基于修正小波變換法對獲得的回波數據做兩次插值的操作,接下來再做降水粒子分類,將得到的降水粒子分類結果與美國國家海洋和大氣管理局官網提供的高分辨率數據(360×1 200)的分類結果進行對比,計算分類誤差。為了說明所提方法的可靠性,此處選取兩組數據進行實驗。

根據圖7和圖8,TAN方法對插值后數據的分類結果接近于NOAA提供的標準數據分類結果,在此基礎上對分類結果誤差進行了統計,經統計分類誤差均處于較低水平,2017年為1.37%,2018年為0.92%。根據上述實驗分析可以看出，本文方法在數據量僅為標準數據25%的條件下能夠實現準確的數據插值處理結果以及降水粒子分類。

圖7 2017年8月17日標準數據與插值重構后數據的 TAN降水粒子分類結果對比圖Fig.7 Comparison diagram of precipitation particle classification results between standard data and interpolated reconstructed data on August 17, 2017

圖8 2018年6月24日標準數據與插值重構后數據的 TAN降水粒子分類結果對比圖Fig.8 Comparison diagram of precipitation particle classification results between standard data and interpolated reconstructed data on June 24, 2018

為了進一步說明TAN降水粒子分類方法的有效性,仍以上述兩組數據為樣本進行降水粒子的分類研究,圖9和圖10分別給出了使用模糊邏輯算法和TAN算法得到的降水粒子分類結果以及NOAA給出的降水粒子分類結果。

圖9 2017年8月17日數據的不同降水粒子分類方法結果對比圖Fig.9 Comparison diagram of the results of different precipitation particle classification methods based on the data on August 17, 2017

分析圖9與圖10,對在不同時間段獲取的雷達偏振參量數據進行降水粒子分類,模糊邏輯法的分類結果均存在一定誤差。圖9(a)中,除了沒有準確識別地雜波外,對其他降水粒子也沒有做到準確分類,將很多非雨的降水類型誤識別為雨。而圖10(a)中,模糊邏輯方法對濕雪和霰完全未能識別出來。而本文所提的TAN降水粒子分類方法針對不同時間段獲得的數據分類效果均較好,與NOAA提供的結果更為相近。綜上可知,對不同時間段的降水粒子進行分類,TAN方法的分類效果均比較穩定,不會產生較大誤差。

圖10 2018年6月24日數據的不同降水粒子分類方法結果對比圖Fig.10 Comparison diagram of the results of different precipitation particle classification methods based on the data on June 24, 2018

4 結 論

對于雙偏振雷達體制下數據分辨率較低時的降水粒子分類問題,本文給出了基于修正小波變換插值-TAN的雷達降水粒子分類方法。該方法首先對原始低分辨率雷達偏振參量數據進行修正小波變換插值以得到插值后高分辨率雷達偏振參量數據。接著，把雷達偏振參量變量當作屬性節點,把表示降水粒子類型的變量作為類節點去構造初始TAN網絡結構,并在此基礎上引入互信息理論,利用離散化的高分辨率偏振參量數據重新進行TAN網絡結構和參數訓練以完成最優TAN網絡的構建。最后,將插值后高分辨率偏振參量數據帶入到訓練好的TAN網絡中實現降水粒子的分類。該方法在數據插值處理中可以更多地保留原始數據的高頻信息,另外通過數據樣本訓練方式計算TAN網絡分類參數的方法比以經驗值作為分類參數的傳統模糊邏輯方法的分類穩定性更好。上述實驗結果表明,該方法可以獲得較為準確的數據插值與降水粒子分類結果。