袁博，張杰林

（1.中國(guó)礦業(yè)大學(xué)（北京），地球科學(xué)與測(cè)繪工程學(xué)院，北京100083；2.核工業(yè)北京地質(zhì)研究院，遙感信息與圖像分析技術(shù)國(guó)家級(jí)重點(diǎn)實(shí)驗(yàn)室，北京100029）

一種改進(jìn)的頂點(diǎn)成分分析端元提取算法

袁博1，張杰林2

（1.中國(guó)礦業(yè)大學(xué)（北京），地球科學(xué)與測(cè)繪工程學(xué)院，北京100083；2.核工業(yè)北京地質(zhì)研究院，遙感信息與圖像分析技術(shù)國(guó)家級(jí)重點(diǎn)實(shí)驗(yàn)室，北京100029）

頂點(diǎn)成分分析算法需要預(yù)先提供端元數(shù)目，端元數(shù)目正確與否對(duì)結(jié)果會(huì)產(chǎn)生較大影響，其算法在實(shí)際應(yīng)用中多次運(yùn)行的結(jié)果不穩(wěn)定。針對(duì)上述缺點(diǎn)提出了一種改進(jìn)的頂點(diǎn)成分分析端元提取算法。該方法在n維光譜空間中生成n個(gè)彼此正交的單位向量，在此基礎(chǔ)上生成與之具有一定夾角的單位向量，將光譜空間中的像元點(diǎn)分別投影在單位向量上以獲取端元。結(jié)果表明改進(jìn)的頂點(diǎn)成分分析端元算法提高了端元提取結(jié)果的穩(wěn)定性。

端元提取；頂點(diǎn)成分分析；純凈像元指數(shù)；正交向量；蝕變信息

高光譜遙感技術(shù)是20世紀(jì)80年代發(fā)展起來(lái)的首次實(shí)現(xiàn)圖譜合一的對(duì)地觀測(cè)技術(shù)，其光譜分辨率可以達(dá)到5～10 nm［1］，可為每個(gè)像元提供數(shù)十個(gè)或數(shù)百個(gè)波段的光譜信息，可以產(chǎn)生一條完整并且連續(xù)的光譜曲線，使得許多在多光譜遙感中不能分辨的物質(zhì)得以識(shí)別。地球表面并不是由單一的物質(zhì)組成的，當(dāng)遙感圖像中的一個(gè)像元包含了不同波譜特性的物質(zhì)時(shí)就會(huì)產(chǎn)生混合像元，在高光譜遙感中混合像元是廣泛存在的，高光譜遙感圖像中每個(gè)像元的光譜信息是其對(duì)應(yīng)地表物質(zhì)光譜信息的綜合，不同的物質(zhì)具有不同的光譜響應(yīng)特征，而一個(gè)像元僅用一條波譜曲線來(lái)代表這些特征［2］，因此對(duì)高光譜遙感數(shù)據(jù)進(jìn)行分析就必須進(jìn)行混合像元的解混。當(dāng)像元中只有一種地物時(shí)，該像元就被稱為端元，端元光譜就是該區(qū)域的特征光譜，因此端元提取是解混的前提。

像元混合分為線性混合和非線性混合。如果像元內(nèi)物質(zhì)混合的尺度大，那么混合像元被視為線性混合；如果像元內(nèi)物質(zhì)混合尺度小，混合緊密，則混合像元被視為非線性混合［3］。因?yàn)椴煌瑫r(shí)間地點(diǎn)獲得的高光譜數(shù)據(jù)具體的非線性混合模式復(fù)雜多變，通常很難找到與實(shí)際情況相吻合的非線性混合模型，而線性混合模型具有簡(jiǎn)單易于操作并且精度滿足需要的優(yōu)點(diǎn)，因此目前線性混合模型應(yīng)用比較廣泛。頂點(diǎn)成分分析（Vertex Component Analysis，VCA）是在線性混合模型基礎(chǔ)上提出的端元提取算法，該算法需要提前確定端元數(shù)目，并且多次提取端元的結(jié)果不穩(wěn)定，針對(duì)VCA的這一缺點(diǎn)筆者提出了一種改進(jìn)的頂點(diǎn)成分分析端元提取算法。

1 端元提取算法研究現(xiàn)狀

目前比較常用的端元提取算法包括純凈像元指數(shù)法（Pure Pixel Index，PPI）、內(nèi)部體積最大法（N-FINDR）、單形體體積增長(zhǎng)算法（Simplex Growing Algorithm，SGA）、VCA算法等［4］。

