多/高光譜遙感數據的類立體紋理特征

趙 巍，崔淑梅，吳 銳，劉家鋒

(1.哈爾濱工業大學計算機科學與技術學院，150001哈爾濱;2.哈爾濱工業大學電氣工程及自動化學院，150001哈爾濱)

源于灰度共生矩陣(Gray Level Co-occurrence Matrix，GLCM)的紋理特征30多年前已被證明其在分析和辨識遙感影像多種不同地物方面具有有效性［1-2］.近年來研究者［3］發現對衛星圖像分類過程中采用紋理特征是克服光譜分析缺失的有效措施.如今遙感圖像研究者［4-7］正以各種方式使用并改善著紋理這一經典的、可塑性高的遙感數據基元描述工具.

目前針對遙感數據的處理通常是對原始圖像提取邊緣、紋理等基本幾何元素或特征.多/高光譜遙感數據的同源同點多幅圖像同時獲取的性質使其具有非常豐富的協同數據可利用;同時，分析遙感圖像數據的有效特征可能就是獲取圖像時所設定的多個光譜的波長或波長范圍［8］，因此研究者［9］對紋理等圖像特征的提取很早就開始向多維化、立體化的方向發展.但目前的解譯方法仍缺乏充分挖掘影像各種空間特征和地物其他屬性特征等方面的探索［10-11］.文獻［9］在進行紋理變換和分析時采用通道信息組合的方式計算3×3像素窗的中心像素點的局部屬性，認為多通道信息融合對邊緣檢測和圖像增強的效果要好于在地物主題映射的效果.Mounir等［12］從16個子波段、3個通道中分別抽取熵、紋理均值、3階矩、4階矩和能量5個參數組成特征向量.這種紋理特征構造方法可以增加信息量并提高分類的精度［13］.此外，采用均值法和主成分分析法，將子圖像各個波段基于GLCM的紋理統計特征向量作為多波段遙感圖像的特征向量［14］，以挖掘同源不同尺度信息，這樣雖增加了紋理分析的復雜性，但分類性能均有不同程度的增益.上述同源遙感數據的紋理研究中較少有直接考慮波段之間像素灰度變化的關聯關系的做法，多是分別提取單波段數據的紋理特征，再在特征層或決策層實現各波段信息融合的方式.比較而言，在多源遙感數據分析中的像素級融合的情況則非常普遍［15-17］.

由于GLCM是對單幅圖像遍歷得到的統計矩陣，這可能限制其在多/高光譜遙感解譯分析能力.本文根據多光譜遙感數據的同源同點多幅同時獲取的特點，提出了一種多光譜遙感圖像底層數據綜合方法，即灰度級差關聯概率矩陣(Gray Level Difference Associated Possibility Matrix，GLDAP)描述.目前，能夠采集立體信息的成像方式大致有兩種:雙目成像和全息成像.對于兩幅同源不同波段的遙感圖像來說，顯然具有比單幅圖像更豐富的視覺信息，當將其中一幅作為參考時，其與另一幅圖像的差異即被顯著化，灰度級差關聯概率矩陣GLDAP以統計的方式記錄了這些差異.GLDAP方法是對兩幅同源圖像的底層數據進行視覺差異統計，是對立體成像過程的模擬.模擬過程中所獲得的信息被稱為類立體信息，在類立體數據的基礎上得到的紋理特征即為類立體紋理.

1 灰度級差關聯概率矩陣

事實上多/高光譜遙感數據除了二維平面像素點的值之外，沒有提供第3維度量信息，并非真正意義上的立體數據.每個波段圖像可看作是透過不同透鏡觀察到的景物，只是無法可靠地指出哪個波段的圖像更能準確、清晰地給出圖像內容的全部細節信息，只是依賴于單波段圖像提供的局部可靠信息來解譯圖像.人們已經根據地物的波譜特性獲得了對特定類型地物有較好可分性的波段范圍，從而有針對性地獲取特定地物在相應單波段圖像上的信息.

本文依據地物譜特征分析，參考區分特定地物的波段信息，從遙感數據中選擇多幅有效波段圖像，統計每兩幅圖像的灰度的協同變化規律，進行視覺差異統計分析，通過計算差異來刻畫和提取有效紋理特征，以克服GLCM方法對單幅圖像統計描述的局限性，實現數據解譯分析的高精度目標.進行同源多幅圖像處理時，對每次選定的兩幅圖像進行像素信息的協同統計，建立兩圖像在同一地理位置上對應像素點之間灰度值變化的關聯關系，計算像素點對的灰度差值改變的概率，構成GLDAP概率矩陣;基于此矩陣計算統計紋理特征，如能量、對比度及熵等.視覺差異分析方法以類似于立體視覺的方式來處理多/高光譜圖像.

