基于分形理論與SVM的河冰高分遙感影像智能識別方法研究

苑希民，韓 超，徐浩田，田福昌

(天津大學水利工程仿真與安全國家重點實驗室，天津 300350)

凌汛災害是黃河影響最為嚴重的自然災害之一，其影響因素眾多，成災機理復雜，并具有突發性和持續性，為黃河凌汛防御工作帶來諸多挑戰，素有“伏汛好搶，凌汛難防”之說[1]。凌汛災害一旦發生，將對沿岸農業生產、基礎設施以及人民生命財產帶來嚴重影響[2]。近年氣候變化影響下，黃河凌汛災害發生頻率呈現升高趨勢，尤其是內蒙古河段，特殊的地理位置、氣候條件及河道走向使其成為黃河凌汛災害最為嚴重的河段之一[3]，1987～2018年間發生不同程度的凌汛災害百余次，由平均1.6年/次上升到0.3年/次[4]，這對我國凌汛災害防御工作提出了更高要求。

眾多學者針對凌汛災害開展了大量研究，苑希民等[5]采用基于遺傳算法的神經網絡方法建立了凌情智能耦合預報模型，對流凌、封河、開河日期進行了預報；劉吉峰等[6]分析了黃河寧蒙河段冰凌2000年以來的變化特點；顧潤源等[7]研究了氣候變化對黃河內蒙古段凌汛期的影響。目前在凌汛研究中，數據分析和數值模擬等主要手段對觀測數據的需求量較大，傳統野外觀測得到的數據在尺度等方面具有局限性，已在水利領域得到廣泛應用的遙感技術具有宏觀性能高、更新周期短、抗人為干擾因素強等優點[8-9]，可有效地彌補傳統觀測方法的不足[10]。

當前遙感監測冰情的技術，大多應用在海洋和湖泊等較為寬闊的水域，主要用于監測海冰面積和邊緣[11]、冰湖數量[12]、冰湖面積[13]等，相比而言，由于大多數河流有限的寬度，河冰的研究較為局限，大多是針對低分辨率大范圍影像，Naira等[14]提出了一種基于閾值的決策樹圖像分類算法來處理MODIS數據并確定河冰范圍的方法；H B Wang等[15]提出了一種基于無人機遙測技術的黃河冰柱定位監測的視頻數據處理方法；Kraatz等[16]提出了一個可供選擇的MODIS河冰監測算法，它可以在無云條件下和通過一些半透明的云識別河冰；Beaton等[17]開發了一種利用MODIS圖像自動檢測駝鹿河等五條河流開河的方法。

以往研究均為對大幅影像的解譯[18]，缺少對河冰細部特征的識別分析以及河冰遙感影像的精細化分類，這也正是本文的先進性所在。

本文針對凌汛災害頻發的黃河內蒙古段，選取典型河段為研究區域對象，對黃河內蒙古包頭段黃牛營子村處彎道的河冰進行分形智能分類識別研究。

1 數據及原理

1.1 技術路線

本文技術路線如圖1所示，具體為：均勻分塊裁剪研究區域的河冰高分遙感影像，在經過灰度化和濾波去噪等預處理措施之后，進行分形特征邊緣檢測，基于該檢測結果和對遙感影像的目視解譯選擇分類的分類樣本，通過支持向量機算法進行智能分類，基于混淆矩陣計算分類精度并對比分析分類結果與分形結果。

圖1 技術路線圖Fig.1 Technology roadmap

1.2 影像預處理

對獲取的遙感影像進行預處理，包括分塊裁剪，灰度化和濾波去噪。將遙感影像均勻分塊，保證每一圖塊僅包含一種河冰類型，便于后續處理區分。

為獲得清晰的、高質量的遙感圖像，本次研究中采用中值濾波進行降噪預處理。中值濾波是當前應用最廣泛的去噪方法之一，它在一定的條件下可以克服線性濾波、均值濾波等帶來的圖像細節模糊，而且對濾除脈沖干擾及圖像掃描噪聲最為有效[19]。

中值濾波器一般采用一個含有若干個點的滑動窗口，將窗口中幾個點灰度值的中值來替代指定點(窗口的中心點)的灰度值。大多選取二維窗口，窗口的尺寸逐漸增大，直到濾波效果滿意為止。