PPI算法是Boardman（1995）提出的在多光譜和高光譜影像中尋找波譜最純凈的像元，從而提取端元的方法［2］。該算法可以從高光譜數(shù)據(jù)中分離出光譜較純凈的像元，減少參與端元選取的像元數(shù)量。之后利用純凈像元生成n維散點(diǎn)圖（n指純凈像元個(gè)數(shù)），進(jìn)行交互式可視化操作選出圖像端元［5］。PPI方法使用隨機(jī)初始化向量進(jìn)行投影，具有隨意性，不可避免地產(chǎn)生冗余向量，不能有效區(qū)分各個(gè)端元的投影，并且該方法用到大量的投影變換，計(jì)算效率較低，消耗時(shí)間較長(zhǎng)。

Winter（1999）提出了內(nèi)部體積最大法（NFINDR），通過(guò)尋找具有最大體積的單形體自動(dòng)地逐步獲取影像中的端元［6］。該方法首先隨機(jī)選擇一定數(shù)量的像元作為端元，并計(jì)算以這些端元為頂點(diǎn)構(gòu)成的單形體的體積，再將新的像元依次代替端元，重新計(jì)算單形體的體積，如果體積增大，則該像元替代之前端元，直到?jīng)]有像元能使單形體的體積增大時(shí)可以得到的所有的端元［6］。但是當(dāng)數(shù)據(jù)存在噪聲時(shí)，該算法使得一些點(diǎn)將不可避免地出現(xiàn)在單形體之外，因此數(shù)據(jù)含有噪聲時(shí)結(jié)果準(zhǔn)確性將降低。

SGA算法是對(duì)N-FINDR算法的改進(jìn)，通過(guò)順序搜索方法來(lái)尋找構(gòu)成最大單形體體積的端元［7］。該方法需要重復(fù)的體積計(jì)算，計(jì)算復(fù)雜度比較高，因此效率較低。

VCA算法理論基礎(chǔ)為線性混合模型，一個(gè)像元矢量的光譜信號(hào)用線性混合模型描述為［8］:

式（1）:B—L×1維的觀測(cè)波譜，其中L—高光譜圖像的波段數(shù)；M＝［m1，m2，…，me］—一個(gè)L×e的端元集矩陣，其中mi（i＝1，2…e）—端元的光譜，其中e—端元個(gè)數(shù)；k—形態(tài)學(xué)比例因子；A＝［a1，a2，…，ae］T，—e×1維的豐度值向量；n—存在于觀測(cè)波譜B各譜段中的L×1維高斯白噪聲。

在線性混合模型下，觀測(cè)向量集合組成一個(gè)凸錐C＝{B∈RL:R＝kMA，1TA＝1，A≧0，k≧0}。凸錐C投影到超平面RTu＝1形成一個(gè)單形體S＝{y∈RL:y＝R/（RTu），R∈CR}。VCA首先隨機(jī)生成一個(gè)初始向量，找到單形體S中在初始向量上具有最大投影長(zhǎng)度的像元作為第一個(gè)端元，然后計(jì)算與該端元正交的向量作為第二個(gè)投影向量，將數(shù)據(jù)投影找到第二個(gè)端元。然后依次生成新的與已提取端元所張成的空間相正交投影向量，提取相應(yīng)的端元，當(dāng)端元數(shù)目滿足預(yù)先確定的數(shù)目e時(shí)算法結(jié)束［9］。以圖1中兩個(gè)端元為例，第一次迭代中，單形體數(shù)據(jù)被投影到第一個(gè)方向f1，最大的投影向量對(duì)應(yīng)著端元m1。接下來(lái)的迭代中，將數(shù)據(jù)投影到與端元m1正交的方向f2，此時(shí)最大的投影向量便對(duì)應(yīng)著端元m2。

VCA算法需要事先確定要提取的端元數(shù)目，端元數(shù)目正確與否對(duì)豐度值求解具有較大影響，但在實(shí)際很多情況下不能預(yù)先確定研究區(qū)域內(nèi)端元個(gè)數(shù)，此時(shí)VCA方法不具有適用性［8］。另外VCA方法在沒(méi)有誤差和噪聲的理論前提下經(jīng)過(guò)n次投影確定出n個(gè)端元，但是高光譜數(shù)據(jù)中不可避免地含有噪聲，實(shí)際投影所得結(jié)果可能是噪聲，因此只憑一次投影就確定一個(gè)端元會(huì)存在很大的誤差。通過(guò)VCA方法多次確定出的端元集差異性較大，其多次運(yùn)行結(jié)果不穩(wěn)定，該方法的準(zhǔn)確率較低。

