梯度提升決策樹在雷達定量降水估測中的應用研究

陳訓來 陳元昭 趙春陽 張珂

（1 深圳市氣象局，深圳 518040；2 深圳南方強天氣研究重點實驗室，深圳 518040；3 阿里云計算有限公司，杭州 310099）

0 引言

定量降水估測（QPE）是現代天氣業務體系的重要組成部分，一直以來是業務預報的難點和重點，也是定量降水預報、強降水短臨預警的基礎，對城市內澇、地質災害精細化預報起重要支撐作用。隨著新一代多普勒天氣雷達觀測網的逐步建成和業務運行，應用多普勒天氣雷達聯合自動氣象站雨量開展定量降水估測已成為一種較為成熟可靠的實現手段[1]。

利用多普勒天氣雷達定量降水估測，相對于用自動氣象站雨量觀測，具有測量范圍廣、時空分辨率高、可及時取得大面積定量降水資料的優點。由于雷達回波（反射率因子）反映的是空中水滴、冰晶等對雷達電磁波的反射或散射，因此雷達回波強度（Z）與降水強度（I）之間存在密切的指數型關系，即Z=A×Ib[2]。通過Z與I之間的對應關系就可以估測出雷達覆蓋范圍內的降水強度。降水定量估測精準度在很大程度上取決于Z-I關系式中系數A和b，降雨常具有時空特性強、變化快的特點，且不同季節、不同地點的降水性質不同，其Z-I關系式也不同，常用來建立Z-I關系的方法有最優化法[3]、概率密度法[4-5]、動態關系法[6]、遺傳算法[7]和動態分級法[8]等。為了得到更高精度的降水估計值，通常利用自動氣象站雨量資料對雷達降水估測值做進一步訂正，包括應用平均校準法[9]、變分同化法[10]、客觀插值分析法[11]和卡爾曼濾波校準法[12]等。

隨著人工神經網絡（ANN）等深度學習技術的發展，很多專家學者利用深度學習的算法開展雷達定量降水估測的研究，由于深度學習是一種非參數運算，且具有自適應學習方法，能有效減少誤差，提高了降水估測的準確性[13-16]。Liu等[13]研究發現基于人工神經網絡的定量降水估測比用Z-I關系更精確。Chiang等[14]利用臺灣六次臺風降雨過程研究表明，采用動態人工神經網絡比Z-I關系有更高的精度和穩定性。邵月紅等[15]通過BP神經網絡方法估測定量降雨，精度要優于Z-I關系式，估測的降雨量與站點實測雨量吻合性較好。傅德勝等[16]基于徑向基函數神經網絡（RBF）建立的雷達定量估測降水模型，其對降水估測精度和穩定性高于以變分法得到的Z-I關系。

目前針對華南地區開展基于深度學習的定量降水估測算法的研究還比較少，本文主要應用三維網格雷達回波資料和自動氣象站數據開展基于梯度提升決策樹（GBDT）算法在雷達定量估測降水中的研究和應用，提高雷達定量降水估測的精度。

1 數據和預處理

1.1 數據

本研究使用的雷達資料為廣東省境內廣州、深圳、韶關、清遠、陽江、河源、汕尾、梅州、湛江和肇慶等10部S波段多普勒雷達在VCP21探測模式下，每隔6 min對9個仰角進行三維體掃描，為確保所有雷達能同步掃描，采用了北京敏視達公司的雷達同步控制和數據采集系統時間同步軟件進行時間同步。為計算方便，在雷達三維體掃描原始數據基礎上，利用Cressman插值法將極坐標格式的雷達探測資料插值到三維直角坐標系中，并采用了蘭紅平等[17]提出的“膨脹-侵蝕”算法對回波缺陷進行質量控制。由于在雷達組網拼圖時，若低于1.5 km，任意兩部雷達之間幾乎沒有共同的覆蓋區域，觀測的范圍太小；若高于5 km，雷達覆蓋的范圍擴大了，然而卻難以形成降水[18]。因此，本研究中用到經過雷達回波質量控制和拼圖處理后得到1.5、2.5、3.5和4.5 km共四層高度上CAPPI（constant altitude plan position indicator）數據，空間分辨率為0.01°×0.01°，水平網格點為約為500 km×300 km，數據覆蓋為2015—2017年的每年4—9月的降水過程。

另外，本研究中還用到該范圍內區域自動氣象站逐小時雨量觀測資料，自動氣象站分布如圖1。

圖1 自動氣象站分布圖（黑點） Fig. 1 Distribution of automatic weather stations(AWSs) (black dots)

1.2 預處理

