分塊低秩圖的遙感影像半監督分類應用*

祖寶開，夏克文+，牛文佳，姜曉慶

1.河北工業大學 電子信息工程學院，天津 300401

2.河北工業大學 河北省大數據計算重點實驗室，天津 300401

3.濟南大學 信息科學與工程學院，濟南 250022

1 引言

近年來，隨著科學技術的飛速發展，高光譜、高分辨率等遙感影像也逐步發展。它已經廣泛應用于礦產資源勘探[1-2]、軍事[3-4]、城市土地利用[5-6]、作物生長監測和生態環境檢測[7-8]等方面。高分辨率圖像具有較高的空間分辨率和豐富的空間信息。然而，數據量大、紋理特征復雜等問題很難解決。高分辨率圖像數據的處理和信息提取一直是困擾研究者的難題[9]。

隨著社會的發展，城市土地利用研究顯得尤為重要[9]。城市土地利用變化對人類社會活動、城市經濟發展和城市生態環境都是有影響的。機器學習方法是遙感圖像分類中最有效的方法之一。主要分為無監督和監督分類方法[10]。監督分類利用先驗知識來訓練分類模型。與此相反，無監督分類不利用先驗信息，通常也被稱為聚類[10]。常用的無監督分類有k-均值和迭代自組織數據分析算法（iterative selforganizing data analysis techniques algorithm，ISODATA）[11]。監督分類方法有統計模型的分類方法、基于神經網絡的分類方法[12]、支持向量機（support vector machine，SVM）[13]和相關向量機（relevancevectormachine，RVM）[14]分類方法。決策樹方法包括基于信息熵增益的ID3算法[15]和利用信息熵增益比的C4.5算法[16]。

然而，大規模遙感數據難以獲取大量標記數據。因為標記樣本需要大量的人力和物力，且時間長，效率低。如果只利用少量的標記樣本，訓練出的學習系統很難有好的泛化能力[17]。而且，如果只使用標記樣本來訓練分類模型，則浪費了未標記樣本所提供的有用信息。為了解決上述問題，半監督學習應運而生。半監督學習利用少量的標記樣本和大量的未標記樣本來學習分類模型。與監督方法相比，半監督學習可以有效地緩解“小樣本”存在的問題[18-20]。Cai等人在線性判別分析的基礎上提出了同時使用標記和未標記樣本的方法，稱為半監督判別分析（semisupervised discriminant analysis，SDA）[21]。但是基于圖的半監督學習算法的性能在很大程度上依賴于圖的構建過程。傳統的圖構造方法K近鄰（K-nearest neighbor，KNN）和局部線性嵌入（locally linear embedding，LLE）等方法在很大程度上依賴參數的選擇，同時利用歐幾里德距離作為測度，容易受到噪聲的影響。

Candès等人為了克服維數災難，提出了稀疏表示與低秩理論[22-23]。稀疏表示和低秩理論的提出是為了解決信號重構的問題。低秩表示（low rank representation，LRR）法構造了一個無向圖（這里稱為低秩（low rank，LR）圖），它根據低秩約束條件聯合獲取整個數據表示[24-25]，對噪聲來說是魯棒的。

低秩表示在遙感圖像分類上的應用還很少，技術還不是很成熟，尤其是針對大規模遙感圖像的分類。隨著數據點（像素）的增加，LRR的計算復雜度呈指數增長。為了解決上述問題，將處理后的圖像按像素劃分為塊（子集），并在每個塊上進行低秩表示。然后，將特征表示子集組合成一個完整的特征圖。應用K近鄰來處理整個低秩圖使之轉換為對稱矩陣，同時，K近鄰最大程度地保留了圖像的局部信息[26]。最后，應用半監督判別分析進行特征提取，可以充分利用標記樣本的標簽信息以及標記和未標記的數據點數據固有的幾何結構信息。

對武漢市區WorldView-2遙感影像進行仿真分析。結果表明，基于塊低秩圖的特征提取方法是一種有效的方法。可以大大提高遙感圖像分類的性能。在簡單的最近鄰分類器以及少量標簽樣本的情況下，也可以獲得顯著的分類性能。

2 基本理論和相關概念

2.1 低秩表示模型

LRR是將數據矩陣X=[x1,x2,…,xn]∈Rm×n表示成字典矩陣A=[a1,a2,…,am]（也稱為基矩陣）下的線性組合，即：

其中，Z=[z1,z2,…,zn]∈Rl×n為低秩表示矩陣，zi為xi的低秩表示系數。為此，求解下列優化問題：

直接求解式（2）較困難。通過矩陣分解方法，式（2）轉化為核范數問題：