何曉寧等在2013年對(duì)VCA算法進(jìn)行了改進(jìn)。該改進(jìn)算法利用圖像空間信息排除噪聲，確定候選端元，再通過(guò)病態(tài)矩陣分析候選端元矩陣的有效性，經(jīng)過(guò)迭代后獲得端元［8］。該方法需要在計(jì)算前提供端元數(shù)目，但是在實(shí)際應(yīng)用中并不總能確定端元數(shù)目，并且該算法具有計(jì)算量大的缺點(diǎn)。

圖1 VCA原理示意圖Fig.1Schematic diagram of VCA

2 對(duì)頂點(diǎn)成分分析（VCA）的改進(jìn)算法

高光譜數(shù)據(jù)可以被視為在n維空間中的點(diǎn)，點(diǎn)的坐標(biāo)是各個(gè)波段的反射率值。改進(jìn)的VCA端元提取算法首先對(duì)數(shù)據(jù)進(jìn)行最小噪聲分離（Minimum Noise Fraction，MNF）變換，去除數(shù)據(jù)中的噪聲，確定含有有用信息的最小波段數(shù)n，然后在變換后的n維空間中生成一個(gè)不與原始數(shù)據(jù)集正交的單位投影向量a1，再計(jì)算出與向量a1正交的單位向量a2作為第二投影向量，此時(shí)單位向量a1和單位向量a2組成一個(gè)二維向量空間。生成第三個(gè)單位向量a3，使之與向量a1和a2所張成的空間正交。以此類推直到第n個(gè)單位向量an，使an與a1、a2、…、an-1構(gòu)成的n-1維向量空間正交，這里n指波段數(shù)，此時(shí)得到n個(gè)彼此之間相互正交的單位向量，即找到了n維空間中的一組基向量，n維空間中的任何像元點(diǎn)都可以有效地投影到單位向量上。

為了保證有效的提取端元，需要適當(dāng)增加投影向量的數(shù)量。計(jì)算與已知n個(gè)向量中每?jī)蓚€(gè)向量ai和aj均成π/4夾角的中間單位向量aij，i＝1…n，j＝i+1…n，經(jīng)計(jì)算可得到n（n-1）/2個(gè)中間單位向量aij，于是最終得到n＋n（n－1）/2個(gè)投影向量，n為選取的波段數(shù)。最后將n維空間中的所有像元點(diǎn)投影到單位向量上，記錄到每個(gè)點(diǎn)投影后為最大投影的次數(shù)，然后設(shè)定閾值，選擇最大投影次數(shù)大于閾值的像元為端元。

圖2 改進(jìn)的VCA原理示意圖Fig.2Schematic diagram of improved VCA

改進(jìn)的VCA端元提取算法原理如圖2所示，ai為生成的單位向量，aj為與ai正交的單位向量，aij為與ai、aj向量夾角均為π/4的單位向量，圖中的點(diǎn)1、點(diǎn)2、點(diǎn)3代表3個(gè)不同像元在n維波譜空間的位置，虛線表示投影方向。記錄像元點(diǎn)在各個(gè)向量上投影后具有最大投影的次數(shù)，點(diǎn)1被記錄為最大投影的次數(shù)為1，點(diǎn)2被記錄為最大投影的次數(shù)為2，點(diǎn)3被記錄為最大投影的次數(shù)為3。

在該算法計(jì)算之前，高光譜數(shù)據(jù)首先經(jīng)過(guò)MNF處理，選取特征值較大數(shù)十個(gè)波段，總的投影向量個(gè)數(shù)控制在數(shù)百個(gè)。該方法原理簡(jiǎn)單，確保了n維空間中不同位置頂點(diǎn)有同等機(jī)率投影到端點(diǎn)位置，最終選取數(shù)百個(gè)投影向量，相比VCA方法，適當(dāng)增加了投影向量，避免了投影向量產(chǎn)生的隨意性，具有一定的抗噪聲能力，使提取的端元信息更準(zhǔn)確，保證了端元提取結(jié)果的穩(wěn)定性。

3 實(shí)驗(yàn)分析

3.1 數(shù)據(jù)來(lái)源

在IDL8.0環(huán)境下實(shí)現(xiàn)改進(jìn)的VCA端元提取算法，并采用核工業(yè)北京地質(zhì)研究院采集的SASI航空高光譜遙感測(cè)量數(shù)據(jù)驗(yàn)證算法的實(shí)用性。實(shí)驗(yàn)區(qū)域的面積約為200 km2。

SASI是加拿大生產(chǎn)的機(jī)載短波紅外（SWIR）高光譜成像系統(tǒng)，其主要技術(shù)指標(biāo)見(jiàn)表1。

