雙流機場雷暴天氣預報方法研究

李典南，徐 海，許東蓓

(1.貴州省氣象臺，貴州 貴陽 550002；2.中國民用航空西南地區空中交通管理局氣象中心，四川 成都 610202；3.成都信息工程大學大氣科學學院高原大氣與環境四川省重點實驗室，四川 成都 610225)

0 引 言

雷暴是一種常見的帶有雷電現象的局地性強對流天氣[1]，其水平發展的范圍為幾千米至幾十千米，垂直發展的高度可達8～15 km，持續時間為幾分鐘到幾小時不等[2]。雷暴的發生通常伴有雷暴大風、短時強降水、冰雹甚至龍卷等一種或多種天氣現象[3]，給航空飛行帶來了十分惡劣的影響。據統計，全球每天約有4.4萬個雷暴發生，并且在任意一個時刻都有0.2～0.4萬個雷暴在活動，其活動面積占比全球面積高達1%[4]，嚴重威脅航空飛行的安全。近年來國內外相關統計顯示，每年由于氣象因素造成的飛行事故超過事故總數的1/3[5]，而在這部分由于氣象因素造成的飛行事故中，又有至少1/3是受到雷暴天氣的影響[6-8]。由此可見，提高雷暴天氣的預報水平，可在一定程度上減少航空飛行事故的發生以及降低飛行事故帶來的損失。

雷暴的形成通常需要3個必要條件：深厚且明顯不穩定的氣層、充沛的水汽以及合適的觸發條件[9]。在雷暴發生發展過程中，大氣各氣象要素在不斷發生變化，當某些氣象要素的強度達到一定程度時，雷暴便會發生。為了表征某個時刻大氣的性質特征，可選用某些參數的數值來衡量，這就引入了對流參數的概念：將基于氣塊法分析對流條件而得出的參數，統稱為對流參數。在實際預報中，可依據對流參數的數值特征來判斷大氣層結的不穩定性、水汽條件和動力條件等的強弱變化情況，并基于此得出未來強對流天氣的發生時間、落區和強度等信息[10]。利用對流參數進行預報的方法有許多，比如指標疊套法[11]，即選取多個因子并計算其指標，通過對這些指標進行分析疊套來進行潛勢預報。還可將多個因子相結合為一個新的預報因子，比如將垂直螺旋度與對流有效位能相結合而成的“垂直能量螺旋度指數(VEH)”[12-13]，VEH能有效體現強對流天氣的動力因子和熱力因子的共同貢獻，楊晶軼等[14]通過實際驗證，發現VEH對雙流機場的雷暴天氣的預報具有一定指示性。此外，建立潛勢預報模型也是一種常用的預報方法，比如：主分量旋轉法、二級邏輯回歸法、Bayes判別法、Logistic回歸判別法、神經網絡法等[15-20]。其大體思路是：篩選出與雷暴發生關系最密切的幾個對流參數或物理量，再將這些參數或物理量利用不同方法來建立預報模型。

成都雙流國際機場位于103°57′02″E、30°34′47″N，地處四川盆地西部的平原腹地，海拔約為504.3 m，是中國大陸的八大航空樞紐之一，是西南地區重要的客貨集散地。雙流機場地勢西高東低，自西北向東南傾斜[21]，受青藏高原大地形和盆地氣候特征的影響，容易在近地面形成較強的垂直風切變以及氣流輻合線、低渦等中小尺度天氣系統[22]。在夏季，雙流機場常處于太平洋副熱帶高壓的西南邊緣，加之機場東面河流縱橫、水汽充沛[23-25]，使得雷暴天氣頻發。

為研究雙流機場的雷暴天氣并為其空中管制提供技術支持，本文參考了常用預報方法的思路來建立潛勢預報模型，并在其基礎上做了一定改進。全文綜合利用了多種氣象資料，首先運用統計學方法分析雙流機場雷暴天氣的時間變化特征，再利用相關性分析篩選出對雙流機場雷暴天氣指示性較好的物理量因子以及主要雷達特征量，最后在此基礎上基于二級邏輯回歸法建立潛勢預報模型。

1 資料與方法

