基于遺傳算法的冰晶冰水含量估計方法

李 海 孫 研 孟凡旺 白 錦 韓雁飛

(1.中國民航大學天津市智能信號與圖像處理重點實驗室 天津 300300；2.中國航空工業集團公司雷華電子技術研究所 江蘇無錫 214063；3.中國民用航空華北地區空中交通管理局天津分局 天津 300300)

0 引言

飛機在巡航階段(一般在22000ft以上)飛行時，高空冰晶是對飛機飛行安全威脅較大的天氣現象。在冰晶高聚集區域，發動機葉片表面的冰晶堆積會造成發動機失速與喘振，燃燒室熄火和損壞等。同時，冰晶還會進入大氣數據傳感器，造成堵塞使得相關傳感器數據異常，造成飛機失控。近年來因遭遇高空冰晶而發生的航空飛行事故時有發生。美國聯邦航空管理局對近二十年來發生的發動機空中事故進行統計，統計結果表明這些事故大部分與高空冰晶相關。這些真實的事件表明高空冰晶對航空運輸業和旅客飛行安全會產生嚴重影響，開展高空冰晶檢測的研究是十分重要的，而高空冰晶冰水含量估計是整個檢測流程的基礎，直接影響后續危險級別判定的準確性。因此針對高空冰晶冰水含量估計的研究具備重要的意義。

機載氣象雷達能探測飛機飛行前方空域的危險氣象條件，可以用作飛行預警或規避危險氣象的重要參考。一般情況下，機載氣象雷達探測降水目標時是根據雷達回波的強度來區分的，并且認為回波強度越大，含水成分越多，對飛行安全的影響越大。而高空冰晶檢測及其危險級別量化可以通過所探測空間的冰水含量(Ice Water Content，IWC)，即單位立體空間的冰晶密度予以指示。IWC是指單位體積空氣中包含的固態水的質量，單位為g/m。

國外學者提出可根據雷達反射率因子Z與IWC之間存在的統計關系，結合經驗公式通過雷達反射率因子反演估計IWC值。1954年空中劍橋研究中心的Atlas根據飛機實測的譜參數數據，獲得了雷達反射率因子-粒子分布-IWC之間的經驗關系，1992年美國的Liao等研究了Ka波段雷達和W波段雷達的Z-IWC關系，1995年Brown用毫米波機載雷達探測了冰晶云的Z-IWC經驗公式，2001年Austin等人通過采用數值模擬的方法，提出了一種基于毫米波雷達反射率和可見光層云微物理參數反演冰水含量的方法，2004年Matrosov等人研究了基于雷達反射率因子的估測海洋性層云含水量的方法，并比較了取不同的降水云與非降水云強度界限時云內含水量的變化，2015年大氣介質空間觀測實驗室的Protat分別給出了35GHz 和95GHz雷達下的Z-IWC關系，并統計對比了實驗誤差。但由于經驗公式形式過于簡單具有較大局限性，難以準確表征復雜的Z-IWC關系模型，特別是使用不同數據進行冰水含量估計時會存在大誤差，魯棒性較差，因此如何獲得精確的Z-IWC關系模型是進行高空冰晶冰水含量準確估計亟待解決的問題。

針對數學模型相對復雜，所需確定參數較多的Z-IWC關系模型，利用一些智能算法進行自適應的參數估計從而得到準確的模型是解決問題的新方案。遺傳算法是一種在處理各種復雜優化問題時效果顯著，且具備較強魯棒性的算法。從1995年開始，國內對遺傳算法的研究也日益增加，作為啟發式搜索算法的一員，遺傳算法在最近幾十年內在信號處理、最優控制、生產規劃等各個領域都有著長足的進步與發展。遺傳算法在求解函數優化問題時，相比其它優化算法，最大的特點是處理對象不是參數本身，而是對參數集編碼獲得基因個體，因此該算法不會受到函數連續性的限制。另外該算法不需要問題的梯度信息等輔助信息，僅根據適應度函數的數值大小就可以選出最優個體，且適應度函數也不會受到連續可微的條件限制，所以適用范圍更廣。針對高空冰晶冰水含量估計問題，如能設計遺傳算法進行參數估計，發揮其魯棒性較好的優點，確定較為精確的Z-IWC關系模型將具有很好的應用價值。

基于此，本文提出了一種基于遺傳算法的機載氣象雷達高空冰晶冰水含量估計方法。該方法對雷達反射率因子與高空冰晶冰水含量關系模型中待定參數進行精細尋優得到最優解，然后使用估計的參數得到確定的模型后利用測試數據的雷達反射率因子進行反演，估計高空冰晶的冰水含量。針對2016年3月17日美國俄克拉荷馬州地區一次真實氣象場景的數值天氣預報模式輸出數據進行實驗，利用遺傳算法估計冰水含量，對比遺傳算法與經驗公式法冰水含量估計結果，并給出不同時間點、不同高度層的冰水含量估計結果，與真實值的對比結果表明，該方法可有效地估計高空冰晶的冰水含量，估計結果較為準確，統計分析的指標參量(包括相關性系數、RMSE、誤差統計占比等)較為理想。