表1 SASI航空成像光譜儀主要技術(shù)指標(biāo)Table 1Main technical performances of SASI airborne imaging spectroradiometer

因?yàn)镾ASI航空高光譜數(shù)據(jù)中1 370～1 415 nm及1 835～1 925 nm波譜范圍內(nèi)共有11個(gè)波段受大氣水汽影響嚴(yán)重，故將其去除。最終選擇950～1 355 nm、1 430～1 820 nm、1 940～2 450 nm共90個(gè)波段參與計(jì)算。

3.2 數(shù)據(jù)處理與結(jié)果

首先對(duì)數(shù)據(jù)進(jìn)行預(yù)處理，采用統(tǒng)計(jì)學(xué)模型中的經(jīng)驗(yàn)線性法（Empirical line），利用已測(cè)得的白板數(shù)據(jù)作為已知點(diǎn)的地面光譜對(duì)高光譜影像進(jìn)行反射率的計(jì)算。之后對(duì)數(shù)據(jù)進(jìn)行MNF操作去除噪聲，得到特征值分布（圖3）。

圖3 MNF特征值分布曲線Fig.3Curve of the eigenvalues of MNF

根據(jù)圖3中的特征值曲線并觀察MNF圖像可知處于20右側(cè)的波段趨于一致為噪聲波段，因此選擇前20個(gè)波段參與之后的運(yùn)算。將選取的波段帶入改進(jìn)的VCA提取端元算法，經(jīng)過(guò)210次迭代運(yùn)算，程序自動(dòng)選擇出端元，生成基于空間位置的端元分布圖，為了驗(yàn)證本算法的穩(wěn)定性，經(jīng)過(guò)多次計(jì)算生成的端元空間位置分布見(jiàn)圖4。

圖4顯示改進(jìn)的VCA算法多次提取結(jié)果的端元分布圖中端元位置保持一致，驗(yàn)證了改進(jìn)算法的結(jié)果具有穩(wěn)定性。而VCA算法多次提取端元的結(jié)果不一致，并且不能將圖像中的端元全部提取出來(lái)。

通過(guò)改進(jìn)的VCA端元提取算法可以獲得端元波譜曲線，再利用ENVI軟件中自帶的波譜分析功能將端元波譜與USGS波譜庫(kù)中數(shù)據(jù)進(jìn)行對(duì)比分析可以識(shí)別出各個(gè)端元。因?yàn)檫x取的SASI數(shù)據(jù)為90個(gè)波段，光譜分辨率為15 nm，所以該步驟需要事先將USGS波譜庫(kù)數(shù)據(jù)重采樣到同樣波段數(shù)和光譜分辨率才可以進(jìn)行。對(duì)比USGS波譜庫(kù)中光譜曲線，經(jīng)分析可以得到圖像中的端元包括白云母、片沸石、黃鉀鐵礬和深綠輝石4種蝕變礦物。

白云母（K2Al4（Si6Al2O20）（OH，F(xiàn)）4）具有440 nm、900 nm Fe2＋吸收譜帶；1 400 nm-OH倍頻、合頻譜帶；在2 200 nm、2 350 nm、2 450 nm存在吸收特征［10］，其診斷性光譜特征是Al-OH鍵在2 180～2 228 nm之間的尖銳而深的單一吸收特征，以及在2 340 nm和2 440 nm附近的較弱Al-OH特征［11］；片沸石（（Ca，Na）2～3Al3（Al，Si）2Si13O36·12H2O），具有965 nm、1 160 nm、1 430 nm、1 940 nm水引起的吸收特征；黃鉀鐵礬（KFe3（OH）6（SO4）2），具有1 400 nm、1 900 nm、2 500 nm-OH振動(dòng)譜帶，900 nm Fe3＋吸收波譜特征［10］；深綠輝石（Ca（Mg，F(xiàn)e3＋，Al）（Si，Al）2O6），具有1 050 nm Fe2＋譜帶和850 nm、450 nm Fe3＋譜帶［10］。

由改進(jìn)算法提取的端元波譜與VCA提取的端元波譜及波譜庫(kù)波譜對(duì)比（圖5）可見(jiàn)改進(jìn)算法比VCA算法提取的波譜更接近于波譜庫(kù)。

以所提取的4種端元的光譜曲線為基礎(chǔ)利用混合調(diào)制匹配濾波方法進(jìn)行蝕變信息的填圖，經(jīng)運(yùn)算可得蝕變信息提取結(jié)果（圖6）。結(jié)果表明以改進(jìn)算法提取的端元光譜為依據(jù)能更準(zhǔn)確更詳細(xì)地提取出蝕變信息。

