王 巍,鄭新奇,原智遠,張路路
(中國地質大學(北京)信息工程學院,北京100083)
地物光譜的多樣性和復雜性,導致傳統分類方法得到的分類圖像不能達到應用目的,必須進行分類后處理工作,從而提高遙感影像分類的精度。初始分類圖像中常見的缺陷是分布于分類圖像上的孤立點、孔洞等,即所謂的“椒鹽現象”或“椒鹽噪聲”,不僅影響遙感圖像分類精度,而且在柵格轉矢量過程中大幅增加了計算量。商業遙感影像處理軟件為用戶提供了基于低通濾波、基于光譜特征或是基于形態學算子的聚類處理方法[1],可以將臨近的、同種類型的小塊分區進行聚類并合并,但存在聚類結果不理想的情況。國內部分學者針對分類圖像中的椒鹽噪聲開展了相關研究,主要解決有兩條途徑:一是針對中、低分辨率影像,多結合地統計學的相關理論,如武文波等提出基于數學形態學的遙感影像分類后處理方法,并在ENDAS平臺上進行了試驗[2];靳志賓等在k-NN分類方法中加入地理權重,能一定程度減少椒鹽噪聲的產生[3];姚娜等將地統計學理論引入遙感臨域,利用臨近像元所蘊含的空間依賴性來改善初始分類結果[4];基于地統計學的方法一般都需要訓練樣本去構建分類后處理的模型,因此其效果必然會受到樣本容量的影響。另外一種解決途徑是采取面向對象分類,如殷亞秋、李家國等借助易康軟件,通過面向對象分類的方式進行水體信息提取,從而在很大程度上抑制了椒鹽噪聲[5],但面向對象方法更多適用于高分辨率影像,針對中、低分辨率影像則偏重于像素級別的處理。
遙感分類圖像的本質是對地物信息的一種數值模擬,各像素間具備空間依賴性,即對于圖像某一區域而言,分散在一定距離內的像素點在數值特征上通常會呈現空間依賴性[6]。因此,針對中分辨率遙感影像,筆者提出一種基于鄰域規則的遙感分類后處理方法,該方法基于地理學第一定律的思想,借鑒元胞自動機中Moore(摩爾)鄰域這一概念,構建由中心像素及其周圍相鄰8個像素組成的“像素元胞”,嘗試對“元胞”鄰域制定聚類規則來實現對分類圖像的分析和調整,力求去除分類圖像的椒鹽噪聲。
地理學第一定律來自于美籍瑞士地理學家Waldo Tobler的觀點:任何事物都與其他事物相聯系,但臨近的事物比較遠的事物聯系更為緊密[7-8]。地理學第一定律表明地理事物或屬性在空間分布上具有相關性,遙感分類圖像是對現實世界的數值模擬。根據地理學第一定律,理論上可以對分類圖像采用基于鄰域規則的聚類方式進行分類后處理,其重點在于規則的制定。
元胞自動機(CA)是一種時間、空間、狀態都離散,空間相互作用和時間因果關系都為局部的網格動力學模型[9]。摩爾鄰域的概念來自于元胞自動機。標準CA中,摩爾鄰域由中心元胞周圍相鄰8個元胞組成。如果說標準摩爾鄰域的擴展半徑為1,那么擴展型摩爾鄰域則以中心元胞開始,往外層進行擴展,例如半徑為2的擴展型摩爾鄰域中會存在5×5的結構,即存在25個“元胞”(如圖1所示)。將這種結構放到分類圖像上,表現為一種3×3或5×5像素大小的一塊區域,聚類規則便基于圖像的摩爾鄰域這種結構而制定。需要注意的是,以中心像素向外擴展過程中,在圖像的邊緣無法獲得完整的摩爾鄰域。例如標準摩爾鄰域由中心像素向周圍擴展一次,所獲取的摩爾鄰域可能不足9個像素,而是僅有4個或6個像素。

