曾 夢(mèng)，寧 彬，蔡之華，谷 瓊

（1.湖北文理學(xué)院計(jì)算機(jī)工程學(xué)院，湖北襄陽(yáng)441053；2.中國(guó)地質(zhì)大學(xué)（武漢）計(jì)算機(jī)學(xué)院，武漢430074）

0 引言

高光譜圖像（HyperSpectral Image，HSI）由上百個(gè)光譜波段組成，包含豐富的空間信息和光譜信息，通常利用分類(lèi)技術(shù)實(shí)現(xiàn)具有細(xì)微光譜差異地物的區(qū)分。目前，高光譜圖像分類(lèi)已經(jīng)被廣泛用于軍事目標(biāo)識(shí)別［1］、海洋環(huán)境監(jiān)測(cè)［2］和城市管理［3］等多個(gè)方面。然而，由于高光譜圖像光譜波段之間相關(guān)性強(qiáng)且具有較高的冗余度，導(dǎo)致分類(lèi)復(fù)雜性高、“維度災(zāi)難”等［4］問(wèn)題的發(fā)生。高光譜波段選擇旨在選取原始波段集中較為重要的波段子集來(lái)降低數(shù)據(jù)的維度，從而降低原始數(shù)據(jù)的計(jì)算量并保證分類(lèi)具有較高精度［5］。

近些年，國(guó)內(nèi)外學(xué)者提出了很多波段選擇方法，大致可以分為三類(lèi)：基于排序策略、基于搜索策略和基于聚類(lèi)策略?；谂判虿呗缘姆椒ㄍㄟ^(guò)評(píng)估每個(gè)波段的重要性來(lái)分配權(quán)重進(jìn)行排序，例如最大方差主成分分析（Maximum Variance Principal Component Analysis，MVPCA）［6］、拉 普 拉 斯 評(píng) 分（Laplacian score，Lap-score）［7］以及基于稀疏表示的波段選擇（Sparsity-based Band Selection，SpaBS）［8］等。基于搜索策略的方法使用啟發(fā)式搜索算法優(yōu)化給定的度量，如演化多目標(biāo)優(yōu)化（evolutionary multiobjective optimization）算法［9］等，但是基于搜索的方法通常耗時(shí)較多?；诰垲?lèi)的方法假設(shè)波段根據(jù)它們的相關(guān)性可分離，通過(guò)聚類(lèi)結(jié)果進(jìn)行波段選擇，例如改進(jìn)的稀疏子空間聚類(lèi)算法（Improved Sparse Subspace Clustering，ISSC）［10］、稀 疏 非 負(fù) 矩 陣 分 解（Sparse Nonnegative Matrix Factorization，SNMF）［11］等。由于基于聚類(lèi)的方法充分考慮了波段之間的相關(guān)性，能夠得到較好的波段選擇結(jié)果，因此得到了越來(lái)越多的關(guān)注。

上述基于聚類(lèi)的算法雖然能得到較好的結(jié)果，但是很多基于聚類(lèi)的方法忽略了光譜波段之間的相關(guān)性。子空間聚類(lèi)算法是廣泛應(yīng)用的聚類(lèi)算法，但很多子空間聚類(lèi)算法僅對(duì)于線性子空間進(jìn)行聚類(lèi)，這不適用于具有非線性關(guān)系的高光譜圖像的波段選擇。最近提出的一種深度子空間聚類(lèi)（Deep Subspace Clustering，DSC）［12］網(wǎng)絡(luò)將深度卷積自編碼器和子空間聚類(lèi)相結(jié)合，使用該算法進(jìn)行波段選擇［13］很好地解決了上述問(wèn)題并取得了較好的效果；然而，DSC算法仍然將自表達(dá)作為監(jiān)督，對(duì)于某些內(nèi)在子空間不獨(dú)立或具有明顯交集的樣本可能表現(xiàn)不佳。

本文利用深度對(duì)抗子空間聚類(lèi)（Deep Adversarial Subspace Clustering，DASC）［14］網(wǎng)絡(luò)進(jìn)行高光譜波段選擇。DASC網(wǎng)絡(luò)是一種新的無(wú)監(jiān)督深度子空間聚類(lèi)模型，該網(wǎng)絡(luò)在DSC 網(wǎng)絡(luò)的基礎(chǔ)上加以改進(jìn)，將DSC 和生成對(duì)抗網(wǎng)絡(luò)（Generative Adversarial Network，GAN）［15］相結(jié)合，引入對(duì)抗學(xué)習(xí)來(lái)監(jiān)督自編碼器的樣本表示和子空間聚類(lèi)，使得子空間聚類(lèi)具有更好的自表達(dá)性能。此外，引入拉普拉斯正則化來(lái)考慮反映圖像幾何信息的局部流形結(jié)構(gòu)，使網(wǎng)絡(luò)性能更優(yōu)。