圖2為河冰原始影像圖塊及其灰度圖塊，對其進行模板為7×7，9×9和11×11的中值濾波后的灰度圖像及二值化后的圖像。從圖中可以看出，采用7×7的模板對原圖像進行中值濾波后，去除了部分噪聲，但圖像中噪聲還是很明顯；采用11×11模板的中值濾波雖然噪聲濾除噪聲能力加強，但是圖像出現模糊和斷續現象;采用9×9模板濾波后噪聲得到了很好的抑制，同時圖像特征得到了很好的保存，因此最終采用9×9中值濾波模板。

圖2 河冰影像及其灰度化和二值圖像Fig.2 River ice image and its grayscale and binary image

中值濾波有三方面優點：1)降低噪聲能力較強；2)在灰度值變化較小的情況下可以得到很好的平滑效果；3)不會使圖像的邊界部分過分模糊。

1.3 河冰影像的分形特征邊緣檢測

圖像中不規則的對象無法用傳統的歐幾里德幾何學來描述[20]。Benoit B. Mandelbrot于1975年創立了分形幾何學，用分形(Fractal)一詞來表述那些沒有特征長度，具有無限精細結構的圖形、構造及現象[21]。分形幾何圖形具有自相似性和遞歸性，易于計算機迭代，擅長描述自然界存在的景物[22-23]。

在實際應用中，常用Richardson定律來計算分形維數[24]。

M(ε)=Kεd-D.

(1)

式中ε=1，2，3……為尺度因子，M(ε)約是尺度ε下的度量特征值，D是分形維數，d是拓撲維數，K是分形系數。對于二維灰度圖像，M(ε)約取為圖像表面積測度A(ε)，則有

A(ε)=Kε2-D.

(2)

一幅圖像可以看成高度正比于其灰度值的山丘，這個曲面的上下ε構成厚度為2ε的“毯子”。對于不同的ε，“毯子”的面積即圖像的表面積A(ε)可以由“毯子”的體積除以2ε得到。

設f(i,j)代表圖像的灰度值，U(i,j,ε)，B(i,j,ε)分別表示上表面和下表面的灰度值，令

U(i,j,0)=B(i,j,0)=f(i,j).

(3)

上下兩張“毯子”分別以如下方法變化：

U(i,j,ε+1)=max {U(i,j,ε)+1,maxm,n∈η[U(m,n,ε)]}.

(4)

B(i,j,ε+1)=min{B(i,j,ε)-1,minm,n∈η[B(m,n,ε)]}.

(5)

式中η={(m,n)|d[(m,n),(i,j)]≤1}，d(?)表示兩點之間的距離。

于是“毯子”的體積為V(i,j,ε)=∑(i,j)∈R(U(i,j,ε)-B(i,j,ε)) ，R表示圖像上以(i,j)為中心的取值區域，故表面積A(i,j,ε)=V(i,j,ε)/2ε。

計算不同尺度下的A(i,j,ε)，由式logA(ε)=(2-D)logε+logK可知利用點對 [logA(i,j,ε),logε]，采用線性最小二乘擬合的方法，由擬合直線的斜率即可求出分形維數D。

黃河河冰的形狀是自然產生的，同一種河冰范圍內的形狀具有一定的自相似性，因此可以利用其分形特征，來實現對河冰影像的分割，從而更加清晰地區分不同種類的河冰。

(6)

K反映了圖像灰度表面積隨尺度變化的空間變化率。對式(2)兩端取對數，得

logA(ε)=(2-D)logε+logK.

(7)

上式表示在logA(ε)-logε坐標下的一條直線，logK為該直線在縱坐標軸logA(ε)上的截距，K相當于該尺度下的灰度曲面面積。當光滑的曲面或灰度變化緩慢的灰度曲面，K值較小；而當起伏較大的灰度曲面或灰度變化較為劇烈的曲面，K值較大。不同紋理灰度表面之間，灰度起伏變化相對無論哪一種紋理圖像來說都更加顯著，因此K值也能夠反映圖形灰度表面的粗糙程度。大多數紋理圖像可以用分形模型進行描述[25]，河冰影像即由不同紋理區域組成，在不同紋理灰度表面之間(即圖像的邊緣處)，灰度變化比較大，即K值較大。所以我們可以用K值作為分形特征，對河冰影像進行分割。算法如下：

(1)以M×M的窗口作為局部處理區域，從河冰影像的起始點開始，從左到右，從上到下，按照ε-毯子法依次計算每個窗口中心像素的分形特征K，從而將河冰影像的灰度空間映射為分形特征K空間；

