基于頂點成分分析的高光譜圖像端元提取算法?

方凌江 粘永健 雷樹濤 倪志揚

（濟南軍區聯勤部后勤信息中心 濟南 250022）

1 引言

高光譜遙感是指光譜分辨率高達10－2λ數量級的光譜成像技術，在地質勘探、軍事偵察等領域得到了廣泛的應用［1］。由于傳感器空間分辨率的限制，以及地物在地表空間分布的復雜性，使得所拍攝的圖像中存在大量的混合像元。端元提取就是從存在混合像元的高光譜圖像中提取出場景的基本組成成分，它是高光譜解混的重要步驟。當物質之間的混合方式屬于非致密混合時，光子在物質內部之間的相互作用可以忽略，可用線性混合模型來描述混合像元的形成機理［2］。近十年來，學者們在線性混合模型的基礎上提出了多種端元提取算法，例如純像元指數（Pure Pixel Index，PPI）［3］、頂點成分分析算法（Vertex Component Analysis，VCA）［4］、N－FINDR［5～6］，迭 代 誤 差 分 析 （Iterative Error Analysis，IEA）［7］，自動形態學端元提取方法（Automated Morphological Endmember Extrac－tion，AMEE）［8］、最小體積單形體分析（Minimum Volumn Simplex Analysis，MVSA）［9］，其中，VCA算法是一種端元逐次提取算法，前面提取的端元不依賴后面提取的端元，有效避免了在端元數目估計不正確的情況下，所提取的端元之間相互干擾的問題，而且該算法提取速度快，計算復雜度低。但是VCA算法是將高光譜圖像中的極值像元作為端元，易受異常像元的干擾，而且它以高光譜圖像中的某個像元作為端元，受噪聲影響也較大。

本文在VCA算法的基礎上，利用地物在地表空間成片分布這一特點，采取空間鄰域迭代計算來提取空間純區域，并將這些區域內的所有像元的均值作為待提取的端元，該算法可以有效地消除異常像元的干擾，并降低噪聲對端元的影響，快速有效地提取出端元。

2 算法描述

2.1 線性混合模型

在線性光譜混合模型中，每個像元的光譜響應由各個端元的光譜響應按照其對應的豐度線性混合而成。用x∈RL×1表示波段數為L的高光譜圖像中某像元的光譜響應，x可以表示為

其中M表示端元數量，N 表示像元數量，E＝［e1，…，eM］∈RL×M表示構成整幅高光譜圖像的端元矩陣，ei表示第i個端元，各個端元仿射獨立。s［n］＝［s1［n］，…，sM［n］］T∈RM×1表示第n個像元的豐度向量。w［n］∈RL×1為誤差向量，表示該像元被噪聲污染的程度。假定噪聲為零均值加性高斯噪聲，則w［n］～N（0，diag（σ2j）），其中diag（σ2j）為L×L的對角矩陣，σ2j（j＝1，2，…，L）為各波段的噪聲方差。根據豐度的物理含義，豐度還要滿足非負約束（Abundance Non－negative Constraint，ANC）及和為 一 約 束 （Abundance Sum－to－one Constraint，ASC），可以表示為

在沒有噪聲的干擾下，所有滿足式（1）～（3）的像元構成一個以各端元為頂點的凸面單形體，VCA正是利用正交投影的方法尋找單形體的頂點來提取端元，詳細步驟可以參考文獻［4］。此時，所有的純像元也都位于單形體的頂點，純像元就是端元，但是在有噪聲的干擾下，單形體體積膨脹，純像元偏離端元的位置，位于端元的周圍，此時端元不能由單個純像元表示。

用單形體來近似擬合有噪聲干擾的高光譜數據，靠近單形體的頂點的像元純度較高，PPI算法正是利用像元純度指數來提取端元的，可以認為，由VCA算法所得到的像元是純像元，但是由于VCA所提取的純像元只有一個，受噪聲影響較大。對某種典型地物而言，通常它在高光譜圖像中分布有一定數量的純像元，如果可以找到更多該地物的純像元，并以這些純像元的均值作為該地物端元提取結果，根據統計學知識，此時的提取精度要比以單個像元作為提取結果的精度要高。

2.2 空譜聯合約束下的端元優化

典型地物在空間上分布連續，純像元的鄰近像元仍然為該類地物的純像元的可能性較大?？兆V聯合約束要求VCA提取的純像元位于光譜相似、且空間連續的一致區域。下面詳細描述空譜聯合約束的實現過程：以VCA提取的單個純像元為種子，采用鄰近像元搜索，逐步擴大純像元的數量，其過程如圖1所示。

圖1 鄰近像元搜索過程