1 相關(guān)工作

1.1 子空間聚類(lèi)

對(duì)于一個(gè)原始數(shù)據(jù)集X=[x1，x2，…，xN]∈?d×N，假設(shè)該數(shù)據(jù)集屬于n 個(gè)線性子空間S=S1，S2，…，Sn，線性子空間的維度分別為d1，d2，…，dn，有如下公式：

s.t.X=XC，diag(C)=0

其中：C ∈?N×N表示自表達(dá)系數(shù)矩陣，diag(C)=0 為了避免平凡解，‖ ? ‖p可以表示為l1或l2正則化。對(duì)于稀疏子空間聚類(lèi)（Sparse Subspace Clustering，SSC）［16］，通常有如下公式：

s.t.di ag( C )=0

其中，λ 表示平衡系數(shù)。使用系數(shù)矩陣C 構(gòu)造相似矩陣M=|C|+|C|T，通常情況下，相似矩陣M 為塊對(duì)角矩陣。將矩陣M用于譜聚類(lèi)進(jìn)行聚類(lèi)，得到最終的聚類(lèi)結(jié)果。

1.2 卷積自編碼器

卷積自編碼器（Convolutional AutoEncoder，CAE）［17］由一個(gè)對(duì)稱(chēng)的編碼器和解碼器組成，編碼器和解碼器都由卷積神經(jīng)網(wǎng)絡(luò)構(gòu)成，用φ=E(X；αe)表示編碼器，其中：X 表示輸入，αe表示參數(shù)，φ 表示編碼器的輸出。對(duì)于解碼器，采用和編碼器相同的結(jié)構(gòu)，不同的是解碼器操作用于重建輸入的X 表示為最后，使用均方誤差作為損失函數(shù)進(jìn)行訓(xùn)練，損失函數(shù)如下：

1.3 深度子空間聚類(lèi)

深度子空間聚類(lèi)（Deep Subspace Clustering，DSC）網(wǎng)絡(luò)將自表達(dá)層嵌入在卷積自編碼器中，輸入數(shù)據(jù)X 經(jīng)過(guò)編碼器得到特征表示，U然后經(jīng)過(guò)自表達(dá)使用UC進(jìn)行全連接，再將UC經(jīng)過(guò)解碼器得到重構(gòu)的DSC的損失函數(shù)L如下：

s.t.diag(C)=0

2 深度對(duì)抗子空間聚類(lèi)實(shí)現(xiàn)波段選擇

2.1 深度對(duì)抗子空間聚類(lèi)

深度對(duì)抗子空間聚類(lèi)（DASC）網(wǎng)絡(luò)由子空間生成器和判別器兩部分組成：在生成器中，使用一個(gè)卷積自編碼器來(lái)學(xué)習(xí)樣本表示，自表達(dá)層嵌入在編碼器和解碼器中用來(lái)得到樣本的相似矩陣并使用Ncut 方法進(jìn)行聚類(lèi)，然后以聚類(lèi)生成的簇為條件生成“真樣本”和“假樣本”；對(duì)于判別器，接收生成的真假樣本，并且區(qū)分出真樣本和假樣本。通過(guò)DASC 網(wǎng)絡(luò)訓(xùn)練，最終得到系數(shù)矩陣C，通過(guò)譜聚類(lèi)得到聚類(lèi)結(jié)果。

2.1.1 生成器

對(duì)于真樣本的獲取，通過(guò)自表達(dá)層學(xué)習(xí)的相似矩陣運(yùn)用Ncut 算 法 得 到k 個(gè) 簇A={A1，A2，…，Ak}及 特 征 表 示Z=[Z1，Z2，…，Zk]，然后計(jì)算每一類(lèi)Ai中每個(gè)樣本到相對(duì)應(yīng)的子空間Si的投影殘差，投影殘差的求解在判別器中執(zhí)行，然后選擇具有較小殘差的mi個(gè)數(shù)據(jù)作為真樣本。投影殘差的公式如下：

其中，Vi表示投影矩陣。由此生成器的對(duì)抗性損失為：