其中，‖?‖*表示矩陣的核范數。

當用樣本X本身作為字典，式（3）可轉換為：

考慮到實際應用中數據的損壞或丟失和噪聲等，目標函數轉化為：

其中，噪聲項選取l2,1范數。λ(λ＞0)平衡目標函數的兩部分。

將式（5）轉化為：

式（6）可以通過求解非精確增廣拉格朗日乘子問題計算，即：

式中，Y1、Y2為拉格朗日乘子，μ＞0為懲罰參數。取得最優解（Z*,E*）后，利用Z*構造一個無向圖的相似度矩陣，其中樣本xi和xj為圖的頂點，max([z*]ij,[z*]ji)為兩個樣本之間的權值。

2.2 半監督判別分析

線性判別分析算法（linear discriminant analysis,LDA）目標是找到一個投影子空間，使得在該子空間中，屬于不同類別的點離得越遠越好，而相同類別的點離得越近越好，即使類間離散度和類內離散度的比值達到最大[21]。

假設給定樣本{x1,x2,…,xl}∈Rn，屬于c個類別。LDA的目的是找到一個投影向量a，使aTSba和aTSta之間的比例最大化。

其中，μ是總樣本的均值向量，是第k類樣本的平均向量，lk是k類樣本個數，是在第k類中的第i個樣本，Sw和Sb分別為類內和類間散度矩陣。總散度矩陣,有St=Sw+Sb，那么式（8）中的目標函數為：

最優解a是與本征值問題的非零特征值對應的特征向量：

當訓練樣本不足時，過擬合將發生。防止過擬合的典型方法是加入正則化約束。Cai提出了半監督判別分析模型[21]。SDA的目的是找到一個投影來代表從標記的數據點中推斷出判別信息，以及從標記和未標記的數據點中推斷出的數據固有幾何結構。具體而言，結合標記樣本和無標記樣本建立一個包含數據集鄰域信息的圖。利用拉普拉斯圖的概念，該圖提供了一個相對于數據流形局部幾何結構的離散的近似值。因此，SDA優化算法可以保留流形結構[21]。SDA優化問題為：

其中，正則化項J(a)控制模型的學習復雜度，而系數α控制模型復雜度與正則化項間的平衡[21]。

J(a)是樣本的鄰接圖矩陣，充分利用了先驗知識。半監督學習算法的關鍵是一致性的先驗假設[21]。傳統的SDA用一個K近鄰的圖G模擬附近的數據點之間的關系。對于樣本點xi和xj，在節點i和j之間設置一個邊。相應的權重矩陣S定義為：

其中，Nk(xi)表示K近鄰的集合。一般而言，映射函數在圖上應該是盡可能光滑的。因此，那些稠密相連的子圖可能擁有相同的標簽。正則化矩陣可以定義如下：

其中，D是一個對角矩陣，它是S的列（或行，因為S是對稱的）的總和。Dii=∑jSij，L=D-S是拉普拉斯矩陣。得到SDA的目標函數：

最優投影向量a由廣義特征值問題求得：

其中Sb=XWl×lXT。

3 分塊低秩圖的遙感影像半監督分類

3.1 分塊低秩圖的半監督判別分析

在上面的章節中，研究了經典半監督判別分析和可以有效捕捉數據全局結構的低秩表示算法。基于上述啟發，提出基于塊低秩表示的正則化圖。

對于特征向量X=[x1,x2,…,xM]∈RN×M,N是波段數（數據維數），M是像素數（樣本數）。隨著樣本數量的增加，低秩表示計算時間呈指數增長。因此，通過分塊理論來研究低秩表示問題。

讓{g1,g2,…,gm}表示m個塊圖像的集合索引，其中每個gi存在S個像素。將同一塊的向量放在一起，即X={Xg1,Xg2,…,Xgm},A={Ag1,Ag2,…,Agm}和E={Eg1,Eg2,…,Egm}。然后每個塊的LRR優化問題轉換為以下形式：

選擇l2,1-norm作為誤差項范數‖·‖l，它是

通過不精確增廣拉格朗日乘子（augmented Lagrange multiplier，ALM）法[25]將增廣拉格朗日函子最小化得到塊低秩（block low rank，BLR）圖然后將這些特征子集合并得到特征圖矩陣Z={Zg1,Zg2,…,Zgm}=[z1,z2,…,zM]∈RS×M。