圖1摩爾鄰域(左)和擴展摩爾鄰域(右)
由于分類圖像的DN值代表地物類型,在研究中,對于分類的圖像主要根據其DN值進行處理,可以忽略決定顏色顯示的RGB值。同時,為了提高運算效率,在對中心像素DN值的檢測過程中,檢測區域設置為標準摩爾鄰域。在進行分類圖像調整時,主要遵從以下幾條既定的鄰域規則。
(1)眾數規則
在一個摩爾鄰域內,某一DN值的頻數較高,則將此DN值設置為該摩爾鄰域的眾數。眾數規則對中心像素DN值的影響表現為中心像素DN值傾向于等于該摩爾鄰域的眾數。
(2)最鄰近規則
中心像素與其臨近像素DN值相等的可能性高,而與距離較遠的像素DN值相等的可能性低。最鄰近規則對中心像素DN值的影響表現為中心像素DN值傾向于等于與其臨近像素的DN值。
(3)像素權重規則
研究中對分類圖像的所有像素采用循環遍歷的方式進行處理。在循環過程中,一部分像素已經經過分析處理,而另一部分像素則還未進行分析處理。因此,像素權重規則在這里主要指進行規則聚類時,經過分析處理的像素對中心像素的影響程度要高于未經過分析處理的像素。在摩爾鄰域內,此規則體現出兩個方面的含義:一是前一像素DN值權重高于后一像素;二是上一行像素DN值權重高于下一行像素。像素權重規則對中心像素DN值的影響表現為中心像素DN值傾向于等于與其臨近的已經過分析處理的像素DN值。
(4)規則的差異性
前面提到過,圖像邊緣無法獲取完整摩爾鄰域。因此,在對分類圖像進行循環遍歷時,根據像素在圖像上的不同部位,上述規則需要做一定程度的修改。規則的差異性主要有3類:一是在圖像的4個角點處,此時摩爾鄰域僅有4個像素;二是在圖像邊緣處,此時摩爾鄰域有6個像素;三是循環至圖像的內部,才能出現完整摩爾鄰域,即9個像素。3種摩爾鄰域情況如圖2所示,圖中黑色部分表示中心像素,灰色部分表示摩爾鄰域。

圖2 圖像上不同位置的摩爾鄰域示意圖
試驗數據為一景成像于2013年10月的北京地區Landsat8中分辨率遙感影像,在ENVI5.1平臺下對Landsat8影像進行數據預處理,過程包括輻射標定和大氣校正;經過預處理的影像數據占用較大的硬盤空間,處理時間較長,需要進行試驗區的挑選和影像裁剪。通過瀏覽北京市影像發現,密云地區地類種類較全,包含了植被、水體、土壤以及不透水面4類主要地類,適用于本次試驗。通過矩形框對經過預處理的Landsat8影像進行裁剪,矩形裁剪框的范圍是北緯40°30'14″至40°19'46″,東經116°40'50″至116°59'27″。研究區在北京市的地理位置如圖3所示。

圖3 研究區在北京市的地理位置
城市遙感中,許多學者借助一種“軟”分類的思想,提出一種V-I-S分類模型用于城市景觀的分析[10-11]。V-I-S指的是植被(vegetation)、不透水面(impervious surface)以及土壤(soil)。通過對試驗區影像解譯可以發現,水體在影像中占有一定的比例。因此,研究中對V-I-S模型進行了擴展,添加水體類型,構建V-I-S-W分類模型,并在ENVI平臺中進行監督分類,分類方法選擇最大似然法。獲取的分類圖像包含4類地物信息:植被、土壤、水體和不透水面。
借助ENVI平臺實現初始分類圖像的獲取,主要步驟如下:
(1)訓練樣本的選擇
使用ENVI提供的ROI工具進行訓練樣本的選取。表1提供的數據表明,樣本可區分度均在1.80以上,滿足軟件計算要求。