對(duì)于假樣本的獲取，在每個(gè)簇Ai中，采樣層從估計(jì)的子空間Si中進(jìn)行隨機(jī)采樣，假設(shè)簇Ai具有mi個(gè)樣本，采樣層從(0，1]均勻分布中隨機(jī)采樣mi個(gè)隨機(jī)向量，然后將αt的每一項(xiàng)和簇Ai中對(duì)應(yīng)樣本相乘并求和，生成mi個(gè)假數(shù)Zˉt=，其中θtj表示θt的第j項(xiàng)。

為了使生成的假數(shù)據(jù)更接近判別器學(xué)習(xí)的子空間，使聚類(lèi)具有更好的效果，將對(duì)抗損失加入到生成器的損失函數(shù)中，公式如下：

s.t.diag(C)=0

其中，λ 表示平衡系數(shù)。式（7）右邊的四項(xiàng)分別表示對(duì)抗性損失、重構(gòu)損失、自表達(dá)損失和自表達(dá)系數(shù)損失。

2.1.2 判別器

在判別器中，通過(guò)投影殘差Lr來(lái)區(qū)分真樣本和假樣本，通常來(lái)說(shuō)，若簇不準(zhǔn)確則真樣本將比假樣本更接近子空間。對(duì)于簇Ai的訓(xùn)練，使用映射殘差定義一個(gè)樣本屬于子空間的概率，損失函數(shù)如下：

其中：[?]+=max(0，?)為附加的邊緣損失，ε 為正邊緣參數(shù)。對(duì)于k個(gè)簇來(lái)說(shuō)，判別器的損失函數(shù)為：

2.2 拉普拉斯正則化DASC

為了更好地利用圖像的局部流形結(jié)構(gòu)信息，在生成器的損失函數(shù)中加入拉普拉斯正則項(xiàng)［18］。將n 個(gè)樣本作為n 個(gè)頂點(diǎn)X=[x1，x2，…，xN]用邊連接起來(lái)構(gòu)成圖，其中xi∈?d。定義W ∈?n×n作為權(quán)重矩陣來(lái)記錄每條邊的權(quán)重，則頂點(diǎn)i 和j之間的權(quán)重可以表示為：

其中：Nk表示具有n 個(gè)頂點(diǎn)的K 最近鄰（K Nearest Neighbors，KNN）的構(gòu)造，σ是熱核參數(shù)。

對(duì)于非線性流形結(jié)構(gòu)，可以通過(guò)線性子空間局部近似，因此定義二次能量函數(shù)Q來(lái)計(jì)算兩點(diǎn)之間的相似性：

其中q(x)=xTw+b為線性表示。

綜上所述，最終得到生成器的損失函數(shù)為：

s.t.diag(C)=0

其中：tr(?)表示矩陣的跡，γ是平衡系數(shù)。

2.3 拉普拉斯正則化DASC實(shí)現(xiàn)波段選擇的流程

利用深度對(duì)抗子空間聚類(lèi)得到最終聚類(lèi)結(jié)果，根據(jù)聚類(lèi)結(jié)果求出每一類(lèi)中波段的平均值作為聚類(lèi)中心，然后計(jì)算每個(gè)波段到該類(lèi)聚類(lèi)中心的距離，選出每個(gè)類(lèi)中最接近聚類(lèi)中心的波段得到波段子集，即為波段選擇的結(jié)果，具體流程如下：

1）輸入三維高光譜數(shù)據(jù)集X，設(shè)置波段選擇數(shù)量k。

2）構(gòu)建DASC網(wǎng)絡(luò)，初始化網(wǎng)絡(luò)參數(shù)。

3）使用深度子空間聚類(lèi)（DSC）網(wǎng)絡(luò)對(duì)數(shù)據(jù)進(jìn)行訓(xùn)練得到自表達(dá)系數(shù)矩陣和特征表示。

4）使用自表達(dá)層得到的數(shù)據(jù)，采用Ncut 方法進(jìn)行聚類(lèi)，得到“真樣本”，并根據(jù)聚類(lèi)個(gè)數(shù)生成m個(gè)假樣本；并判斷生成器訓(xùn)練是否大于5 次，如果是則進(jìn)行下一步，否則返回到步驟3）。

5）根據(jù)生成器生成的真假樣本，對(duì)判別器進(jìn)行訓(xùn)練。判斷如果總迭代次數(shù)大于300，則進(jìn)行下一步，否則返回到步驟3）。