1 雷達反射率因子與冰水含量關系模型

雷達反射率因子與粒子直徑的6次方成正比，一般用來反映單位體積內云、冰晶、雨滴等粒子群的濃度情況。冰水含量(IWC)反映了單位立體空間的冰晶密度。對雷達反射率因子與冰水含量關系進行建模是利用雷達反射率因子反演冰水含量的關鍵。

根據Heymsfield等人的研究，冰晶粒子群的等效雷達反射率因子和冰水含量表達式分別如公式(1)和公式(2)所示。

(1)

(2)

由于公式(1)與公式(2)均含有斜率參量成分，因此可利用共有成分將兩公式進行聯立從而得到Z-IWC關系模型。對公式(1)進行變換，提出斜率參量成分可得

(3)

對公式(2)進行變換，提出斜率參量成分可得

(4)

將公式(3)與公式(4)進行聯立可得到

(5)

對公式(5)進行變換可得

(6)

將確定的、、等參數帶入公式(6)可得

(7)

公式(7)即為求得的Z-IWC關系模型，模型中的待定參數為、、。在Z-IWC模型中，只要知道了、、三個參數，就可以得到準確的Z-IWC模型，即可以通過雷達反射率因子準確反演冰水含量。其中，是冰晶質量冪律分布函數(D)中表征冰晶粒子成分的參數；是冰晶粒子數濃度分布函數()中表征冰晶粒子尺寸的參數；是冰晶粒子數濃度分布函數()中決定分布形狀的參數。

2 基于遺傳算法的Z-IWC參數估計

本文使用遺傳算法對Z-IWC模型中的三個未知參數、、進行參數估計，以確定具體的Z-IWC關系模型，從而完成冰水含量的估計。遺傳算法的核心部分包括：

1)種群初始化與編碼；

2)適應度函數的選取；

3)遺傳選擇；

4)交叉和變異；

5)迭代終止。

利用遺傳算法進行參數估計的具體實現流程如圖1所示。下面對利用遺傳算法進行Z-IWC模型參數估計的主要操作過程進行具體論述。

圖1 基于遺傳算法的Z-IWC模型參數估計流程圖

1)編碼與種群初始化。遺傳算法進行參數估計，首先對要優化問題的解進行編碼，編碼的目的主要是為了改變優化問題解的表現形式，使其方便遺傳運算。本文采用二進制編碼的方式對Z-IWC關系模型中的三個參數、、分別進行編碼，將每一個染色體(解)編碼為基因型。完成編碼操作后，接著需要設置初始種群。對于模型中每一個待定參數，本文通過隨機生成的形式，來創造一個初始種群，文中是隨機生成具有均勻分布的200個個體，作為初始種群。

2)適應度函數的選取。遺傳算法的收斂速度以及能否找到最優解，直接依賴于適應度函數的選取。為了使估計的冰水含量值更加接近冰水含量真值，本文構造滿足下列優化準則的目標函數(,,)作為適應度函數如式(8)所示。

(8)

式(8)中，IWC(,,)為利用氣象場景數據空氣密度、水汽混合比等參數計算出來的冰水含量值，作為IWC的真實值，IWC′(,,)為利用式(7)的Z-IWC關系模型求得的冰水含量估計值，為使用的樣本個數。

4)遺傳交叉與變異。遺傳算法的收斂性，主要取決于變異和交叉兩大核心操作。變異的實質是在群體中發掘個體的多樣性，提高算法的局部隨機搜索能力，而交叉操作決定產生新個體的速度，對實現算法的全局搜索性能具有關鍵性的作用。為了更好地完成交叉和變異操作，以盡可能得到最優解，本文設置交叉概率為06，設置變異概率為008。

5)迭代終止。本文設置進化代數計數器=0，設置最大進化代數=100，運行終止條件設置為若優化準則目標函數的目標值滿足()≤=001，或者運行迭代次數=100次，則停止遺傳進化過程，選取最后一代的最優個體作為最優解，以得到確定的Z-IWC模型。利用確定的模型對雷達數據進行處理，得到高空冰晶冰水含量估計結果。

3 冰水含量估計方法流程與步驟

基于遺傳算法的機載氣象雷達高空冰晶冰水含量估計方法具體實現步驟如下：

1)步驟1：對雷達數據進行質量控制；

2)步驟2：建立高空冰晶冰水含量與雷達反射率關系模型(Z-IWC關系模型)；

3)步驟3：利用遺傳算法對Z-IWC關系模型中的未知參數進行精細尋優；