表1 訓練樣本的可區分度
(2)基于最大似然法的分類
使用ENVI中Maximum Likelihood工具進行最大似然法分類,設定保存路徑和文件名,其他參數采取默認設置。運行工具后,軟件自動進行計算,并將初始分類圖像保存到設置路徑中。
本文在IDL平臺下實現文中提出的聚類規則。選擇IDL的原因,是因為IDL不僅提供了功能強大圖像處理函數庫,而且與ENVI軟件高度集成,使其功能得到進一步加強。通過IDL實現的程序包含下列功能:
(1)分類數據的輸入和輸出
ENVI標準格式的分類圖像,其DN值、地物類型及顯示顏色三者間的關聯信息都寫在圖像的頭文件中,而分析處理重點是僅包含DN值的二維圖像(矩陣),可以在IDL批處理模式中,調用ENVI內置函數進行相關數據的讀取操作;寫出數據時,ENVI和IDL也提供了多種寫出方式,可以將數據文件和頭文件分開輸出到本地儲存,也可以通過ENVI函數將數據文件和頭文件一并寫入本地存儲。表2中列舉了IDL代碼中用到的數據讀寫函數及其功能。

表2 IDL批處理模式下用到的數據讀寫函數
(2)眾數的求取
分類圖像DN值的數值類型為長整形。IDL提供的Histogram()函數來實現直方圖統計的功能,根據直方圖就可以計算數組的眾數。
在IDL編譯器中,求取眾數的實現方式如下面代碼所示:
FUNCTION CAL_MODE,data
ht=HISTOGRAM(data,OMIN=omin)
in=WHERE(ht NE 0)
maxv=MAX(ht[in],index)
RETURN,omin+WHERE(ht EQ maxv)
END
(3)構建摩爾鄰域
對IDL二維數組中元素的訪問方式除了行列序號外,借助IDL數組的位置索引機制也是一種非常便利的訪問手段。IDL中規定:索引從0開始計數,由數組左上角至右下角,按列逐漸遞增,因此二維數組最后一位元素的索引為二維數組元素個數減1。構建摩爾鄰域時,借助索引機制,可以方便地獲取摩爾鄰域內像素的位置和DN值信息。
下面以圖像的內部的摩爾鄰域為例說明(如圖4所示)。首先為摩爾鄰域進行編碼,序號從1—8,9號位是中心像素,用編碼表示以中心像素向外進行摩爾擴展時可獲取的像素位置信息;然后根據圖像矩陣的行列號設定轉換規則,將位置編碼映射到圖像矩陣像素的位置索引上;最后通過位置索引信息便可以獲取摩爾鄰域內所有像素的DN值。通過這種映射方式,只需獲取位置編碼、索引或DN值三者其中之一,就能查詢其他兩者的信息。

圖4 摩爾鄰域構建及映射關系
(4)情景推導
程序載入分類圖像后,通過對圖像所有像素循環遍歷的方式來調整圖像,根據既定的規則逐一進行處理,每次僅調整一個像素,筆者將這一過程命名為情景推導。表3是通過偽代碼形式展示情景推導過程,其書寫方式基于IDL語法規則,因篇幅原因,更加復雜的推導過程以Rule_Set替代。部分推導過程中,需要將待檢測中心像素周圍的像素作為新的中心像素進行規則判別,這一過程類似于函數遞歸調用,因此以RECURSION替代。另外,IDL中沒有In運算符,查詢功能在IDL代碼中可以通過Where()函數實現;而or在IDL中是位運算符,實現條件選擇功能可以通過if-else組合來完成。