由于雷達回波為空中水滴、冰晶等對雷達電磁波的反射或散射，自動氣象站記錄的是落到地面的水滴，在選取雷達的回波強度數據時，需要綜合考慮多種因素的影響。梁建茵等[19]認為，對低緯度地區，自動氣象站點上降水強度與其上空約5 km×5 km區域內的平均雷達回波強度關系較好，使用平均雷達回波強度的方法可以在一定程度上減少降水估測的誤差。因此，為了更有效的利用三維雷達回波強度數據進行定量降水估測，本文建立以自動氣象站為中心、1 h內逐6 min的5 km×5 km格點范圍的雷達回波強度數據（圖2）。研究中用到的自動氣象站降水資料為整點的1 h雨量，與雷達資料的匹配方式為：T時刻記錄的過去1 h雨量（T-1到T的1 h累計雨量）與T-1到T的逐6 min、10個時次5 km×5 km格點的雷達回波反射率相匹配，形成一條數據樣本，作為降雨定量估測模型的輸入數據。本文共建立100萬組訓練數據集和10萬組測試數據集。

2 基于GBDT的雷達定量估測降水算法

梯度提升決策樹是一個應用廣泛、效果較好的監督式機器學習模型，它于2001年由Friedman[20]提出的組合決策樹算法，由決策樹（decision tree）[21]和梯度提升（gradient boosting）[22]組合而成，提高回歸預測的精度。

GBDT是一種迭代的決策樹算法，其迭代的目的是通過計算上一次模型的負梯度來改進模型，每一次新模型的建立都減少了上一個模型的殘差，然后在殘差減少的梯度方向上建立新的決策樹；每個新模型的輸入都是上一個模型輸出，輸出是下一個模型的輸入，預測模型的最終結果為之前所有模型結果的相加。該算法由多個弱預測模型（通常是決策樹）通過集合的形式組合產生預測模型，所有樹的結論加起來作為最終答案。該方法的優點在于可以盡量避免過擬合現象，每迭代一棵樹，殘差代表的都是全局最優的絕對方向。由于GBDT算法考慮方差和偏差之間的權衡，具有較好的泛化能力，幾乎可以解決所有的數據回歸問題并具有高檢測精度的特性，且在非線性數據處理方面具有良好的效果，被廣泛應用消息檢測[23]、徑流預測[24]等研究領域。

本文將每一個數據樣本是具有一組特征屬性的雷達回波三維網格數據和自動氣象站小時雨量數據組成，定義訓練數據集：

圖2 以自動氣象站為中心的雷達回波反射率數據 Fig. 2 Three dimensional data of radar reflectivity for GBDT as input (The gauge location is at the center of the grid)

假設通過遞歸迭代分割將特征空間（即輸入變量）劃分為m個子集，則GBDT模型表示為m個模型函數相加來預測輸出：

基于GBDT的定量降水估測QPE模型算法的詳細步驟如下：

1）使用常量值初始化模型，估計使損失函數極小化的常數值，它是只有一個根節點的樹，即β是一個常數值：

迭代訓練模型，在上一次模型損失函數的梯度下降方向上建立模型，從m=1到M（M為迭代次數）

計算損失函數的負梯度在當前模型的值，將它作為殘差ri的估計值：

2）將式（5）中得到的估計殘差ri作為輸入，擬合到一顆回歸樹，求得回歸樹的葉節點區

3 個例分析

利用建立的訓練數據集進行GBDT算法的建模和訓練，得到基于GBDT算法的雷達定量降水估測模型，并用測試數據集進行檢驗和誤差分析。作為比較，本文還利用測試數據集建立了固定Z-I關系法和動態Z-I關系法進行定量降水估測。其中固定Z-I關系法應用汪瑛等[8]通過統計歷史上強降水日的雷達反射率因子和1 h降水資料，得到適合廣東地區的系數A和b的值，即Z=220.8I2.23，這個統計結果可直接用于雷達定量降水估測，且計算速度快。動態Z-I關系法主要利用自動氣象站雨量計降水資料，基于Z-I關系算法，實時計算當前時刻最優系數A和b的值，使雷達定量估測降水與區域自動站觀測降水的誤差最小，作為下一時次的Z-I關系計算多普勒雷達估測降水，該方法的優點是前期不需要收集大量的資料樣本做統計準備，只需要統計上一個小時的最優Z-I關系。

圖3給出了2017年4月21日07時的影響廣東區域一條颮線雷達定量估測降水和06—07時的1 h雨量實況。該小時降水的主要特點是有一條東北—西南向覆蓋廣東中部地區、小時雨量大于10 mm的降雨帶，在廣州和茂名存在兩個強降水中心（小時雨量在40 mm以上），其中最大小時雨量出現在廣州，為65.3 mm（圖3d）。圖3a為固定Z-I關系法的定量降水估測，與自動氣象站實況雨量相比，降水分布和落區較為一致，但是降水量級估測偏小很多，對兩個強降水中心的降水估測較實況比偏小20～30 mm。動態Z-I關系法對降水的估測比固定Z-I關系法有所改善（圖3b），對位于廣州的強降水中心估測在40 mm以上，雖然仍比實況偏小，但在降水量級上比固定Z-I關系法已有10 mm的提高。圖3c是基于GBDT算法反演的雷達定量估測降水，不同等級降水量的面積和強度相對固定Z-I關系法和動態Z-I關系法都有所改進，尤其對位于廣州的強降水中心的降水量在50～55 mm，與實況更接近，改善了降水量偏小的問題。

4 雷達定量降水結果和誤差分析