記v為VCA在第i步所提取的純像元，搜索鄰近像元記為FindNV＿Pixels，實現的流程如表1所示。

表1 鄰像元搜索的流程圖

在表2中，V表示純像元集合，Nb表示取鄰域，本文中采用8鄰域搜索，mean為取均值，dist為衡量兩個像元相似度的準則，τ為門限。由FindNV＿Pixels得到的V為純像元區域，記V＝｛v1，…，vK｝，K 為V 中純像元的數量，V 具有兩個性質：

1）?vk1，vk2∈V（k1≠k2），都有dist（vk1，vk2）≤τ成立，體現了對純像元區域具有光譜相似的約束。

2）?vk1，vk2∈V（k1≠k2），存在一條路徑U＝｛vk，vk∈V｝，使得vk1、vk2連通，體現了對地物在空間上具有連續性的約束。

然而，還有兩個尚待解決的問題：一是dist的選擇；二是τ的選擇。下面根據統計檢測理論中的Neyman－Pearson（NP）準則來解決這兩個問題。vi、vj這兩個像元的差ε為

其中μ是由vi、vj兩個像元成分不同所引起的非噪聲差值。定義隨機變量r

由式（4）與式（5）可知，如果xj為該類地物的純像元，則ε～N（0，2diag（）），r為自由度為L 的chi平方分布；否則，ε～N（μ，2diag（）），r為自由度為L，非中心參量為‖μ‖的非中心chi平方分布。定義兩個假設：

為提高純像元區域的純度，則應使vj誤判為純像元的概率最小，也就是使P（H0｜H1）最小。預設虛警概率PFA（PFA＝P（H1｜H0）），根據NP準則得到：

其中τ為判決門限，可通過χ2L分布及PFA得到。與SPA、AMEE僅利用光譜信息來選取相似像元不同，本文算法利用了空間上下文信息。

在高光譜圖像中，通常會有大量的異常點存在，而降維及去噪等預處理措施并不能完全去除異常點，所以在降維及去噪后的高光譜圖像中仍然存在一定數量的異常點。由于異常點大多都屬于數據點云的極值點，VCA會受異常點的干擾，甚至在提取結果中直接將異常像元作為端元提取的結果，從而導致算法錯誤。而異常點通常都是孤立點，具有數量少、在空間上不連續以及與周圍點的差異比較明顯等特征，可以利用這些特征來消除端元提取過程中異常點的影響。設定純像元區域大小門限為η，當K≤η時就判定為異常像元。綜上，算法的整體流程為

表2 算法的整體流程圖

3 實驗結果

將本文提出的算法記為SS－VCA，利用 Hysime［10］來估計噪聲，端元數目已知。與 AMEE、SSEE、SPEE等算法一樣，本文算法也利用了空間上下文信息，但是上述算法的性能依賴于多個參數的設置，而且參數的設置與場景的空間復雜程度有關，很難用統一的標準來比較本文算法與上述三種算法的優劣，而且本文的重點在于強調對VCA的改進作用，所以在實驗中并不將上述三種算法作為比較對象。比較準則為各個端元與真實端元之間的角距離，用ei、＾ei分別表示參考端元與算法提取的端元，其定義為

3.1 仿真數據實驗

從美國地質調查局（USGS）光譜庫中任意選擇五種地物的光譜作為仿真數據的端元，該光譜庫中地物光譜的波長范圍是0.38μm～2.5μm，光譜分辨率為10nm，共有224個波段。豐度的生成方法是，先利用 MRF（Markov Random Field）生成空間分辨率為400×400的含有三種類別地物的分類圖，再由線性混合模型模仿由于傳感器分辨率不足而引起的地物混合現象，生成50×50的高光譜圖像。下面用兩組實驗來測試算法的性能。在所有的試驗中，虛警概率PFA設置為10－2，異常像元大小η為1。

圖2 不同信噪比條件下的算法性能比較

實驗1（測試算法的抗噪能力）。在仿真數據中加入加性零均值高斯白噪聲，并逐漸改變信噪比（SNR）的大小，由10dB逐漸增加到30dB，步長為5dB。圖2給出了不同算法在不同信噪比條件下的平均角距離，記為A＿SAD。從圖2可以得出，與VCA相比，SS－VCA具有更強的抗噪能力。

實驗2（測試異常像元對算法的影響）。在USGS光譜庫中另外選擇一種地物，并在空間中隨機選擇五個不相鄰的位置，將原來數據中的像元的光譜替換為當前地物的光譜，以模仿高光譜圖像中的異常像元。此時SNR設為10dB。由圖3可以直觀看出，VCA算法提取的端元中包含了異常像元，而SS－VCA則正常提取端元。

圖3 異常像元存在條件下的算法性能比較