Haralick等［18］提出經典的GLCM紋理分析方法設置窗口遍歷一幅圖像，統計差異并記錄在灰度共生矩陣中，在灰度共生矩陣上提取包括能量、熵在內的14種紋理特征.一幅圖像在0°，45°，90°和135°等4個方向上的灰度共生矩陣的計算為

式中:P為灰度共生矩陣;#為集合中元素的數目;(k，l)、(m，n)分別為像素點坐標;i，j分別為灰度級;d為兩像素之間的距離.

本文在GLCM方法的基礎上提出的GLDAP紋理分析方法，統計同源兩目標波段圖像F和F'中在相同方向、相同像素點距離d上灰度級差分別為Δi與Δi'的像素對的數目;以此來反映兩圖像在相同方向和像素距離上的灰度變化相關性.GLDAP矩陣的定義為

式中:~P為相應方向的灰度級差關聯概率矩陣;#為集合中元素數目;(k，l)、(m，n)和(k'，l')、(m'，n')分別為F和F'中的相應兩點坐標;Δi、Δi'分別為F和F'中相應兩點灰度差值;d為窗口內兩像素間距.

按照GLDAP的定義，應當選定方向和點間距d，根據式(1)同時遍歷兩個波段的圖像并歸一化便可獲得GLDAP矩陣.

GLDAP矩陣是非對稱的，其行列值分別是兩幅圖像的灰度差的絕對值，很好地反映了在兩幅同源圖像的相同距離和相同方向上像素灰度差的相關變化的統計信息.當兩幅圖像灰度分布一致時，主對角線和為1，各非主對角線和為0;隨著灰度分布差異增大和紋理復雜程度增加，主對角線和將隨之減小，各非主對角線和相應增加.因此，對角線和也是一種重要的統計特征.GLDAP矩陣主對角線元素表示兩幅圖的相鄰像素灰度差相同的概率、各非主對角線元素表示兩幅圖的相鄰像素灰度差不同的概率，計算公式為

式中:GLDAP矩陣大小為K×K;i，j分別為矩陣的行、列數;t為所求取的對角線;t=0為主對角線.

2 實驗結果與分析

用GLDAP進行多光譜遙感圖像中的水體、城區、山區3種類型的地物進行分類，并同GLCM的分類結果進行比較.為使實驗具有可比性，實驗過程中為GLCM矩陣方法提供與GLDAP方法相同的兩幅單波段圖像.設提取D種紋理特征，則GLCM方法采用文獻［12］的常規處理方式，分別獲得兩幅圖的紋理特征再構成特征向量，維數為2D;GLDAP方法除了這些紋理特征外還使用了主對角線和特征來反映兩波段具有相同灰度差的統計量，其特征維數為D+1.

設圖像的大小為N×N，預設窗口尺寸、灰度級及特征數等參數為常量，兩種方法的時間復雜度均為O(N2).

實驗系統的決策模塊采用支持向量機(Support Vector Machine，SVM)方法，選擇徑向基函數(Radical Basis Function，RBF)作為核函數.軟件在裝有VC、Open CV環境的PC機上實現.

2.1 樣本與前期實驗準備

采用Landsat-7的增強型專題制圖儀ETM+多光譜數據.第1～5及7波段的圖像如圖1所示，其中B為波段;圖中不包括波段6(熱紅外)數據.單波段圖像灰度級為256，為減少計算量和避免信息損失，將灰度級壓縮為64級.借助ERDAS Imagine遙感圖像處理系統為ETM+多光譜數據準備訓練和測試樣本.

圖1 ETM＋多光譜波段數據圖

水、陸兩類分類的預備性實驗中，進行了第1～5、7波段兩兩組合以及與GLDAP矩陣的不同計算窗口尺寸的組合測試;d取值為1，即只比較相鄰兩像素點的灰度變化;選取能量、對比度、均一性和相關性4個特征;選擇了近1 500個樣本點對支持向量機SVM進行訓練，對近350 km2的水陸共存的遙感圖像進行分類實驗測試.兩類分類的實驗結果顯示，GLDAP有著與GLCM相當的分類能力，將GLDAP應用于遙感圖像的分類是完全可行的;選用第4、7波段組合和第1、7波段組合時，與5×5大小的窗口配合使用，地物邊界效果理想.

基于上述實驗參數條件和實驗平臺，用6個波段的10 973個樣本對SVM分類器進行訓練，截選若干區域進行了水體、城區和山區3類地物樣本的分類測試.選擇兩幅有以上3種地物的區域作為本文實驗數據處理的演示圖像，其合成ETM+數據的第3、2、1波段組合為模擬真彩色圖像如圖2，3所示.肉眼觀測該圖可分辨城區紋理細密，山區紋理相對粗糙，水體灰度分布均勻、紋理單一.