1.1 資料數據

所用的資料主要為：①2013—2018年雙流機場的逐小時觀測資料。該資料由雙流機場例行每小時記錄一次，要素包括風向、風速、氣溫、露點、相對濕度、修正海平面氣壓、云、能見度、天氣現象等等。②歐洲中期天氣預報中心(ECMWF)2013—2018年的ERA-interim逐6 h再分析資料(水平分辨率為0.125°×0.125°)，主要為高空資料。③成都市氣象局多普勒天氣雷達S波段SC型號的雷達產品資料，主要包含回波頂高、1.5°仰角基本反射率、3.4°仰角基本反射率、反射率最大回波出現的高度、組合反射率、垂直累積液態水含量和3 h累積降水量，所選取的數值是以雙流機場為中心、半徑20 km范圍內的最大值。

1.2 研究方法

1.2.1 相關系數的計算 相關系數是最早由統計學家卡爾·皮爾遜設計的統計指標，是研究變量之間線性相關程度的量，一般用字母r表示。記某一天雷暴過程的發生與否為Y，發生為Y=1，不發生為Y=0；記某因子為X，X的樣本序列為Xi(i=1,2,3,…)。為判斷X因子與雷暴發生與否Y的相關性，可運用點雙序列公式[26]來計算相關系數：

(1)

1.2.2 多元回歸方程的建立 利用二級邏輯回歸的思路建立多元回歸方程組(第一級預報方程、第二級消空方程)[27]，該方法的基本思路是把一個預報事件與多個預報因子均看成隨機事件，將隨機事件Y的出現用“1”記錄，未出現用“0”記錄，將事件Y與n個預報因子之間的關系看成是n個因子已經出現的條件下事件Y出現與否的關系。

在邏輯回歸中，X和Y均為“0”“1”化后的矩陣，具體方法為：X表示因子與閾值之間的關系，當因子的值大于其閾值，則記為X=1，反之記為X=0；Y表示當日有無雷暴過程，若當日有雷暴過程發生則記為Y=1，反之記為Y=0。通過擬合計算，可建立如下形式的多元回歸方程(b0,b1,b2, ……,bn為回歸參數)：

Y=b0+b1×X1+b2×X2+b3×X3…+bn×Xn

(2)

1.2.3 預報因子的檢驗方法 在確定好預報因子和預報模型之后，將一份樣本數據帶入檢驗，以評估其預報效果，具體步驟為：

首先，將預報因子X1、X2、X3……Xn的指標(閾值)帶入對應的第一級預報方程，用k1代表計算出的Y值；再將樣本數據帶入同一方程，用Y(預報)代表利用樣本數據計算出的Y值。若有Y(預報)≥k1，則預報為有雷暴過程，反之為無雷暴過程。

第二步，將預報因子Xn+1、Xn+2、Xn+3……X2n的指標(閾值)帶入對應的第二級消空方程，用k2代表計算出的Y值，同樣再將樣本數據帶入第二步的方程得出Y(消空)，預報有無雷暴的方法與第一步相同。

第三步，結合第一級預報方程、第二級消空方程的預報結論來綜合判定：若某日Y(預報)預報為“無”，則預報結論為“無”。若某日Y(預報)預報為“有”，則進一步利用第二級消空方程進行判定——若Y(消空)也預報為“有”，則結論為“有”；若Y(消空)預報為“無”，則結論為“無”。

第四步，通過計算和統計，將預報結論帶入王名才等[28]公式來分析預報效果：

(3)

(4)

(5)

其中，N為實況發生雷暴的總次數，n為預報有雷暴的總次數，n1為正確預報出雷暴的次數，n′為預報有雷暴但實況無雷暴的次數。

2 雙流機場雷暴變化特征

根據雙流機場的逐小時天氣現象資料，將每一次雷暴天氣自開始時間至結束時間的整個過程記為1次(除1.2節外均按此標準統計)。通過統計(表1)，2013—2018年雙流機場共有221 d出現雷暴過程，累計發生283次。其中伴隨降水過程的雷暴共發生了187 d，累計218次，占雷暴總次數的77.03%；其余為干雷暴過程，共發生了62 d，累計65次，占雷暴總次數的22.97%。下面對這6 a的雷暴天氣進行月變化和日變化的特征討論。