6）將DASC 網(wǎng)絡(luò)訓(xùn)練得到的系數(shù)矩陣C 用于譜聚類(lèi)得到聚類(lèi)結(jié)果，計(jì)算得到每個(gè)類(lèi)的聚類(lèi)中心，并計(jì)算每個(gè)類(lèi)中離聚類(lèi)中心最近的波段作為代表波段，得到波段子集，即為波段選擇結(jié)果。

3 實(shí)驗(yàn)結(jié)果與分析

3.1 實(shí)驗(yàn)數(shù)據(jù)和預(yù)處理

Indian Pines 數(shù)據(jù)集使用AVRIS 傳感器在美國(guó)印第安納州西北部的一塊小區(qū)域收集，包含145 像素×145 像素，預(yù)處理后波段數(shù)為200，主要地物為16類(lèi)。為方便實(shí)驗(yàn)，選取像素區(qū)間在[50 ～120，50 ～120]的數(shù)據(jù)子集進(jìn)行實(shí)驗(yàn)，共包含8類(lèi)，如圖1所示。

圖1 Indian Pines 數(shù)據(jù)集Fig.1 Indian Pines dataset

Pavia Center數(shù)據(jù)集使用ROSIS傳感器在意大利北部帕維亞大學(xué)收集，包含1 096 像素×715 像素，預(yù)處理后波段數(shù)為102 個(gè)，主要地物為9 類(lèi)。為方便實(shí)驗(yàn)，選取像素區(qū)間在[100 ～200，200 ～400]的數(shù)據(jù)子集進(jìn)行實(shí)驗(yàn)，共包含6 類(lèi)，如圖2所示。

圖2 Pavia Center 數(shù)據(jù)集Fig.2 Pavia Center dataset

詳細(xì)的數(shù)據(jù)總結(jié)在表1 中給出。所有方法均使用python 3.5進(jìn)行實(shí)驗(yàn)評(píng)價(jià)，DASC算法采用tensorflow框架。

表1 Indian Pines和Pavia Center數(shù)據(jù)集參數(shù)設(shè)置Tab.1 Parameters setting of Indian Pines and Pavia Center datasets

3.2 實(shí)驗(yàn)結(jié)果對(duì)比

使用Indian Pines 數(shù)據(jù)集和Pavia Center 數(shù)據(jù)集設(shè)計(jì)實(shí)驗(yàn)驗(yàn)證DASC 算法進(jìn)行波段選擇的效果。實(shí)驗(yàn)使用支持向量機(jī)（Support Vector Machine，SVM）［19］分類(lèi)器進(jìn)行波段選擇后的數(shù)據(jù)分類(lèi)，采用總體精度（Overall Accuracy，OA）、卡帕系數(shù)（Kappa Coefficient）以及平均精度（Average Accuracy，AA）來(lái)評(píng)價(jià)分類(lèi)效果。在實(shí)驗(yàn)中，設(shè)置波段數(shù)為5～40，間隔為5個(gè)波段。對(duì)于Indian Pines 數(shù)據(jù)集，設(shè)置式（14）中的λ=0.4，α=15，β=1，γ=10；對(duì)于Pavia Center 數(shù)據(jù)集，設(shè)置式（14）中的λ=0.2，α=20，β=1，γ=15。對(duì)于自編碼器的通道數(shù)設(shè)置為10-50-100-C-100-50-10，卷積核大小為3×3×3。將DASC 和Lap-score、SpaBS、ISSC、SNMF、MVPCA、DSC 等6種方法進(jìn)行比較，得到的實(shí)驗(yàn)結(jié)果如圖3、圖4 所示。從圖中可以看出，采用DASC方法進(jìn)行波段選擇的效果明顯比其他方法更優(yōu)。DSC和ISSC方法的效果比SpaBS、SNMF、MVPCA方法好，特別是DSC 方法表現(xiàn)出較好的實(shí)驗(yàn)效果，這表明將深度學(xué)習(xí)和子空間聚類(lèi)算法相結(jié)合，對(duì)高光譜圖像波段選擇具有很好的應(yīng)用價(jià)值，而使用DASC 方法，解決了DSC 方法的不足，在OA、AA和Kappa系數(shù)的比較中都得到了比DSC更好的結(jié)果。

圖3 Indian Pines數(shù)據(jù)集實(shí)驗(yàn)結(jié)果比較Fig.3 Experimental result comparison on Indian Pines dataset