表3 情景推導過程
為了對比分析,研究中對初始分類圖像應用ENVI軟件的Clump Classes工具進行了聚類處理,形態學算子大小設置為3×3,與摩爾鄰域大小相同。圖5顯示了影像不同區域上原始圖像、初始分類圖、Clump Classes工具處理結果以及本研究提出方法處理結果的對比,初始分類圖像經過基于鄰域規則的聚類處理后,椒鹽噪聲的抑制效果較為顯著,且視覺上較Clump Classes工具的處理結果而言更加接近于原始影像呈現的狀態。
分類精度方面,因缺少地表真實信息,在統計分類精度時采用人工解譯獲取地表真實ROI的方式,將ROI作為地表真實信息。用于精度驗證的ROI數量分別是:植被668像素、土壤653像素、水體666像素以及不透水面693像素。所有操作均基于ENVI平臺,精度統計結果見表4。通過結果可以看出,本次試驗中,基于鄰域規則的聚類方法處理的分類圖像精度最高,達到95.63%,Kappa系數為0.941 8。

表4 分類精度統計

圖5 部分區域影像及分類結果對比
通過試驗結果可以看出,應用基于鄰域規則的聚類方法對初始分類圖像進行分類后處理,可以在很大程度上抑制椒鹽噪聲,提高分類精度,滿足了實際應用的需求。同時,使用地表真實ROI的精度驗證結果表明,基于鄰域規則的聚類結果較之ENVI軟件提供的Clump Classes工具的處理結果而言精度要更高,形態上更加接近于影像所呈現的狀態;但是在地塊邊緣處理上,Clump Classes工具給出的結果更加平滑,能更好消除地塊的鋸齒狀邊緣。條件所限,兩者都缺乏真實地物信息的檢驗,因而結果還有待進一步驗證。
基于鄰域規則的聚類方法關鍵點在于對規則的制定,規則越合理、推導越深入、復雜度越高,處理的結果應該更加精確;如果加入對地塊邊緣的處理規則,也能在一定程度上消除鋸齒狀邊緣。在本研究中,規則設計是建立在對圖像所有像素進行循環遍歷處理的基礎之上,隨著規則復雜度的提高、計算量的增大,處理效率必然會受到影響。因此,如何設計行之有效、效率更高、適用性更強的規則將是未來的研究重點。
[1] 鄧書斌.ENVI遙感影像處理方法[M].北京:科學出版社,2011.
[2] 武文波,楊貴軍,陳步尚.基于數學形態學遙感影像分類后優化處理[J].遼寧工程技術大學學報(自然科學版),2003,22(2):172-175.
[3] 靳志賓,蒲英霞,陳剛,等.基于地理加權的k-NN高分辨率遙感影像分類算法改進[J].遙感技術與應用,2013(1):97-102.
[4] 姚娜,林宗堅,張錦雄,等.遙感影像分類后處理的地統計方法[J].武漢大學學報(信息科學版),2013,38(1):15-18.
[5] 殷亞秋,李家國,余濤,等.基于高分辨率遙感影像的面向對象水體提取方法研究[J].測繪通報,2015(1):81-85.
[6] 張錦雄,GOODCHILD M F.野外空間采樣的漸進式策略[J].武漢大學學報(信息科學版),2008,33(5):441-445.
[7] TOBLER W.On the First Law of Geography:A Reply[J].Annals of the Association of American Geographers,2011,94(2):304-310.
[8] 孫俊,潘玉君,和瑞芳,等.地理學第一定律之爭及其對地理學理論建設的啟示[J].地理研究,2012(10):1749-1763.
[9] 黎夏,葉嘉安,劉小平,等.地理模擬系統:元胞自動機與多智能體[M].北京:科學出版社,2007.
[10] SETIAWANA H,MATHIEU R,et al.Assessing the Applicability of the V-I-SModel to Map Urban Land Use in the Developing World:Case Study of Yogyakarta,Indonesia[J].Computers,Environment and Urban Systems,2006,30(4):503-522.
[11] WENG Qihao,LU Dengsheng.Landscape as a Continuum:an Examination of the Urban Landscape Structures and Dynamics of Indianapolis City,1991—2000,by Using Satellite Images[J].International Journal of Remote Sensing,2009,30(10):2547-2577.