圖2 ETM＋多光譜波段數據第3、2、1波段合成示意圖

圖3 ETM＋多光譜波段數據第3、2、1波段合成示意圖

2.2 分類實驗結果樣例與分析

對遙感圖像進行水體、城區和山區的3類分類，選用第4、7波段組合和第1、7波段組合，配合5×5的滑動窗口.分類結果用黑、白、灰3色標記:將水體標記為白色，城區(人工跡地)為黑色，山區為灰色.

2.2.1 第4、7波段圖像分類的結果

對圖2的第4、7波段數據處理結果如圖4所示.圖4(a)是用GLCM方法的類別標記結果，能夠大略看出山區的輪廓，但誤識多.如:將部分山區地貌誤識為城區，部分區域因紋理較單一被誤識為水體;城區部分綠地與山區紋理很相近而被誤識為山區;水體中間部分河床較高而水體很淺，被誤識為山區.

圖4(b)是用GLDAP方法對山區附近圖像的分類結果，可以清晰地看出山區輪廓，山溝之間的居民區能正確地標記為城區，水體中間的淺灘部分也能正確識別.城區的綠地等部分由于紋理較粗糙，也存在誤分類情況，但明顯少于GLCM方法的誤分量.這兩種方法對水體部分的識別性能差異不大.

圖4 圖2的第4、7波段采用GLCM方法與GLDAP方法分類處理結果

對圖3的第4、7波段數據處理結果如圖5所示.圖5(a)為GLCM方法的分類結果，山區和城區兩類地物都互相存在大量的誤分類，幾乎無法分辨;水體兩岸較為接近的區域誤分成山區，而右下角灰度值較單一的山區部分則被誤分為水體.使用GLCM的誤識情況更為嚴重.

圖5(b)為GLDAP分類結果，從圖中可明顯看到，山區和城區能夠很清晰地分辨，且分類準確性很高.

圖5 圖3的第4、7波段采用GLCM方法與GLDAP方法分類處理結果

2.2.2 第1、7波段圖像分類的結果

采用第1、7波段圖像分類的結果分別如圖6、7所示.從圖6、7可以看出，使用GLCM紋理分類時水體和山區的區分度很差，水體的大部分被誤分成了山區，而山區的大部分識別成了水體.造成山區和水體互有誤識的情況主要是:1)水體的反射主要在第1、2波段，由于水體的反射率較高，第1波段圖像水體的顏色較淺，而且受水體的深淺和懸浮物質的影響第1波段的紋理較第4或第7波段的復雜;水體對其他的波段吸收都很強，特別是第4波段，因此第4波段的水體圖像呈黑色，且紋理單一;2)山區地表主要是綠色植被，綠色植被在藍光區域(第1波段)有一個吸收谷，因此山區呈現的顏色也較深，紋理也較單一，與水體較接近;導致互相誤識的情況大量發生.而由于GLDAP方法統計了兩個波段間的灰度變化的關系，因此受單波段圖像的影響較小，第1、7波段組合3種類型的區域基本上能夠正確識別，山區的正確識別率較GLCM方法有大幅提高.

可見，選用GLDAP處理多波段圖像時能夠避免因單波段中某兩類地物的可分性差而造成的誤識.

圖6 圖2的第1、7波段采用GLCM方法與GLDAP方法分類處理結果

圖7 圖3的第1、7波段采用GLCM方法與GLDAP方法分類處理結果

2.3 統計數據及性能分析

兩種方法分類效果的定量比較及平均分類正確率統計如表1所示.平均分類正確率是從樣本圖像中隨機選取地物的若干典型樣本區，計算每類地物正確分類像素數與實際該類像素數的百分比，將各類正確分類率取平均得到的.

表1 地物識別平均正確率%

實驗表明，采用第1、7波段組合和第4、7波段組合，GLDAP對多光譜遙感圖像進行分類比傳統的GLCM方法更準確;第4、7波段組合的GLCM和GLDAP紋理方法優于第1、7波段組合的結果.

實驗還分別統計了3種地物互有誤識時所對應的6種平均誤識率(如水體誤識為城區、城區誤識為山區的比率等).表1中第1、7波段組合實驗GLCM方法的水體(14.4%)、山區(13.3%)的識別率極低，而城區識別率卻高達98.3%.表2給出了兩種方法的被誤識為城區的平均統計量，從中可看出，GLCM將水體和山區誤識為城區的比率之和高達80.8%，而GLDAP方法在誤識為城區的比率則較該方法其他類型的誤識差異不大.

表2 第1、7波段組合誤識為城區的統計量比較%