圖4 Pavia Center數(shù)據(jù)集實(shí)驗(yàn)結(jié)果比較Fig.4 Experimental result comparison on Pavia Center dataset

3.3 波段選擇分析

為了更好地分析波段選擇的結(jié)果，使用Indian Pines 數(shù)據(jù)集，對(duì)波段選擇數(shù)為20 時(shí)的波段號(hào)進(jìn)行分析，分別得到不同方法選擇20 個(gè)波段時(shí)的波段號(hào)。為了直觀表示波段選擇的結(jié)果，畫(huà)出Indian Pines數(shù)據(jù)集200個(gè)波段的信息熵，信息熵越大表明該波段的信息量越多，波段更具有代表性。將選擇的波段號(hào)在圖中畫(huà)出，如圖5、圖6所示。

圖5 選擇波段數(shù)為20時(shí)不同方法的波段分布比較Fig.5 Comparison of band distribution of different methods when the number of bands is 20

圖6 選擇波段數(shù)為20時(shí)光譜波段信息熵Fig.6 Spectral band information entropy when the number of bands is 20

從圖6 中可以看出，信息熵的值比較平均，僅在某些地方呈下降趨勢(shì)，因此選擇波段分布較為均勻的方法更好。對(duì)應(yīng)到圖5 中，對(duì)于不同方法得到的波段選擇結(jié)果有很大的差別，其中Lap-score 和MVPCA 方法得到的波段號(hào)較為集中，而且僅出現(xiàn)在某一區(qū)間內(nèi)，波段之間具有很強(qiáng)的相關(guān)性；SpaBS、ISSC 和SNMF 方法選擇的波段在前面較為集中，而后面比較分散，分布不均勻；DSC 和DASC 方法得到的波段號(hào)分布均勻，且選出的相鄰波段較少，具有很好的代表性，因此能得到更好的分類(lèi)效果。通過(guò)3.2 節(jié)的實(shí)驗(yàn)結(jié)果比較，可以證實(shí)DASC方法選出的波段比其他方法更優(yōu)。

3.4 計(jì)算時(shí)間對(duì)比

該實(shí)驗(yàn)將DASC 方法和其他6 種方法在波段選擇數(shù)量均為20 時(shí)的計(jì)算效率進(jìn)行比較。不同方法在Indian Pines 和Pavia Center 數(shù)據(jù)集上的運(yùn)行時(shí)間如表2所示。實(shí)驗(yàn)方法在配置為Inter Xeon E5-2620 CPU 2.10 GHz 和32 GB 內(nèi)存的服務(wù)器上進(jìn)行。

表2 波段數(shù)為20時(shí)不同方法的運(yùn)行時(shí)間對(duì)比 單位：sTab.2 Running time comparison of different methods when the number of bands is 20 unit：s

由表2 可以看出，ISSC 和MVPCA 方法所用時(shí)間較少，SpaBS 和SNMF 方法計(jì)算效率較低，DSC 和DASC 的計(jì)算效率適中。同時(shí)結(jié)合分類(lèi)精度和實(shí)驗(yàn)效果來(lái)看，DASC的綜合表現(xiàn)最好，具有很好的應(yīng)用價(jià)值。

4 結(jié)語(yǔ)

本文提出使用DASC 方法來(lái)實(shí)現(xiàn)高光譜波段選擇。通過(guò)在編碼器和解碼器中引入自表達(dá)層來(lái)模仿傳統(tǒng)子空間聚類(lèi)的“自表達(dá)”屬性，充分運(yùn)用光譜信息和非線性特征轉(zhuǎn)換得到波段之間的相互關(guān)系，并且將深度子空間聚類(lèi)和生成對(duì)抗網(wǎng)絡(luò)相結(jié)合，引入對(duì)抗學(xué)習(xí)來(lái)監(jiān)督自編碼器的樣本表示和子空間聚類(lèi)，使得子空間聚類(lèi)具有更好的自表達(dá)性能。同時(shí)，加入拉普拉斯正則化來(lái)考慮反映圖像幾何信息的局部流形結(jié)構(gòu)，使模型更加魯棒。實(shí)驗(yàn)結(jié)果表明，DASC方法在波段選擇上具有更好的效果，綜合來(lái)看具有較高的應(yīng)用價(jià)值。但是，由于該方法無(wú)法使用批梯度下降方法進(jìn)行訓(xùn)練，因此在計(jì)算效率方面存在不足，這將是后續(xù)研究和改進(jìn)的重要方向。