(2)采用自適應門限對K空間進行二值化。找出K空間的最大值Kmax，并選取適當的系數a (0

1.4 河冰影像的智能分類及評價

支持向量機(SVM)是一種基于統計學習理論的模式識別算法，能夠根據有限的樣本信息尋找最優分類結果。SVM非線性分類的思想是把輸入樣本經過非線性變換映射到高維核空間，在高維核空間中尋找最優分類超平面，將樣本集分為兩類，并使兩類間隔最大[26]。

設有一個最優超平面[27]能將兩個類別分開，可用下式(8)表示：

f(x)=wTφ(x)+b.

(8)

其中w、b分別為超平面的法向量和截距，φ(x)表示將x映射到更高維空間后的特征向量。

為使間隔最大化，有

(9)

(10)

利用核函數K(xi,yi)進行求解，如下

K(xi,xj)=<Φ(xi)，Φ(xj)>=Φ(xi)TΦ(xj).

(11)

即可得到，

(12)

其中αi指樣本。

SVM中有3種非線性核函數:多項式核函數、徑向基核函數和Sigmoid核函數。本研究選擇徑向基核函數作為分類函數。與其他核函數相比，徑向基核函數的參數更少，且性能穩定[28]。

本研究選擇SVM基于以下幾點考慮：

(1)SVM是基于小樣本統計理論的學習方法，可以較為快速地處理較大數量的訓練樣本[29]；

(2)SVM利用核函數將本研究中低維線性不可分問題轉化為高維線性可分問題，同時利用核函數簡化了高維空間問題的求解難度，提高了分類效率；

(3)SVM是基于結構風險最優化原則，與其他傳統學習方式相比，避免了過學習問題，泛化推廣能力強[30]。

在分類型模型評判的指標中，常見的方法有如下3種，分別是混淆矩陣、ROC曲線、AUC面積。本研究選用混淆矩陣進行分類后的評價。混淆矩陣就是分別統計分類模型歸錯類，歸對類的觀測值個數，然后把結果放在如表1中展示出來。

表1 混淆矩陣精度報表Table 1 Confusion matrix accuracy report

驗證樣本(即真實值)是a，模型認為是a的數量為TP，驗證樣本是a，模型認為是b的數量為FN； 驗證樣本是b，模型認為是a的數量為FP，驗證樣本是b，模型認為是b的數量為TN。

將這4個指標一起呈現在表格中，就能得到混淆矩陣。

主要的評價指標為準確率ACC和Kappa系數，ACC代表分類模型所有判斷正確的結果占總觀測值的比重[31]，其計算如下式

(13)

Kappa系數是主要用于一致性檢驗的分類精度指標，起計算如下式

(14)

2 河冰高精度分類識別過程與結果分析

2.1 研究區域及實驗區域概況

黃河內蒙古段位于黃河干流最北端，緯度較高，相比黃河上游蘭州段北進約4.5°，從伊克昭盟的拉僧廟進入內蒙古境內，至榆樹灣出境，全長840 km，而總落差僅162.5 m[32]，研究區域圖如圖3所示。內蒙古的冬天嚴寒而漫長，黃河流經之處的冬季平均氣溫較上游驟降4～6℃[33]。由于氣溫的大幅降低以及河道比降小等因素，黃河內蒙古段下游先封河，上游先開河，流凌封凍通常早蘭州20余天，而解凍開河要晚一個多月。同時該河段彎道和河面束窄遍布，極易發生卡冰結壩[34]，尤其是在每年的開河期，大量冰水混合下泄，會形成極大的凌峰，使得黃河內蒙古段凌汛災害頻發且嚴重，據1951～2011年寧蒙河段各斷面凌情及凌災統計，平均每年發生冰情5.77次，每兩年就會發生一次較嚴重的凌災[35]，這對兩岸人民群眾的生命財產安全形成了巨大的威脅。

圖3 研究區域圖Fig.3 Study area map

此次實驗區域如圖4所示，位于黃河內蒙古包頭段黃牛營子村處彎道，所選區域尺寸為3 000 m×3 000 m。影像區域內河冰種類豐富，包括冰塞(下面帶有冰礁的冰蓋)，透明冰，清溝以及平滑冰蓋，實驗中按照這四類冰對影像進行分類。

圖4 實驗區域Fig.4 Experimental area

2.2 數據介紹

本研究采用的是Worldview-3遙感數據，成像時間為2018年2月24日，分辨率0.31m每像素，包括紅綠藍3個波段。研究區域影像有以下特點：

(1)清溝顏色呈黑色，透明冰呈深藍色，目視解譯較為清楚明了；

(2)清溝上下游有冰蓋及冰礁等河冰，均呈白色，冰礁在影像的表現形式上表現為圓狀斑點；

(3)由于較薄的平滑冰蓋與遍布冰礁的冰塞顏色相仿，不容易區分，但是后者較前者導致凌汛的風險大大增加，因此需要進一步的處理對二者進行細致的區分。

基于以上這些圖像特征，傳統的單一的算法很難實現準確的區分和識別，因此融合多種處理技術，采用有針對性的處理方法是很有必要的。

2.3 影像的預處理

將區域圖均勻分割成900個100 m×100 m的圖塊，剔除實驗區域中的裸地、森林等與河冰無關的圖塊，保留涉及到河冰的圖塊，如圖5所示。

圖5 均勻分割后圖塊Fig.5 Evenly segmented graph block

保證每一個圖塊僅包含一種河冰的種類，這樣在分形程序處理之后，通過得到的不同處理結果能夠有效地區分河冰的種類。

對圖塊進行灰度化，避免影像失真，易于后續處理。

濾波去噪。使用中值濾波方法進行去噪，經調試比較，采用9×9模板的中值濾波去噪能力最強，同時也能夠較好地保留河冰的特征。

2.4 河冰影像的分形特征邊緣檢測實驗

采用上述基于河冰分形特征的影像邊緣檢測與分割的方法進行實驗，選擇2×2窗口作為局部處理區域，系數a取0.8，選取尺度為2<ε<8。

通常認為，黃河河曲段所結冰蓋是由岸冰、流凌冰塊、水內冰花結合生成[36]。在水流有湍急流動的區域，水體產生大量的水內冰，水內冰向下游流動、積聚，與流凌冰塊或岸冰粘結附著，另外黃河泥沙含量高，在結冰的過程中，還會出現異類成核的現象，即水中圍繞上游下泄冰塊以及水體中的雜質首先結冰，在冰蓋形成時后也處于冰蓋中冰厚最大的位置，極易形成冰礁和錨冰等阻水，增加凌汛風險。

從實驗結果看，邊緣檢測能夠很好地識別到冰礁，如圖6(a)中的冰礁位置，分形處理后如圖6(b)所示，每一處冰礁被識別成團狀，1986年在國際冰凌學術會議上曾給冰塞定義為：阻礙水流的碎冰或冰花的固定集合體。從這個意義上來說，大的錨冰或阻水的冰礁都屬于冰塞，只是規模大小和部位不同。因此本研究將識別出冰礁的冰蓋統一歸類為冰塞。

圖6 冰塞影像及其分形處理結果Fig.6 Ice jam image and its fractal processing results

而對于厚度較薄的平滑冰蓋來說，通常是由兩岸岸冰向河心延伸而成，另外清溝的出現減少了向下游輸送的流冰，使得清溝下游的冰蓋厚度較薄，表面平滑，分形紋理不夠明顯，因此分形程序無法檢測到邊緣，其處理結果的像素灰度值全部為0，表現為黑色。

通過不同的分形結果，可以有效地區分冰塞和平滑冰蓋，全部圖塊識別結果如圖7所示。從分形結果圖中可以看出，首先對于清溝的識別十分準確，可以看出分形識別的清溝輪廓與原圖相符；另外冰塞主要分布在圖中紅框的位置，而平滑冰蓋主要分布于圖中黃框位置。

圖7 圖塊分形結果圖Fig.7 Graph block fractal result graph

分析清溝產生的原因，首先，上游彎道處形成冰塞卡冰，使下游的來冰量減少，無法形成冰蓋；其次，整個彎道呈幾乎180°拐彎，流凌和冰塊在離心力的作用下堆積在下游彎道的外側，使彎道內側形成清溝。該影像時間為二月末，黃河沿岸氣溫有所回升，此時清溝的出現和擴大也預示著開河的臨近，有關部門應該加強警惕，謹防開河期凌汛災害的發生。

2.5 河冰影像智能分類

對2017年11月27日(距所研究高分影像的時間最近)的遙感影像進行水體識別，獲得了較為精確的水體范圍，由于黃河內蒙古段支流流入補給的流量較少，冬季水位變化小[37]，因此可以認為該水體范圍與河冰范圍相差無幾，可以用水體范圍圖層對研究區域進行裁剪，剔除非河冰部分的像素，減少其他地物對河冰識別精度的干擾，裁剪之后影像如圖8所示。

圖8 裁剪后的實驗區域 圖9 分類樣本選取圖 圖10 支持向量機分類結果圖Fig.8 Trimmed experimental area Fig.9 Classification sample selection diagram Fig.10 SVM classification results

根據分形特征邊緣檢測結果得到的河冰類別，結合目視解譯，對原始影像進行分類樣本的選取。將影像中河段分為前兩公里和后一公里兩段，在前兩公里河段進行分類樣本的選取，后一公里河段用于驗證。共選取4種分類樣本，共246個，如圖9所示(其中紅色是冰塞63個，綠色是平滑冰蓋39個，藍色是透明冰59個，黃色是清溝85個)。

本研究實現SVM多分類的方法分兩步，第一步是依次把某個類別的樣本歸為一類，其他剩余的樣本歸為另一類，根據4個類別樣本構造出4個SVM分類器，如下：

(1)冰塞為一類，平滑冰蓋、透明冰、清溝為另一類；

(2)平滑冰蓋為一類，冰塞、透明冰、清溝為另一類；

(3)透明冰為一類，平滑冰蓋、冰塞、清溝為另一類；

(4)清溝為一類，冰塞、透明冰、平滑冰蓋為另一類。

第二步是將上述4種分類結果進行測試，每個測試都有一個結果f1(x)，f2(x)，f3(x)，f4(x)，分類時將未知樣本分類為具有最大分類函數值的那類，于是最終的結果便是這4個值中最大的一個作為分類結果。利用支持向量機進行分類，分類結果如圖10所示。

2.6 分類結果評價

基于混淆矩陣對遙感影像的分類精度進行評價，通過對研究區域的目視解譯結果隨機取樣，在后一公里河段隨機生成 125個樣本作為精度驗證樣本分類結果的混淆矩陣如表2所示。

表2 分類結果混淆矩陣Table 2 Confusion matrix of classification results %

從表1可得，與驗證樣本對比，該方法將94.21%的平滑冰蓋像元正確分類為平滑冰蓋，5.76%的平滑冰蓋像元分類為冰塞，0.03%的平滑冰蓋像元無法分類，分類為其他類別；將96.68%的清溝像元正確分類為清溝，3.32%的清溝像元無法分類，分類為其他類別；將97.71%的透明冰像元正確分類為透明冰，2.24%的透明冰像元分類為冰塞，0.05%的透明冰像元無法分類，分類為其他類別；將79.87%的冰塞像元正確分類為冰塞，18.69%的冰塞像元分類為平滑冰蓋，1.43%的冰塞像元分類為透明冰，0.01%的冰塞像元無法分類，分類為其他類別。即本方法準確預測分類了94.21%的平滑冰蓋像素，96.68%的清溝像素，97.71%的透明冰像素以及79.87%的冰塞像素，根據式(13)計算準確率達93.05%，即研究區域影像中共有93.05%的像素被準確分類。另外Kappa系數達 0.9026，說明本分類方法在考慮錯分漏分情況下依然表現出較好的分類效果，說明本分類方法較為準確可信。

3 結語

針對黃河河冰高分遙感影像的幾何特征、空間特征和紋理特征等復雜信息，本文選取黃河內蒙古段為研究對象，提出一種高分遙感影像凌汛災害信息識別提取方法，開展河冰影像精細化分類與智能識別研究，具體研究結論如下：

(1)對河冰的高分遙感影像進行分形特征邊緣檢測，精細化提取河冰信息，將河冰影像分為冰塞、透明冰、清溝及平滑冰蓋四類，準確地識別出河冰的種類及不同類別所在的位置；

(2)基于河冰分形特征的邊緣檢測結果選取河段分類樣本，利用支持向量機進行智能分類識別，總體分類結果達到93.05%，智能方法分類結果精度較高，且與分形特征邊緣檢測的結果一致，本研究方法相對準確可靠；

(3)隨著高分遙感影像的飛速發展，通過對河冰影像智能分類分析研究，可對凌汛災害進行有效的識別及信息提取，可為凌汛災害風險的預測評價提供理論方法支撐，對提升凌汛災害的綜合防控能力具有重要意義。