同時還看到，與第4、7波段組合相比，第1、7波段組合時水體被誤識得最為嚴重.這與特定波段反映地物譜信息的能力不同直接相關，導致了所構建的波段組合形成的視覺差別和基于此視差獲得的統計信息在反映地物能力上的強弱有別.在對多/高光譜數據進行研究時，配合波段相關性分析，采用GLDAP方法充分測試不同波段組合，將對高效、完整獲取類立體數據的可利用信息起到至關重要的作用.

3 結論

1)本文從各波段數據的內在關聯出發，根據地物的波譜特征，選擇雙波段統計灰度協同變化規律，建立視覺差異統計方法，基于GLDAP矩陣提取類立體紋理特征，以克服GLCM方法對單波段圖像統計描述的局限性;用類立體視覺方式處理多/高光譜圖像，以提高遙感數據解譯分析精度.方法同樣可用于高光譜數據分析.

2)實驗表明，GLDAP紋理提取方法具有良好的分類性能，3種地物分類效果明顯優于GLCM方法，并能在一定程度上避免GLCM因單波段中兩類地物可分性差而導致誤識的情況.兩種方法的時間復雜度相同，為多/高光譜遙感數據解譯人員在相同時間條件下獲得更優良的性能提供了新的可操作性強的分析途徑.

［1］WESZKA J S，DYER C R，ROSENFELD A.A comparative study of texture measures for terrain classification［J］.Comparative and General Pharmaclology，1976，1(4):269-285.

［2］SHANMUGAN K S，NARAYANAN V，FROST V S，et al.Textural features for dadar image analysis［J］.IEEE Transactions on Geoscience and Remote Sensing，1981，19(3):153-156.

［3］CARLEER A P，WOLFF E.Urban land cover multilevel region-based classification of VHR data by selecting relevant features［J］.International Journal of Remote Sensing，2006，27(6):1035-1051.

［4］OHANIAN P P，DUBES R C.Performance evaluation for four class of texture features［J］.Pattern Recognition，1992，25(8):819-833.

［5］NYOUNGUI A N，TONYE E，AKONO A.Evaluation of speckle filtering and texture analysis methods for land cover classification from SAR images［J］.International Journal of Remote Sensing，2002，23(9):1895-1925.

［6］RAJESH K，JAWAHAR C V，AENGUPTA S，et al.Performance analysis of textural features for characterization and classification of SAR images［J］.International Journal of Remote Sensing，2001，22(8):1555-1569.

［7］舒寧.關于多光譜和高光譜影像的紋理問題［J］.武漢大學學報(信息科學版)，2004，29(4):292-295.

［8］RICHARDS J A，JIA X P.Remote sensing digital image analysis，an introduction［M］.Berlin，Heidelberg:Springer-Verlag，2006.

［9］JAMES R I，GARY W P.Texture transforms of remote sensing data［J］.Remote Sensing of Environment，1981，11:359-370.

［10］舒寧.關于遙感影像處理分析的理論與方法之若干問題［J］.武漢大學學報(信息科學版)，2007，32(11):1007-1015.

［11］駱劍承.遙感影像智能圖解及其地學認知問題探索［J］.地理科學進展，2000，19(4):289-296.

［12］MOUNIR A K，AHMED H，DRISS A.Textural feature selection by joint mutual information based on gaussian mixture model for multispectral image classification［J］.Pattern Recognition Letters，2010，31(10):1168-1174.

［13］TANG X，STEWART W K.Optical and sonar image classification:wavelet packet transform vs fourier transform［J］.Computer Vision and Image Understanding，2000，79(1):25-46.

［14］田艷琴，郭平，盧漢清.基于灰度共生矩陣的多波段遙感圖像紋理特征的提取［J］.計算機科學，2004，31(12):162-163.

［15］盧焱.多源遙感影像融合方法研究［J］.吉林大學學報(地球科學版)，2007，37(11):176-179.

［16］賈永紅，李德仁.多源遙感影像像素級融合分類與決策級分類融合法的研究［J］.武漢大學學報(信息科學版)，2001，26(5):430-434.

［17］許君一，盧秀山，卿熙宏，等.基于熱傳導模型的像素級遙感圖像融合［J］.測繪學報，2010，39(6):566-571.

［18］HARALICK R M，SHANMUGAM K，DINSTEIN I.Textural features for image classification［J］.IEEE Transactions on System，Man and Cybernetics，1973，SMC(3):610-621.

［19］BARALDI A，PARMIGGIANI F.An investigation of the textural characteristics associated with gray level Cooccurrence matrix statistical parameters［J］.IEEE Transaction on Geosciences and Remote Sensing，1995，33(2):293-304.