由于圖鄰接矩陣Z是不對稱的，需要對稱化處理才能滿足SDA正則圖的對稱性要求。經典的對稱化過程是矩陣轉置和其本身的均值。然而，鄰接圖矩陣Z不是一個方陣。之前的研究表明，低秩K近鄰（low rankK-nearest neighbor，LRKNN）圖可以顯著提高SDA[26]的性能。因此，K近鄰算法的處理不僅解決了鄰接矩陣對稱化的要求，而且最大限度地保留了圖像局部信息。因此，引入基于塊低秩K近鄰正則圖對遙感影像進行半監督特征提取。

如果zi屬于zj的k個最鄰近點，或者zj是zi的k鄰近點，那么zi和zj互為近鄰點。這里使用熱核方法為S分配權重，如下所示：

其中，Nk(zi)表示zi中的k近鄰。

3.2 基于塊低秩圖的遙感影像半監督分類

給定一個屬于c類的標記數據集和一個未標記數據集其中，第k類有lk個樣本，假設{xi,…,xl}中的數據點根據自身的標簽來排序。

基于塊低秩圖的半監督判別分析的步驟如下：

步驟1構造鄰接圖：構建X的低秩近鄰圖矩陣S，計算該圖的拉普拉斯矩陣L=D-S。

步驟2構建標記圖：為標記圖構建權重矩陣Wl×l∈Rm×m，如下：

步驟3特征值問題：計算廣義特征向量問題及其對應的非零特征值。

其中，X=[x1,x2,…,xl,xl+1,…,xm]，W的秩為d，得到d個特征向量a1,a2,…,ad。

步驟4SDA特征映射：讓A=[a1,a2,…,ad]，A是一個n×c的變換矩陣。通過X→Z=ATX，數據可以映射到d維子空間。

由式（18）可以得到投影矩陣A：

其中，S是KNN變換之后的圖。L=D-S是拉普拉斯矩陣。Sij是代表附近數據點關系模型的圖形。Dii=∑jSij是對角矩陣D。

將W通過映射矩陣Φ映射到d維子空間中。

步驟5最后，執行K近鄰分類器進行分類。

4 實驗仿真與分析

實驗的硬件平臺是處理器為Intel Core CPU 2.60 GHz和8 GB RAM的機器，在Matlab 2016b上進行實驗分析。

4.1 研究區概況

實驗區為湖北武漢市轄區，地理位置在30°52′N～30°48′N與114°27′E～114°32′E之間（圖1）。實驗數據為WorldView-2影像數據，圖像大小為2 500×2 500個像素的tif格式文件，影像質量良好，分辨率為2 m。WorldView-2影像是由美國Digital Globe公司發射的衛星獲取的，該數據具有很高的分辨率。該衛星系統能夠采集海量的影像數據，存儲速度極快，而且定位精度極高，存儲量極大。數據實驗區內地物類型豐富，主要有河流、湖泊、建筑物、草地、農業用地、林地等。在武漢市中心城區附近，景觀分布較為破碎，裸地和建筑用地、草地和林地等在空間上分布混亂。

Fig.1 Study area圖1 研究區域

4.2 精度評價指標

分類精度評價非常重要，通過精度評價，能夠直觀、有效地獲取分類結果中的信息。下面給出幾種常用的遙感影像分類精度評價指標。

在遙感分類領域，經常使用混淆矩陣（confusion matrix）表示分類精度。對角線為正確分類的樣本數，而非對角線表示錯誤分類的樣本數。對角線的值越高表示結果越好。

分類精度（classification accuracy，CA）是指每個類別中正確分類的像素的百分比。總體準確率（overall accuracy，OA）指正確分類像素的百分比。平均準確率（average accuracy，AA）是計算所有類平均正確分類像素的百分比。為了使測量更加客觀，使用Kappa系數來估計正確分類的像素的百分比。

制圖精度（producer’s accuracy）是指分類器將整個影像的像元正確分為某類的像元數（混淆矩陣中對角線值）與該類真實參考總數（混淆矩陣中該類列的總和）的比率。

用戶精度（user’s accuracy）是指正確分到某類的像元總數（對角線值）與分類器將整個影像的像元分為該類像元總數（混淆矩陣中該類行的總和）比率。

4.3 對比方法

為了驗證本文方法的有效性，給出以下幾個對比方法。為了公平，對比方法均包含KNN正則化處理過程。

塊k近鄰（blockKNN，BKNN）圖：K近鄰圖（KNN圖）每個頂點只與自己的k個近鄰連接[26]。這里近鄰數k設置為4。采用高斯核計算權重。

塊局部線性嵌入（block locally linear embedding，BLLE）圖：從局部線性嵌入方法通過最小化l2重構樣本。將下式求得的局部線性重構系數作為樣本連接權值[27]。

4.4 分類結果

