孫 誠, 沈 芳, 唐儒罡
(華東師范大學 河口海岸學國家重點實驗室, 上海 200241)
我國長期以來飽受洪澇災害的頻繁侵擾, 嚴重威脅了人民群眾的生命和財產安全. 遙感技術不受空間限制, 且可迅速地獲得洪水淹沒信息, 已成為洪澇災害監(jiān)測和評估的常用手段[1-2]. 星載SAR(synthetic aperture radar, 合成孔徑雷達) 工作于微波波段, 具有地表穿透能力, 不受云霧限制, 突破了光學遙感技術受天氣影響大的局限性[3], 在洪水監(jiān)測方面具有獨特的優(yōu)勢.
閾值分割法及機器學習方法常被用于SAR 影像水體信息的提取. 閾值分割法的關鍵在于在圖像灰度的分布中確定合適的閾值, 以區(qū)分水體和非水體. 賈詩超等[4]研究了Sentinel-1SAR 雙極化數(shù)據(jù)之間水體信息提取的關系, 提出了基于閾值分割的SDWI (sentinel-1 dual-polarized water index) 水體指數(shù)法. 閾值分割法的操作邏輯簡單, 計算時間短, 但在地物類型復雜的情況下, SAR 影像上灰度值相等的像元并不能完全對應同一地物類型, 難以選擇最優(yōu)閾值, 使得閾值方法易受到圖像噪聲和強度不均勻性的影響, 因此基于影像單個像元的方法具有局限性[5].
為了更好地獲得SAR 影像中的信息, 有學者采用GLCM (gray-level co-occurrence matrix) 提取SAR 影像中的紋理信息, 建立多特征空間, 并利用機器學習模型進行水體信息的提取[6-7]. Lyu 等[8]提出了結合4 種灰度共生矩陣紋理特征與SVM (support vector machines) 分類器, 在區(qū)分水體目標區(qū)域與其他地物及對灰度共生矩陣提取水體等方面進行了初步探索. 胡德勇等[9]使用Radarsat-1 單波段單極化數(shù)據(jù)對水體和居民地信息進行了提取, 但由于單極化影像信息有限和分類算法的限制, 存在一定的錯分、誤提. 此外, 與單極化SAR 影像相比, 全極化SAR 影像的信息更為豐富, 圖像分類的性能更高[10]. 也有學者采用RF (random forest)[11]、GBDT (gradient boosting decision tree)[12]等機器學習算法, 但傳統(tǒng)的機器學習算法由于內存占用較大, 存在模型訓練時間較長的問題[13]. 近年來, 基于集成學習的LightGBM (light gradient boosting machine, 輕量級梯度提升器) 算法, 因其學習能力強及預測精度高的特點, 越來越多地被應用于各類學科領域[14-16], 相比于SVM、RF 及GBDT 等機器學習算法, 其具有更快的訓練速度和更低的計算代價, 可以快速地處理大數(shù)據(jù)量和高特征維度的數(shù)據(jù)[17]. 因此, LightGBM 算法更適用于對速度與精度有較高要求的洪澇災害淹水信息應急提取.
本文基于Sentinel-1SAR 影像, 旨在建立一種結合SAR 紋理信息和LightGBM 算法的洪水淹沒地區(qū)遙感應急監(jiān)測方法, 并將該方法的水體提取結果與SDWI 水體指數(shù)法、SVM、RF 及GBDT 等其他方法進行定量和定性對比, 以測試該方法的提取精度和運行效率, 最后對淮河流域中的典型洪水淹沒區(qū)域進行淹水范圍提取和分析應用實踐.
淮河流域是我國南北區(qū)域的自然分界線, 降水時空分布不均, 每年的6—9 月為汛期, 10 月至次年5 月為非汛期, 年平均降水量為920 mm, 汛期年降水量最高達1 600 mm, 流域平均每5 年就會發(fā)生一場較為嚴重的洪澇災害. 因此, 本文的研究區(qū)域選取了淮河流域中游主干道附近地區(qū), 該區(qū)域內支流眾多, 源短流急, 在汛期易形成洪澇災害. 研究區(qū)域范圍見圖1, 底圖為Sentinel-1 影像.
圖1 研究區(qū)地理位置Fig. 1 Map of the study area
采用SAR 影像為IW (interferometric wide swath) 模式下Level-1 級別Sentinel-1 的GRD (ground range detected) 影像, 包括VH (垂直發(fā)射水平接收)、VV (垂直發(fā)射垂直接收) 兩種極化模式, 分辨率為10 m, 重訪周期最短為6 d. 利用歐空局開發(fā)的SNAP 軟件, 對該模式下的影像進行了軌道校正、熱噪聲去除、輻射定標、相干斑濾波、分貝化處理及地理編碼等預處理. 軌道校正對從歐空局網(wǎng)站下載的Sentinel-1 原始影像中的原始軌道數(shù)據(jù)進行校正; 熱噪聲是SAR 系統(tǒng)自帶的背景噪聲, 會影響雷達得到后向散射信號的精度, 需要在軟件中進行熱噪聲去除; 輻射定標是將SAR 接收的后向散射信號轉化為后向散射系數(shù); 相干斑是SAR 影像分類時的噪聲, 在軟件中過濾相干斑可提高影像分類的精度; 由于SAR 接受的后向散射信號之間差距不明顯, 使用分貝化處理的方式可指數(shù)化后向散射系數(shù),便于可視化和影像分類; 最后使用地理編碼的方法對SAR 影像進行地理坐標的校正, 完成預處理[18].
研究區(qū)域內選取成像時間為2020 年3 月5 日、2020 年3 月17 日、2020 年5 月16 日及2020 年11 月12 日共4 天的Sentinel-1SAR 影像作為數(shù)據(jù)源. 使用隨機抽樣的方法在研究區(qū)域內裁剪樣本區(qū),并利用成像時間相同的總體云量小于20%且水體周圍無云干擾的Sentinel-2 影像, 對樣本區(qū)進行目視解譯, 獲得樣本區(qū)的水體范圍信息. 最終共得到訓練集552 萬個, 驗證集125 萬個.
本文提出了一種基于SAR 紋理信息和LightGBM 算法的WEGL (water extraction based on GLCM and LightGBM) 方法進行水體提取. WEGL 方法由兩個主要步驟構成: 第一是SAR 影像紋理信息的提取 (詳見2.1 節(jié)); 第二是LightGBM 算法構建, 完成LightGBM 參數(shù)優(yōu)化和紋理特征變量優(yōu)化 (詳見2.2 節(jié)). 首先對研究區(qū)域的Sentinel-1 影像進行預處理和裁剪, 利用Sentinel-2 光學影像與Sentinel-1 SAR 影像進行配準, 得到樣本區(qū)的水體分布, 進一步基于GLCM 計算SAR 紋理信息, 構建訓練樣本數(shù)據(jù), 進行LightGBM 算法訓練, 得到水體提取結果. 具體步驟實現(xiàn)的流程見圖2.
圖2 方法實現(xiàn)流程Fig. 2 Flow chart of the method performance
GLCM 最早由Haralick 等[19]提出, 是一種描述在預定計算窗口內相鄰像元或間隔一定距離像元灰度相關關系的矩陣, 常使用于地物復雜的遙感影像中[20-21]. 本文試驗使用的8 個紋理特征及其描述見表1. 表1 各式中:N是GLCM 窗口的行列數(shù);i,j為圖像灰度級數(shù);u代表均值;σ代表方差;Pi,j代表從圖像灰度為i的像素出發(fā), 灰度為j的像素出現(xiàn)的概率.
表1 紋理特征表達式及描述Tab. 1 Expressions of texture feature and their description
影像的紋理信息與計算灰度共生矩陣選擇的方向、步長、影像的灰度級及窗口的大小有關[22-23].經(jīng)過大量試驗, 選擇0°、45°、90° 及135° 這4 個方向的GLCM 平均值, 步長大小選擇1, 影像灰度級選擇64, 窗口大小選擇7 × 7, 對表1 中的8 個紋理特征信息進行提取.
LightGBM 是一種基于GBDT 框架的分類模型, 其帶有深度限制的Leaf-wise 葉子生長策略及直方圖作差加速等創(chuàng)新技術[24], 可支持高效的并行訓練, 能夠快速處理大數(shù)據(jù)量和高特征維度的數(shù)據(jù)[25].
2.2.1 精度評價指標
選取了總體精度、交并比、F1指標及Kappa 系數(shù)共4 種精度評價指標, 對試驗結果精度進行定量分析.
總體精度 (AO) 常用于衡量圖像分類的性能, 表示預測正確的水體像元和陸地像元占總像元數(shù)量的比重, 其表達式為
式(1)中:nTP表示水體像元被判定為水體像元的數(shù)量;nFP表示陸地像元被判定為水體像元的數(shù)量;nTN表示陸地像元被判定為陸地像元的數(shù)量;nFN表示水體像元被判定為陸地像元的數(shù)量.
交并比 (intersection of union, IoU) 表示預測結果與真實地物的相關度, 其表達式為
F1指標是衡量精準率和召回率的綜合指標, 其表達式為
式(3)中: 精準率 (P) 表示真實水體像元占所有預測為水體像元的比例, 表達式為
召回率 (R) 表示真實水體像元占所有實際為水體像元的比例, 表達式為
Kappa 系數(shù) (K) 為衡量分類精度的經(jīng)典指標, 其表達式為
式(6)中:p0代表水體和陸地像元中正確分類的像元數(shù)量,pe的表達式為
式(7)中:a1和a2分別代表水體和陸地的實際像元數(shù);b1和b2分別代表水體和陸地的預測像元數(shù).
2.2.2 參數(shù)優(yōu)化
在進行訓練之前, 需要為LightGBM 算法尋找最優(yōu)參數(shù), 確定最優(yōu)識別精度下的最優(yōu)參數(shù)組合.在其他參數(shù)都保持默認參數(shù)的情況下, 使用遍歷的方法尋找LightGBM 算法收斂的迭代次數(shù). 經(jīng)過試驗, 200 次迭代之后, LightGBM 算法的水體識別準確率趨于穩(wěn)定, 不再有明顯的提升. 保持LightGBM算法迭代次數(shù)為200 次, 進一步測試算法識別準確率隨樹的最大深度 (max_depth) 變化的情況, 試驗結果如圖3 所示. 當樹的最大深度較小時, LightGBM 算法的準確率曲線隨迭代次數(shù)收斂的速度變慢,且識別準確率較低. 隨著樹的最大深度的提高, LightGBM 算法的準確率曲線隨迭代次數(shù)收斂的速度加快, 且識別準確率明顯升高. 從圖3 中可看出, 樹的最大高度為8 時, 模型的識別準確率表現(xiàn)最好.在確定模型最佳迭代次數(shù)和樹的最大深度后, 使用網(wǎng)格滾動調參的方式, 構建不同參數(shù)的排列組合,代入模型中進行測試, 以得到最優(yōu)識別精度下的最優(yōu)參數(shù)組合. LightGBM 算法的相關參數(shù)和調參的范圍及步長如表2 所示.
表2 LightGBM 算法參數(shù)及調優(yōu)Tab. 2 Parameters and optimization of LightGBM algorithm
圖3 LightGBM 算法識別準確率Fig. 3 Accuracy of LightGBM algorithm for identification
2.2.3 紋理特征變量優(yōu)化
通過采用LightGBM 算法的信息增益函數(shù), 表征紋理特征變量的重要性排序, 以代表特征變量對模型訓練的貢獻程度[14]. 圖4 分別顯示了VH、VV 兩組單極化影像中, 各紋理特征變量對模型訓練的重要性排序及對識別水體的貢獻度. VH 影像中, 特征變量重要性排序前4 位的特征分別為: 均值、均質性、相異度、相關性; VV 影像中, 特征變量重要性排序前4 位的特征分別為: 對比度、均值、相關性、角二階矩. 進一步, 將上述兩組特征重要性排序前4 位的特征組合成8 個紋理特征并進行模型訓練, 得到的LightGBM 算法水體提取精度為98.40%, 與使用16 個紋理特征得到的98.41%精度基本一致. 因此, WEGL 方法剔除了一些特征貢獻度較低的紋理特征, 在保持精度基本一致的情況下, 可縮減一半的數(shù)據(jù)量和內存占用, 大大提高了水體提取效率.
圖4 VH、VV 極化影像特征重要性Fig. 4 Importance degree of VH/VV SAR image features
Sentinel-1 系列衛(wèi)星的特點是有多種成像方式, 可實現(xiàn)單極化、雙極化等不同的極化方式[26]. 為驗證基于雙極化SAR 影像的分類性能是否優(yōu)于單極化影像, 將單極化影像和雙極化影像加入對比試驗中. 為驗證LightGBM 算法的分類性能, 在保持訓練集相同的情況下, 與LightGBM 同時進行了SVM、RF 及GBDT 算法的預測試驗, 并記錄各方法在不同數(shù)據(jù)特征條件下的預測精度及訓練耗時,具體參見表3.
表3 不同算法的對比試驗Tab. 3 Comparative tests of different algorithms
從表3 的試驗結果看出, 對于每種算法而言, 基于VH + VV 雙極化方式的SAR 影像的算法提取精度要優(yōu)于VH、VV 單極化方式的SAR 影像, 說明雙極化SAR 影像中的紋理信息更為豐富, 水體識別的性能更高, 故WEGL 方法選擇雙極化SAR 作為輸入的影像. 在水體提取精度方面, 基于雙極化影像的SVM 的水體提取準確率為95.19%, RF 和GBDT 的準確率均高于SVM, 分別為96.23%、96.41%. LightGBM 算法的水體提取準確率最高, 為98.40%. 在訓練耗時方面, SVM 的用時最長, 共用時超過4 h, 說明在特征維度較大的情況下, SVM 算法擬合時間較長. RF 和GBDT 算法的訓練耗時有進一步優(yōu)化, 為1 h 左右, 基于決策樹思想的算法的訓練效率明顯提升. LightGBM 算法的訓練耗時最短, 僅為3 min, 說明在處理海量數(shù)據(jù)上LightGBM 算法具有明顯優(yōu)勢, 更適用于洪水淹沒地區(qū)應急提取等場景.
為驗證WEGL 方法的合理性及有效性, 在研究區(qū)域內選取了河道區(qū)、湖泊區(qū)及洪水淹沒區(qū)3 類典型區(qū)域進行方法的評價, 并將其與SDWI 水體指數(shù)法、SVM、RF 及GBDT 算法進行對比.
圖5 區(qū)域為河道區(qū), 經(jīng)目視解譯和百度地圖驗證, 1 處為含水量較高的淺灘, 2 處為河道邊的道路.SDWI 水體指數(shù)法對河道提取較為完整, 但存在錯分現(xiàn)象, 其將1 處的淺灘和2 處的道路識別為水體;對比前者, SVM 算法在1、2 處仍存在錯分現(xiàn)象, RF、GBDT 算法在2 處無錯分現(xiàn)象, 但將1 處的淺灘識別為水體; WEGL 方法識別出了完整的河道部分, 在1 處和2 處的錯分現(xiàn)象明顯改善, 表現(xiàn)最好.
圖5 河道區(qū)水體提取結果Fig. 5 Water body extraction results of river
圖6 區(qū)域為湖泊區(qū), 3、4、5 處為湖泊水體的邊緣. 水體邊緣處水深較淺, 與水體內部的像元灰度值存在差異. 基于像元的SDWI 水體指數(shù)法無法準確識別邊緣水體, 存在較多錯分和漏提取現(xiàn)象;SVM、RF、GBDT 算法去除了大量噪聲影響產生的錯分現(xiàn)象, 但在3、4、5 處水體邊緣的漏提取現(xiàn)象仍然存在; WEGL 方法對水體邊緣的提取更完整和平滑, 錯分和漏提取的現(xiàn)象較少.
圖6 湖泊區(qū)水體提取結果Fig. 6 Water body extraction results of lake
圖7 區(qū)域為洪水淹沒區(qū), 6、7 處為泥沙含量較高的洪水水體. SDWI 水體指數(shù)法提取的水體在6、7 處內部存在噪聲現(xiàn)象; SVM 算法在6、7 處改善了水體提取的噪聲, 但漏提取現(xiàn)象仍存在; RF、GBDT 算法的錯分現(xiàn)象較少, 但仍存在水體內部噪聲; WEGL 方法對洪水水體的提取結果表現(xiàn)最好,水體內部噪聲和漏提取現(xiàn)象較少.
圖7 洪水淹沒區(qū)水體提取結果Fig. 7 Water extraction results of flooded area
各方法提取水體精度的定量評價結果如表4 所示. WEGL 方法的相關評價指標最好, 與SDWI 水體指數(shù)法比較, 總體精度、IoU、F1指標及Kappa 系數(shù)分別提高了2.00%、16.89%、0.11 及0.12; 與SVM 比較, 總體精度、IoU、F1指標及Kappa 系數(shù)分別提高了0.60%、5.35%、0.04 及0.04; 與RF 比較, 總體精度、IoU、F1指標及Kappa 系數(shù)分別提高了0.34%、6.10%、0.04 及0.04; 與GBDT 比較, 總體精度、IoU、F1指標及Kappa 系數(shù)分別提高了0.24%、3.30%、0.02 及0.02.
表4 水體提取結果與精度對比Tab. 4 Accuracy comparison of water extraction results
使用WEGL 方法對2020 年7 月淮河流域洪水受災較嚴重的霍邱和潁上段進行洪水淹沒監(jiān)測. 圖8分別展示了霍邱和潁上段災前災中的水體范圍及淹沒范圍的空間分布. 監(jiān)測表明, 2020 年7 月27 日淮河流域霍邱和潁上段北岸受災情況嚴重, 相比于災前(2020 年7 月3 日)的監(jiān)測影像, 水面積由194.92 km2變?yōu)?43.50 km2, 監(jiān)測區(qū)域內的淹沒面積達到148.58 km2. 根據(jù)霍邱縣人民政府網(wǎng)站(http://www.huoqiu.gov.cn) 和潁上縣人民政府網(wǎng)站 (http://www.ahys.gov.cn) 的報道, 2020 年7 月兩縣受強降水影響嚴重, 霍邱縣農作物受災面積為542.27 km2, 潁上縣湖區(qū)受災面積共17.39 km2, 其中農作物為15.83 km2. 上述受災統(tǒng)計結果較為碎片化, 無法精確到具體的受災區(qū)域, WEGL 方法的洪水淹沒監(jiān)測結果更為客觀和完整, 有助于對受災情況的整體評估.
圖8 霍邱和潁上段災前災中水體信息提取情況Fig. 8 Outline drawing of water extraction results in pre-flood and flooding of Huoqiu and Yingshang County
針對目前常用的SAR 影像水體提取方法精度不足、運行效率低等問題, 通過GLCM 提取SAR 影像紋理信息, 并進行LightGBM 算法參數(shù)優(yōu)化和紋理特征變量的優(yōu)化, 建立了WEGL 方法, 大幅提高了水體提取的運行效率, 并在水體提取精度上有進一步優(yōu)化. 經(jīng)試驗分析, 與SDWI 水體指數(shù)法、SVM、RF 及GBDT 算法相比, WEGL 方法提取河道、湖泊和洪水淹沒區(qū)的精度均具有優(yōu)勢, 在一定程度上抑制了道路、裸地等地物的影響, 提取的水體邊緣更加清晰且完整. 除了目標提取精度的優(yōu)勢,WEGL 方法的運行效率也顯著提升, 更加適用于洪澇災害淹沒地區(qū)的應急監(jiān)測. 將WEGL 方法應用于淮河流域霍邱和潁上段的洪澇災害監(jiān)測, 發(fā)現(xiàn)2020 年7 月洪水期間, 水面積由194.92 km2變?yōu)?43.50 km2, 受災面積達到148.58 km2, 結果表明WEGL 方法具有時空可擴展性, 可用于不同時期和區(qū)域的洪澇災害監(jiān)測.
WEGL 方法成功實現(xiàn)了洪澇災害期間淹水范圍的快速監(jiān)測, 為后續(xù)研究中需要大面積水域信息快速提取的場景提供了新的思路. 同時, 可進一步開發(fā)集成WEGL 方法的軟件系統(tǒng), 實現(xiàn)淹水信息和洪澇受災情況的自動獲取. 因受限于Sentinel-1 衛(wèi)星的重訪周期, 對重點地區(qū)進行全天候監(jiān)測的難度較大. 在后續(xù)研究中可嘗試基于WEGL 方法, 加入不同類型的SAR 數(shù)據(jù), 提高衛(wèi)星的觀測頻率, 進一步提高洪澇災害淹沒范圍監(jiān)測的時效性.