3.2 真實數據實驗

實驗中采用了美國加州印度松樹測試地的高光譜數據，該數據由機載可見光及紅外成像光譜儀（Airborne Visible／Infrared Imaging Spectrometer，AVIRIS）拍攝而成，該數據波長范圍是0.4μm～2.5μm，光譜分辨率為10nm，共有224個波段，空間分辨率為3.7m，截取大小為86×83子塊作為實驗對象，該塊主要分布椰菜、玉米、雜草以及不同類別的生菜等自然地物，雖然該高光譜圖像空間分辨率較高，但由于地物分布復雜，仍存在相當數量的混合像元。由于信噪比太低或者水吸收，36個波段被移除（波段號為1～4，104～112，150～167，220～224），僅保留188個波段用于數據處理。

為評價算法的性能，在實驗中根據地表真實標記圖手工從高光譜數據中選擇各類地物的多個純像元，并取各類純像元的均值作為參考光譜。考慮到部分地物內部存在較大的光譜變異性，難以用一個端元來描述，在實驗中取端元數量為10。在提取端元之后，由全約束解混FCLS［11］估計的豐度，由SS－VCA提取的端元所得到的豐度估計如圖4所示。

比較真實標記圖與估計豐度得到端元對應的地物，其中生菜（5wk）對應（d）、（f）、（h）三個端元，說明該地物具有較大的內部變異性，而雜草對應（g）、（j）兩個端元，分別表示開始衰老的雜草及綠色雜草。SS－VCA、VCA的端元提取結果分別與參考光譜的比較如表3所示。

圖4 Salinas－A數據及SS－VCA豐度所提取端元對應的豐度

表3 Salinas－A數據的角距離比較（粗體數字表示最好的結果）

在克服異常像元的影響方面，在VCA提取的端元中，有一個端元為異常像元，而SS－VCA沒有。在提取結果的精度方面，由表3可以得出SS－VCA所提取的端元具有更小的平均角距離，精度更高。

4 結語

本文算法幾乎不受異常像元的影響，并具有較強的抗噪能力。但它只能找出其中一塊屬于某類地物的純區域，不能將空間不相鄰但屬于該類地物的所有純區域同時找出來，抗噪性能的提高還有一定的空間，這是下一步要研究的方向。

［1］Keshava N，Mustard J F.Spectral Unmixing［J］.IEEE Signal Processing Magazine，2002，19（1）：44－57.

［2］Adams J B，Gillespie A R.Spectral－Mixture Analysis［M］.Remote Sensing of Landscapes with Spectral Images：A physical Modeling Approach，Cambridge University，2004：126－167.

［3］Boardman J W，Kruse F A，Green R O.Mapping Target Signatures via Partial Unmixing of AVIRIS Data［C］.1995.

［4］Nascimento J M P，Bioucas－Dias J M.Vertex Component Analysis A Fast Algorithm to Unmix Hyperspectral Data［J］.IEEE Transaction on Geoscience and Remote Sensing，2004，43（8）：898－910.

［5］Winter M E.N－FINDR：An Algorithm for Fast Autonomous Spectral Endmember Determination in Hyperspectral Data［C］／／Proceedings of SPIE，Denver，1999，3753：266－275.

［6］Plaza A，Chang C I.An Improved N－FINDR Algorithm in Implementation.Algorithms and Technologies for Multispectral，Hyperspectral and Ultraspectral Imagery［C］／／Proceedings of SPIE，Bellingham，WA，2005，5806：298－306.

［7］Neville R A，Staenz K，Szeredi T，et al.Automatic Endmember Extraction from Hyperspectral Data for Mineral Exploration［C］／／Proceedings of International Conference on Airborne Remote Sensing，Ottawa，Canada，1999：21－24.

［8］Plaza A，Martinez P，Perez R，et al.Spatial／Spectral Endmember Extraction by Multidimensional Morphological Operations［J］.IEEE Transaction on Geoscience and Remote Sensing，2002，40（9）：2025－2041.

［9］Li J，Bioucas－Dias J E M.Minimum Volume Simplex Analysis：a Fast Algorithm to Unmix Hyperspectral Data［C］／／Proceedings of IEEE International Conference on Geoscience and Remote Sensing，Boston，USA，2008：7－11.

［10］Bioucas－Dias J M，Nascimento J M P.Hyperspectral Subspace Identification［J］.IEEE Transaction on Geoscience and Remote Sensing，2008，46（8）：2435－2445.

［11］Heinz D C，Chang C.Fully Constrained Least Squares Linear Spectral Mixture Analysis Method for Material Quantification in Hyperspectral Imagery［J］.IEEE Transaction on Geoscience Remote Sensing，2001，39（3）：529－545.