表1 雙流機場2013—2018年雷暴發生天數和次數統計Tab.1 Statistics on the number of days and the number of thunderstorms at Shuangliu Airport from 2013 to 2018

2.1 月變化特征

2013—2018年的初雷在2月、3月各出現1次，4月出現4次；終雷在8月出現2次，9月出現3次，10月出現1次。圖1為雙流機場2013—2018年雷暴逐月發生次數統計圖，由圖可見雷暴的發生次數隨月變化呈顯著的單峰型：1—7月雷暴發生次數隨時間增加，7月雷暴發生的次數最多，8月雷暴發生的次數僅次于7月。經計算，每年7、8月雷暴發生的次數之和，均達到各年全年發生總次數的一半以上。從9月開始雷暴發生次數明顯減少，11月—次年1月基本無雷暴發生。可見春季(3—5月)和夏季(6—8月)是雷暴的多發季，其中以夏季7、8月發生次數最多。這是因為在春季和夏季，雙流機場邊界層附近的氣流受日照加熱抬升較強，并且機場東側水域豐富，近地面潮濕，尤其在炎熱的夏季，雙流機場處太平洋副熱帶高壓和青藏高壓的過渡區中，常受西南暖濕氣流和高原波動槽的影響，所以形成了高溫高濕且不穩定的環境，有利于雷暴天氣發生[29]。而秋季和冬季四川盆地層結較穩定，水汽條件和抬升條件較弱，雷暴發生條件不足，故發生次數較少。

圖1 雙流機場2013—2018年雷暴各年的逐月發生次數柱狀圖Fig.1 Histogram chart of the number of monthly thunderstorms occurrence in Shuangliu Airport from 2013 to 2018

2.2 日變化特征

將每小時作為一個時段，若某個時段有雷暴過程發生或有雷暴過程持續，則記該時段發生1次，統計結果如圖2所示。由圖可知，雷暴在上午發生的次數最少，6 a里平均每小時僅發生十幾次。午后雷暴發生次數逐漸增加，21時—次日06時(北京時，下同)是高發時段，6 a共計發生521次，占比總次數的58.54%，這期間平均每小時發生雷暴超過40次。若以每3 h為1個時段，則雷暴發生次數最多的時間段是21時—00時，共計發生193次；其次是00—03時169次，03—06時159次，18—21時98次。綜合來看，雷暴在一天之內的時間分布呈明顯的“夜雷多、日雷少”特征。

圖2 雙流機場2013—2018年雷暴在一天各時段內發生次數折線圖Fig.2 Line chart of the number of thunderstorms occurrence in Shuangliu Airport from 2013 to 2018 during various periods of the day

形成該特征的主要原因有：第一，雙流機場所處地形較閉塞，近地面潮濕，云霧較多，夜間云層上部輻射冷卻速度較快，云層下部由于云霧的保暖作用而冷卻較慢，氣層上冷下暖趨于不穩定，易形成雷暴天氣[30]。第二，夏季太平洋副熱帶高壓伸入我國西南地區且活動頻繁，太平洋副熱帶高壓白天多西伸，使雙流機場受西南暖濕氣流影響，為雷暴的發生發展預備了良好條件；夜間多東退，高原西風槽東移，多在后半夜移動到雙流機場附近引起雷暴天氣[31]。此外，弱冷鋒在經過四川盆地時，由于受到地形抬升和太陽輻射的共同作用，在白天常表現為鋒消；而夜間隨著地表輻射增強，又逐漸表現為鋒生，易形成鋒面雷暴[32]。

3 預報因子的篩選

3.1 物理量因子

基于雷暴形成的3個必要條件[33]，初步選取了850 hPa假相當位溫、850與500 hPa假相當位溫差、相對濕度、比濕等物理量因子以及對流有效位能、K指數、強天氣威脅指數、熱力總指數、抬升指數、風切變指數等對流參數(以下統稱為物理量因子)。借助歐洲中心ERA-interim逐6 h再分析資料，針對雙流機場2013—2018年的雷暴天氣進行物理量因子的相關系數計算(公式(1))。計算樣本為2013—2018年每年的4—8月(2014年8月資料不全，故排除，下同)，共計887 d，其中有208 d發生雷暴天氣。物理量因子Xi的取值為每日的算術平均值。