為了描述定量估測降水的誤差和精度特征，在檢驗和分析時，使用平均誤差（ME）、相對誤差（RE）和均方根誤差（RMSE）等進行統計評價，具體計算公式如下：

圖3 2017年4月21日07時雷達1 h定量降水結果和自動氣象站實測降水（a）固定Z-I關系法；（b）動態Z-I關系法；（c）GBDT算法；（d）雨量實況 Fig. 3 Radar-derived 1h QPE and gauge spatial distributions on 06:00 BT 21 April 2017(a) 1 h QPE of Z-I relations; (b) 1 h QPE of dynamic Z-I relations; (c)1 h QPE of GBDT scheme; (d) observation

式中，n為總樣本數，i為QPE雨量計的序列，Pi和Gi分別為雷達估測雨量和自動氣象站實測雨量。

用平均誤差評估數據的整體差異，其值越接近0則數據整體差異越小等評價方法；使用相對誤差評估數據的可信度，其值越小數據可信度越高；用均方根差評估雷達估測雨量與自動氣象站觀測雨量的離散程度，其值越小分布越集中；從整體上評估不同方法對降水估測精度的高低，所有誤差都在站點上計算。為了使評估的結果具有客觀性，剔掉沒有通過質量控制的自動氣象站觀測雨量數據，將1 h雨量劃分為1～5、5～10、10～20、20～30、30～40、40～50 mm和50 mm以上的7個降水量級，分析這三種方法對降水的反演效果。

圖4 給出了這三種方法的雷達定量降水和自動氣象站觀測雨量散點圖。越接近直線點表明雷達定量降水的效果越好，高于直線說明高估，低于直線點為低估。可以看出，固定Z-I關系法對雷達定量降水的反演存在低估，隨著降水量級的增大，低估程度趨于明顯，與實測降水的相關系數為0.71；動態Z-I關系法對雷達定量降水的估測準確性有所改進，雷達估測降水和實測降水的相關系數略有提高，為0.74；應用GBDT算法定量估測的降水與實況觀測相關性較好，相關系數提高到0.80，散點圖顯示集中性較好，改善了降水估測的精度，但是仍然存在低量級降水高估、高量級降水低估的現象。

表1是對不同量級降水的雷達定量降水誤差分析，可以看出，固定Z-I關系法對10 mm/h以下的降水存在高估，而對10 mm/h以上的降水明顯低估，對50 mm/h以上的強降水低估平均誤差達到34.14 mm，相對誤差為56.95%。動態Z-I關系法對雷達定量降水估測效果優于固定Z-I關系法，各量級的平均誤差和相對誤差較固定Z-I關系法有所減小，但仍然存在對10 mm/h以下量級降水高估，而10 mm/h以上量級降水低估的現象。這可能是由于動態Z-I關系法雖然根據自動站觀測雨量和雷達回波關系不斷調整A和b值，提高雷達定量降水估測能力，但仍然存在對同一時刻的降水類型沒有進行分類，所有降水類型都用同一個A、b值的Z-I關系的不足。

圖4 實測雨量和雷達定量估測降水散點圖：（a）固定Z-I關系法；（b）動態Z-I關系法；（c）GBDT算法 Fig. 4 The scatter diagram of hourly rainfall estimation: (a)Z-I relations; (b) dynamic Z-I relations; (c) GBDT scheme

基于GBDT算法反演的降水平均誤差和相對誤差都小于固定Z-I關系法和動態Z-I關系法的結果，尤其對30 mm/h以上量級降水的估測，相對誤差比動態Z-I關系法都提高了5%以上，對降水反演的系統誤差更小，改善了降水估測的精度。RMSE的含義是雷達估測雨量與自動氣象觀測的降雨量之間的離散程度，對5～20 mm/h量級的降水估測，動態Z-I關系法RMSE最小，但對于20 mm/h以上的強降水，GBDT算法估測降雨量的RMSE都小于固定Z-I關系法和動態Z-I關系法，表明GBDT算法對強降水的雨量估測擬合度較高。

表1 不同量級降水的雷達定量降水誤差分析 Table 1 Deviation statistics of different rainfall intensity categories

5 結論

利用廣東省新一代多普勒氣象雷達三維拼圖數據和地面自動氣象站雨量資料，建立了一種基于梯度提升決策樹GBDT算法的雷達定量估測降雨模型，通過與固定Z-I關系法和動態Z-I關系法進行誤差分析和對比，結果表明：

1）固定Z-I關系法對小于10 mm/h的降水量級存在高估，但對10 mm以上量級的降水存在低估，隨著降水量級的增大，低估現象越明顯；動態Z-I關系法優于固定Z-I關系法，對雷達定量降水的估測準確性有所改進。

2）通過對平均誤差、均方根誤差和相對誤差等的統計結果表明，基于梯度提升決策樹GBDT算法估測的降雨量精度要高于固定Z-I關系法和動態Z-I關系法，特別是對30 mm/h以上的強降水的估測，GBDT算法的強降水中心和分布與自動氣象站觀測較為吻合，對短時強降水低估現象有一定的改進作用。