4)步驟4：根據確定后的參數得到最終確定的Z-IWC關系模型；

5)步驟5：利用確定的Z-IWC關系模型對雷達數據進行計算，反演IWC得到高空冰晶冰水含量估計結果。

總體流程框圖如圖2所示。

圖2 總體流程圖

4 實驗及結果分析

為了驗證所提方法性能，本文將2016年3月17日美國俄克拉荷馬州地區西經95°至100°，北緯30°至34°范圍內的一次真實氣象場景的數值天氣預報模式輸出數據作為樣本進行檢測分析，分別選取俄克拉荷馬州地區2016年3月17日14時與2016年3月17日16時兩個時間點的數據進行實驗。

為了驗證所提方法的有效性和準確性，針對兩個不同時刻點的不同高度層的數據分別給出了本文方法進行冰水含量估計與經典經驗公式法進行冰水含量估計的實驗結果對比圖，如圖3和圖4所示。本文方法與1954年Atlas等人提出的經驗公式IWC=0.064Z，1997年Aydin等人提出的經驗公式IWC=0.104Z以及2002年Mace等人提出的經驗公式IWC=0.104Z進行對比。從圖3和圖4的實驗結果中可以看出，經驗公式法因受經驗取值的影響，冰水含量估計結果均產生了不同程度的偏差，誤差較大，而采用遺傳算法完成冰水含量估計，誤差較小，具有更優越的魯棒性。

圖3 2016年3月17日14時不同高度層Z-IWC參數估計結果對比圖

圖4 2016年3月17日16時不同高度層Z-IWC參數估計結果對比圖

從真值與估值相關性的角度，下面給出了俄克拉荷馬州地區2016年3月17日14時與2016年3月17日16時真值IWC與估值IWC的相關性散點圖，分別如圖5和圖6所示。圖中橫坐標是利用遺傳算法估計得到的IWC值，縱坐標是IWC真值，當估計值完全等于真實值時數據點應該全部分布在對角線上，那么越靠近對角線則說明估計的越準確，從圖5和圖6可以看出數據點都在對角線附近，說明誤差較小，估值與真值之間具有很強的相關性，證明了本文方法進行高空冰晶冰水含量估計的準確性。

圖5 2016年3月17日14時不同高度層真值與估值冰水含量相關性散點圖

圖6 2016年3月17日16時不同高度層真值與估值冰水含量相關性散點圖

在給出真值IWC與估值IWC相關性散點圖的同時，對利用本文方法估計得到的IWC值的誤差大小及出現的概率情況進行統計，分別給出2016年3月17日14時和2016年3月17日16時不同高度層數據的誤差統計直方圖如下圖7和圖8所示。圖7和圖8中的誤差為估計的IWC值減去IWC真值的絕對誤差，由圖可以看出兩個時刻點不同高度層數據的真值與估計值的誤差出現在0附近的概率較大，均占到80%左右，故驗證了算法的準確性。

圖7 2016年3月17日14時不同高度層誤差分布直方圖

圖8 2016年3月17日16時不同高度層誤差分布直方圖

另外，分別統計本文方法與經驗公式法進行冰水含量估計時估值與真值之間的相關性系數與均方根誤差(Root Mean Square Error，RMSE)，2016年3月17日14時數據的統計結果如表1和表2所示，2016年3月17日16時數據的統計結果如表3和表4所示。由此可見本文方法與經驗公式法相比，相關性系數更接近1，表明具有很強的相關性。同時本文方法的RMSE都小于0.1，說明估計值與真值之間偏差較小，而經驗公式法的RMSE偏差較大，對比表明本文算法準確性更高。

表1 IWC估值與真值之間的相關性系數(2016年3月17日14時數據)

表2 IWC估值與真值之間的RMSE(2016年3月17日14時數據)

表3 IWC估值與真值之間的相關性系數(2016年3月17日16時數據)

表4 IWC估值與真值之間的RMSE(2016年3月17日16時數據)

5 結束語

本文提出了一種基于遺傳算法的機載氣象雷達高空冰晶冰水含量估計方法。首先對高空冰晶冰水含量與雷達反射率因子關系進行建模，得到理論推導的基礎模型之后，利用遺傳算法對模型中待定參數進行精細尋優獲得最優解，然后利用估計的參數得到確定的模型后對測試數據的雷達反射率因子進行反演，估計高空冰晶的冰水含量。通過對仿真數據進行試驗驗證，利用遺傳算法估計冰水含量，對比遺傳算法與經驗公式法冰水含量估計結果，并給出不同時間點、不同高度層的冰水含量估計結果，與真實值的對比結果表明，該方法可有效地估計高空冰晶的冰水含量，估計結果較為準確，統計分析的指標參量(包括相關性系數、RMSE、誤差統計占比等)較為理想。