由于分類結果對訓練樣本的選取有很大的依賴性，因此選取的訓練樣本要具有典型性和代表性。本實驗借助ENVI處理軟件選取各個地物特征明顯的區域作為訓練樣本。從數值結果和視覺效果上進行了比較。一般來說，地物像元的一致性只能通過采樣來檢查，因為檢查所有的像元是不現實的。驗證樣本可以通過野外實地調查或在高分辨率圖像上通過目釋解譯的方法獲取各個分類的地表真實感興趣區[28-29]。由于實驗影像分辨率為2 m，結合武漢市區的文字圖片資料，在影像上選取6 270個驗證樣本和4 762個訓練樣本，如表1所示。表1給出每類分類精度、總體分類精度等評價指標。

Table 1 Classification results for different feature graphs表1 不同特征圖的分類結果

表2～表5給出了混淆矩陣、制圖精度和用戶精度等。從表1～表5可以看出，BLR分類方法無論是在總體分類精度上還是Kappa系數上都比BSR、BKNN、BLLE方法有優勢。Kappa系數比其他方法提高了約6個百分點。對比表1中的三種方法分類結果，可以看出基于BLR的影像分類結果比基于BSR的總體準確率提高了4.83%，Kappa系數提高了5.79%，說明基于BLR的分類能夠較好提高高分辨率遙感影像的分類精度。雖然基于BSR方法的總體精度和Kappa系數較高，但從表3可以看出其對耕地、裸地和草地混分較為嚴重，耕地的制圖精度只有70.35%。這是由于耕地和草地光譜較為相似，二者光譜特征在各個波段波動范圍較大，在同一波段上反射率具有一定的重疊。稀疏表示圖缺乏對數據的全局約束，每個樣本求解過程相互獨立。基于分塊低秩圖方法分類結果的總體精度比BKNN和BLLE方法的總體精度分別提高了5.57%和4.94%，Kappa系數分別提高了6.67%和5.92%。這是因為K近鄰和LLE使用固定的全局參數，利用歐幾里德距離來確定圖的權重，其圖矩陣結構對噪聲敏感。

Table 2 Classification results of BLR graph method表2 BLR正則圖方法分類結果

Table 3 Classification results of BSR graph method表3 BSR正則圖方法分類結果

Table 4 Classification results of BKNN graph method表4 BKNN正則圖方法分類結果

基于塊低秩圖的分類方法中耕地的制圖精度達到90.43%，比稀疏表示方法提高了20.08%。這是因為LR圖解決了對比方法的缺點，它根據低秩約束獲得數據的全局信息，并且對噪聲表現出魯棒性。同時使用K近鄰來處理整個低秩圖，保留了數據的局部信息。因此，塊低秩圖方法顯著提高了分類性能，對于大規模遙感影像是非常好的特征提取方法，能夠在一定程度上改善耕地和草地的分類精度。然而，對于耕地、裸地和草地，由于在城市中樣本數比較少，容易受到其他相對較大規模地物的混淆，導致基于塊稀疏表示方法的分類精度沒有較大的提高。由表2看出，基于BLR方法的各類制圖精度都在90%以上，能夠滿足一般遙感影像制圖精度的需求。因此，基于分塊低秩的方法用于大規模遙感影像分類有一定的優勢。

Table 5 Classification results of BLLE graph method表5 BLLE正則圖方法分類結果

從視覺效果上對實驗結果進行比較，如圖2所示，基于分塊低秩的分類方法和其他對比方法相比，有效地提高了影像的分類精度。例如，圖中左上角水體區域以及右下角的裸地等誤分為其他地物的問題得到了有效的改善，在一定程度上提高了影像分類精度，在影像分類上具有一定的可行性。遙感影像的分類受到眾多因素制約，許多問題還有待進一步研究。實驗中缺乏足夠的實地調查信息和受目視判讀的影響，樣本地物類別判斷上可能存在一定的偏差。實驗只選擇了光譜信息作為特征數據，在未來可以考慮將圖像的紋理等信息引入到圖像分類中。

5 結論

Fig.2 Land cover classification of Wuhan municipal district圖2 武漢市轄區土地利用覆蓋分類結果

本文把塊低秩表示方法應用到高分遙感影像半監督分類上。為了降低計算復雜度，整個圖像被分塊并分別實現每個塊上的低秩表示。因此，通過低秩表示保留了影像的全局結構。另外，K近鄰算法保留了影像的局部幾何結構。因此，基于塊低秩正則圖的半監督判別分析特征提取方法可以有效提高圖像分類的性能。對武漢市WorldView-2影像進行實驗表明，本文方法是有效和魯棒的，可以實現大規模遙感影像的精確分類。