圖5 光譜曲線對(duì)比Fig.5Comparison of spectral curve

圖6 蝕變信息提取結(jié)果Fig.6Extracting results of alteration information

（，Continued on page 44）（，Continued from page 38）

4 結(jié)論

提出了改進(jìn)的VCA端元提取算法，計(jì)算n維空間的一組基向量并選擇與之成特定夾角的單位向量共同作為投影向量，在選出的數(shù)百個(gè)單位向量上進(jìn)行像元投影，確定閾值選出端元。改進(jìn)的VCA算法無(wú)需預(yù)先確定圖像的端元數(shù)目，通過(guò)適量投影自動(dòng)選擇端元，算法同時(shí)改進(jìn)了VCA算法不具有抗噪聲能力的缺點(diǎn)，在保證計(jì)算效率的同時(shí)提高了端元提取結(jié)果的穩(wěn)定性。

［1］芶盛.高光譜遙感圖像光譜特征提取與匹配技術(shù)研究［D］.成都：成都理工大學(xué)，2011.

［2］林娜，楊武年，劉漢湖.基于高光譜遙感的巖礦端元識(shí)別及信息提取研究［J］.遙感信息，2011（5）：114-117.

［3］Hassan N，Hashim M.Decomposition of mixed pixels of ASTER satellite data for mapping Chengal（Neobalanocarpusheimiisp）tree［C］//Control System，Computing and Engineering（ICCSCE），2011 IEEE International Conference on.Penang：IEEE，2011：74-79.

［4］夏威.高光譜遙感圖像解混和波段選擇方法研究［D］.上海:復(fù)旦大學(xué)，2013.

［5］Boardman J W，Kruse F.Analysis of Imaging Spectrometer Data Using N-Dimensional Geometry and a Mixture-Tuned Matched Filtering Approach［J］.GeoscienceandRemoteSensing，IEEE Transactions on，2011，49（11）：4 138-4 152.

［6］錢進(jìn)，鄧喀中，劉冬.基于優(yōu)化N-FINDR算法的高光譜遙感影像礦物識(shí)別［J］.金屬礦山，2012，08：84-87，91.

［7］王麗姣.基于單形體體積增長(zhǎng)的高光譜圖像端元提取及快速實(shí)現(xiàn)［D］.杭州:浙江大學(xué)，2015.

［8］何曉寧，曹建農(nóng)，高坡，等.基于頂點(diǎn)成分分析法的端元提取改進(jìn)算法［J］.測(cè)繪通報(bào)，2013（7）：30-34.

［9］Nascimento J M P，Dias J M B.Vertex component analysis：A fast algorithm to unmix hyperspectral data［J］.IEEE Transactions on Geoscience and Re mote Sensing，2005，43（4）：898-910.

［10］燕守勛，張兵，趙永超，等.礦物與巖石的可見(jiàn)-近紅外光譜特性綜述［J］.遙感技術(shù)與應(yīng)用，2003，18（4）：191-201.

［11］劉圣偉，甘甫平，閆柏琨，等.成像光譜技術(shù)在典型蝕變礦物識(shí)別和填圖中的應(yīng)用［J］.中國(guó)地質(zhì)，2006，33（1）：178-186.

An improved algorithm of end member extraction based on Vertex Component Analysis

YUAN Bo1，ZHANG Jielin2
（1.China University of Mining&Technology（Beijing），College of Geoscience and Surveying Engineering，Beijing 100083，China；2.National Key Laboratory of Remote Sensing Information and Image Analysis Technology，Beijing Research Institute of Uranium Geology，Beijing 100029，China）

Vertex Component Analysis algorithm needs to determine the number of end members at the beginning，therefore the determined number has great influence on the extract results，and the algorithm was found not stable in practical applications.To overcome the shortcomings of Vertex Component Analysis，this paper proposed an improved algorithm.This improved algorithm first find orthogonal unit vectors in the n-dimensional spectral space，and generates a number of unit vectors which have certain angle with the original unit vectors，then pixel points in the spectral space are projected on all of the unit vectors to extract end members.The result shows that this improved algorithm improves stability of the extraction result of end members.

extraction of endmembers；Vertex Component Analysis；Pure Pixel Index；orthogonal vectors；alteration information

TP79；P57

A

1672-0636-（2016）01-0033-06

10．3969/j．issn．1672-0636．2016．01．006

2015-08-20

袁博（1989—），男，河北保定人，中國(guó)礦業(yè)大學(xué)（北京）碩士研究生；研究方向:高光譜遙感。E-mail:boziaaaaa@163.com