如表2所示，將通過α=0.05顯著性水平檢驗的因子按相關系數的絕對值從大到小排列，相關系數越大說明因子對雷暴的預報指示性越好。由表2可知，對流有效位能、K指數、850與500 hPa假相當位溫差、850 hPa比濕、850 hPa假相當位溫、強天氣威脅指數對雷暴有較高的指示意義，相關系數絕對值均高于0.4。除此之外，700 hPa比濕、熱力總指數、抬升指數與雷暴的相關性也較好，相關系數絕對值均高于0.3。

表2 物理量因子與雷暴發生情況的相關系數Tab.2 The correlation coefficient between the physical factors and thunderstorms

統計學中，箱線圖(又稱箱圖)能以簡單的組合圖形直觀地反映較大容量樣本的值分布情況和數據批中隱含的結構信息，可通過分析數據批的分布形狀、排除異常數據點，來進行數據批之間的比較[34]。利用相關系數絕對值大于0.3的前9個物理量因子的日均值，繪制雷暴日、非雷暴日的箱線圖(圖3)。由圖可知，每組因子箱線圖的箱體之間存在一定交集(值域交集)，計算各因子箱體交叉部分占雷暴日箱體的比例，得出的結論按從小到達排列為：對流有效位能0%；K指數0%；抬升指數6.56%；850與500 hPa的假相當位溫差10.49%；850 hPa比濕15.27%；850 hPa假相當位溫16.18%；強天氣威脅指數20.21%；熱力總指數29.07%；700 hPa比濕40.94%。當某因子在雷暴日和非雷暴日的值域交集越少，說明該因子對雷暴的預報指示性越好，即對流有效位能、K指數、850與500 hPa的假相當位溫差和抬升指數對雷暴的預報更有利。

圖3 2013—2018年887個樣本有無雷暴情況下各物理量因子箱線圖(上下*號：第99百分位數和第1百分位數；上下短實線：第90百分位數和第10百分位數；矩形內實線：中位數；矩形上下邊線：第75百分位數和第25百分位數)Fig.3 Box plot of each physical factors in 887 samples with or without thunderstorms from 2013 to 2018(Upper and lower *: the 99th percentile and the 1st percentile; upper and lower short solid lines: the 90th percentile and the 10th percentile; solid line within the rectangle: median; upper and lower edges of the rectangle: 75th percentile and 25th percentile)

在樣本容量固定的情況下，因子過多不僅對擬合方程起的貢獻不大，因子本身存在的隨機因素還容易對方程的穩定性產生影響[19]。基于雷暴發生的3個基本條件以及相關系數和箱線圖所顯示結果綜合考慮，最后確定出如下4個物理量預報因子：對流有效位能、K指數、850 hPa比濕、850與500 hPa的假相當位溫差。

3.2 雷達因子

多普勒天氣雷達是監測實況天氣的重要手段，廣泛運用于強對流天氣的監測和預警。許多雷達產品都可用于反映強對流天氣的發展狀態，比如回波頂高、垂直累積液態水含量常被用來評估風暴強度[35]，雷達回波強度、回波頂高等常用來反映閃電和雷暴云發展的高度[36-37]等等。

結合實際業務經驗[37-42]，初步選取了與強對流天氣聯系較緊密的7個雷達因子進行相關系數的計算(表3)。計算樣本為2013—2018年每年4—8月的雷暴日，共計208 d，雷達因子值取逐小時整點時刻的值，共計4 992 h，其中發生雷暴的時間有873 h。表3中雷達因子已按相關系數的絕對值從大到小排列，且均通過α=0.05的顯著性水平檢驗。

表3 多普勒天氣雷達因子與雷暴發生情況的相關系數Tab.3 The correlation coefficient between the doppler weather radar factors and thunderstorms

以相同方法繪制箱線圖(圖4)，并計算各雷達因子箱體交叉部分(值域交集)占雷暴日箱體的比例，結論如下：1.5°仰角基本反射率0%；垂直累積液態水含量0%；3 h累積降水量13.52%；3.4°仰角基本反射率15.79%；回波頂高20.00%；組合反射率52.17%；反射率最大回波出現的高度46.03%。綜合相關系數和箱線圖兩方面考慮，最終確定的4個雷達預報因子為：回波頂高、1.5°仰角基本反射率、3.4°仰角基本反射率、垂直累積液態水含量。

圖4 2013—2018年4992個樣本有無雷暴過程的情況下各雷達因子箱線圖(上下*號：第99百分位數和第1百分位數；上下短實線：第90百分位數和第10百分位數；矩形內實線：中位數；矩形上下邊線：第75百分位數和第25百分位數)Fig.4 Box plot of each radar factors in 887 samples with or without thunderstorms from 2013 to 2018(Upper and lower *: the 99th percentile and the 1st percentile; upper and lower short solid lines: the 90th percentile and the 10th percentile; solid line within the rectangle: median; upper and lower edges of the rectangle: 75th percentile and 25th percentile)

4 潛勢預報模型的建立

4.1 預報因子閾值的計算

利用雷暴日物理量預報因子月平均值以及雷暴發生的整點時刻的雷達預報因子月平均值，繪制如圖5所示折線圖。圖5a顯示，物理量預報因子的月平均值隨時間變化的趨勢基本一致：最小值出現在4、5月，最大值出現在7、8月，具有明顯的月變化特征或季節變化特征，說明物理量預報因子閾值的計算需要分春季、夏季考慮。而圖5b顯示，雷達預報因子的月平均值無明顯的時間特征，說明滿足雷暴發生條件的雷達預報因子受季節影響較小，故閾值的計算無需分季節考慮。

圖5 2013—2018年雙流機場雷暴日物理量預報因子(a)、雷達預報因子(b)的月平均值Fig.5 The monthly average of thunderstorm daily physical predictors (a) and radar predictors (b) of Shuangliu Airport from 2013 to 2018

將物理量預報因子按季節分為春夏兩組(4、5月為春季，6—8月為夏季)，雷達預報因子不分組；再將每組每個因子的值按從大到小排列，選取每個因子的第一、四分位值作為該因子的預報閾值[17]，即雙流機場雷暴天氣的預報指標，計算結果如表4所示。

表4 雙流機場雷暴天氣的預報因子閾值Tab.4 Thresholds of predictors for thunderstorms in Shuangliu Airport

4.2 潛勢預報方程的建立

參照1.2.2節，用X1、X2、X3、X4分別表示850 hPa比濕、K指數、850與500 hPa的假相當位溫差、對流有效位能，利用這些物理量預報因子建立如公式(2)所示的春季Y1、夏季Y2多元回歸預報方程，作為第一級預報方程：

Y1=0.012+0.100×X1+0.081×X2+0.180×X3+0.038×X4

(6)

Y2=0.052+0.142×X1+0.243×X2+0.200×X3+0.031×X4

(7)

若僅使用Y1、Y2方程進行預報，容易造成命中率較高、虛警率也較高的情況(已驗證，結果略)，故進一步利用雷達預報因子建立Y3方程作為第二級消空方程：

Y3=0.050+0.064×X5+0.163×X6+0.063×X7+0.093×X8

(8)

其中，X5、X6、X7、X8分別表示回波頂高、1.5°仰角基本反射率、3.4°仰角基本反射率、垂直累積液態水含量。最后再對求得的上述回歸方程Y1、Y2、Y3進行F檢驗[43]，結果表明這3個方程是顯著的。

4.3 預報效果的檢驗

用于檢驗的樣本為2018年4—8月的相關數據，利用1.2.3節的檢驗方法：首先，將表3中所示的物理量預報因子的閾值按季節帶入對應的第一級預報方程(公式(6)或公式(7))，再將雷達預報因子的閾值帶入第二級消空方程(公式(8))，用k1、k2、k3分別代表計算出的Y1、Y2、Y3的值，即各方程的閾值：k1=Y1=19.00、k2=Y2=85.32、k3=Y3=7.85。再將2018年4—8月每日的物理量預報因子日均值依季節帶入對應的第一級預報方程Y1、Y2，在春季，當Y1(預報)≥k1時，預報為有雷暴過程，反之無雷暴過程；同理，在夏季，當Y2(預報)≥k2時，預報為有雷暴過程，反之無雷暴過程。第三步，運用第二級消空方程，將樣本中每日逐小時的雷達預報因子值帶入Y3，當Y3(消空)≥k3時，預報為有雷暴過程，反之為無雷暴過程。第四步，結合兩級方程的結論來綜合判定：若Y1(或Y2)、Y3均預報為“有”，則最終結論為“有”，其他情況下最終結論均為“無”。最后，將最終預報結果帶入公式(3)、公式(4)和公式(5)進行計算，結論如表5所示。

表5 利用Y1、Y2、Y3方程預報2018年4—8月雷暴天氣Tab.5 The conclusion of using Y1, Y2, Y3 equations to forecast thunderstorm from April to August 2018

由表5可知，預報指標在春季的命中率、虛警率、臨界成功指數分別為100%、50.00%、50.00%，命中率雖然很高，但由于空報數較多，使得虛警率較高、臨界成功指數較低。指標在夏季的命中率、虛警率、臨界成功指數分別為80.49%、8.33%、75.00%，命中率雖不如春季，但因為空報數較少，所以虛警率較低、臨界成功指數較高。因此可以認為這套指標以及建立的兩級預報方程，對雙流機場雷暴天氣的預報具有一定指示意義，且綜合來看在夏季的預報效果更好。在實際運用中，可將第一級預報方程用作24 h潛勢預報，若預報結論顯示第2 d有雷暴天氣，則加強天氣監測，再配合第二級消空方程，以更準確地判斷臨近時間的天氣情況。

5 結論與展望

5.1 結論

本文利用雙流國際機場2013—2018年的逐小時氣象觀測資料、歐洲中心ERA-interim逐6 h再分析資料、多普勒天氣雷達產品資料，對雙流機場的雷暴天氣進行了時間特征分析，并選取預報因子建立潛勢預報模型，結果表明：

①雙流機場2013—2018年共有221 d發生雷暴天氣，累積發生283次，77.03%的雷暴過程伴隨降水。雷暴發生次數隨時間呈單峰型，春夏為多發季，其中以夏季7、8月發生的次數最多。雷暴多于午后發生，21時—次日06時是高發時段，呈現明顯的“夜雷多、日雷少”特征。

②依據雷暴發生的基本條件和實際業務經驗選取物理量和雷達參數，通過相關分析，篩選出預報因子：對流有效位能、K指數、850 hPa比濕、850與500 hPa假相當位溫差、回波頂高、1.5°仰角基本反射率、3.4°仰角基本反射率、垂直累積液態水含量。再基于二級邏輯回歸的思路建立預報方程和消空方程：利用物理量預報因子分春、夏兩季建立兩個第一級預報方程；利用雷達預報因子建立一個第二級消空方程。通過回代檢驗，結果表明，夏季的命中率低于春季，但因為夏季預報的空報數較少，所以虛警率和臨界成功指數相對于春季都表現較好。由此可認為所建立的潛勢預報模型對雙流機場雷暴天氣的預報具有一定指示性，且綜合來看在夏季的預報效果更好。

5.2 反思與展望

論文中仍有許多的內容需要完善和探討，首先是建立預報方程方面，使用的僅僅是二級邏輯回歸法，若選擇多種方法(比如神經元方法、SVM方法等)，就可通過對比檢驗得出最優方案，或形成集成預報。其次是方程的檢驗方面，使用的是2018年的數據進行回代，若資料充足，選取2019年或之后的資料進行檢驗，會更具有說服力。最后，在雷達資料選取上，由于雙流機場的雷達資料缺失較多，故本文使用的雷達資料是成都市氣象局的多普勒天氣雷達產品資料，但若采用雙流機場的雷達資料進行分析，得出的結論會更加有說服力。