一種基于最小體積與優(yōu)化約束條件的高光譜解混方法
【技術(shù)領(lǐng)域】
[0001] 本發(fā)明涉及一種高光譜數(shù)據(jù)解混方法，具體涉及一種基于幾何學(xué)最小體積算法解 混的方法。
【背景技術(shù)】
[0002] 高光譜遙感能利用很窄的電磁波波段從感興趣的物體上獲取有關(guān)數(shù)據(jù)。高光譜遙 感是在電磁波譜的可見光、近紅外、中紅外和熱紅外波段范圍內(nèi)，獲取許多非常窄的連續(xù)光 譜數(shù)據(jù)的技術(shù)，目前先進(jìn)的成像光譜儀可以收集到成百上千個非常窄的光譜波段信息。
[0003] 利用高光譜解混技術(shù)能夠確定高光譜圖像中的端元數(shù)、端元光譜曲線以及每個像 素點的豐度值，進(jìn)而能夠確定某種物質(zhì)在圖像觀察區(qū)域的分布。在解混方法中，高光譜的混 合模型分為兩類，分別為線性混合模型以及非線性混合模型。線性模型假設(shè)物質(zhì)成份在空 間布局上是離散混合，忽略在不同類型物質(zhì)之間的多重散射量，即每個像元點的光譜幅度 僅是這個像元中存在的各物質(zhì)成份光譜信號，以相對貢獻(xiàn)量多少的線性疊加。對于高光譜 線性混合模型來說，解混方法可以分為三大類，分別為：1)基于凸面幾何學(xué)方法的線性解混 模型；2)基于數(shù)據(jù)統(tǒng)計的方法，代表算法有ICA等，其核心思想在于假設(shè)原信號在體積上是 獨立分布的；3)基于稀疏模型的解混方法，代表算法有BP、BTON等。
[0004] 基于幾何學(xué)的解混算法又能夠被分為兩類:一類是基于純像元假設(shè)的解混算法， 另一類是基于最小體積法的解混算法?；诩兿裨僭O(shè)的解混算法需要保證圖像中的每個 不同的端元至少完全占據(jù)一個像素，否則會出現(xiàn)端元錯誤估計的情況。而基于最小體積的 解混算法即使在缺失純像元的情況下也能夠較高地估計出端元，但缺點是其優(yōu)化目標(biāo)函數(shù) 是非線性且非凸的，很多情況下，只能得到近似的局部最優(yōu)解，同時其約束條件數(shù)量非常 多，使其算法效率沒有基于純像元假設(shè)的解混算法高。

【發(fā)明內(nèi)容】

[0005] 為了解決基于最小體積的高光譜解混算法中，約束條件過多所導(dǎo)致的數(shù)據(jù)存儲空 間大，運算時間長，精確性難以提高的問題，本發(fā)明提供了一種基于最小體積與優(yōu)化約束條 件的高光譜解混方法，通過優(yōu)化約束條件改進(jìn)基于最小體積法的高光譜解混方法，得出更 準(zhǔn)確的解混效果。
[0006] 本發(fā)明的目的是通過以下技術(shù)方案實現(xiàn)的：
[0007] -種基于最小體積與優(yōu)化約束條件的高光譜解混方法，對高光譜圖像進(jìn)行線性解 混以求取端元矩陣時，考慮初值對求解端元矩陣的影響，并將求解體積最小問題的非線性 問題轉(zhuǎn)化為線性規(guī)劃問題，近似求解得到非線性優(yōu)化的局部最小解;并考慮求解端元的非 負(fù)物理特性，從而提高求取端元的準(zhǔn)確率。具體包括如下步驟：
[0008] 步驟一:數(shù)據(jù)加載與預(yù)處理。
[0009]步驟二:篩選圖像采樣點，構(gòu)造優(yōu)化的約束條件，尋找滿足初始條件的數(shù)據(jù)。
[0010]步驟三:將非負(fù)非線性規(guī)劃問題轉(zhuǎn)換為線性規(guī)劃問題，對線性規(guī)劃問題的目標(biāo)函 數(shù)結(jié)合優(yōu)化后的約束條件進(jìn)行求解，計算中間變量矩陣Hne3W、gn?。
[0011] 步驟四：根據(jù)變化率檢測終止條件判斷是否終止迭代計算，若不滿足終止條件，則 返回步驟三，繼續(xù)更新中間變量矩陣11_、g ne3W。
[0012] 步驟五:解出滿足非負(fù)性要求的端元矩陣，并計算豐度系數(shù)，完成圖像的解混。
[0013] 本發(fā)明與現(xiàn)有技術(shù)相比具有如下優(yōu)點：
[0014] 1)通過篩選數(shù)據(jù)采樣點，同時利用某一采樣點必然能包含在由其他頂點構(gòu)成的單 形體內(nèi)的性質(zhì)，有效縮減了約束條件的個數(shù)使得計算開銷變小，同時解混效率有較大提高。
[0015] 2)不論是否存在純像元，本方法都能夠經(jīng)過解算，得到高光譜圖像的端元，進(jìn)而得 到豐度系數(shù)，完成高光譜圖像的解混。
【附圖說明】
[0016] 圖1為本發(fā)明的流程圖；
[0017] 圖2為數(shù)據(jù)點集、真實端元集、估計端元集在2(P-1)維上的投影顯示；
[0018]圖3為真實端元光譜圖與估計端元光譜圖在224個波段上的顯示。
【具體實施方式】
[0019] 下面結(jié)合附圖對本發(fā)明的技術(shù)方案作進(jìn)一步的說明，但并不局限于此，凡是對本 發(fā)明技術(shù)方案進(jìn)行修改或者等同替換，而不脫離本發(fā)明技術(shù)方案的精神和范圍，均應(yīng)涵蓋 在本發(fā)明的保護(hù)范圍中。
【具體實施方式】 [0020] 一:本實施方式提供了一種基于最小體積與優(yōu)化約束條件的高光譜 解混方法，假設(shè)光譜圖像數(shù)據(jù)FVR&Ar，L是波段數(shù)，N是圖像的像素點數(shù)量。若在光譜圖像 中有P個端元，首先對光譜信息進(jìn)行降維處理，將數(shù)據(jù)從L維映射到p-1維;接著構(gòu)造一個單 形體，該單形體能夠包含所有降維后的數(shù)據(jù)，通過循環(huán)求解線性規(guī)劃子問題，不斷更新單形 體頂點矩陣，最終使單形體的體積收斂至穩(wěn)定，最后通過單形體的頂點矩陣得到圖像的端 元矩陣，進(jìn)而求出圖像的豐度系數(shù)，完成高光譜圖像的解混任務(wù)。
[0021 ]如圖1所示，共分為五個步驟，具體步驟如下：
[0022]步驟一:數(shù)據(jù)加載與預(yù)處理。
[0023] 1)加載數(shù)據(jù)矩陣Γ0 ，其中L是波段數(shù)，N是圖像的像素點數(shù)量；
[0024] 2)對數(shù)據(jù)矩陣y 行降維。已知光譜圖像中，端元個數(shù)為p，將數(shù)據(jù)矩陣Y降 維至p-Ι維數(shù)據(jù)矩陣f e &具體降維方法為：
[0025] f =QfU (.1) ?
[0026] 其中，U=[yi-d，y2_d，. . .，yN_d]:
.[，數(shù)據(jù)矩陣Y表示為列向量的形式： Y=[yi,y2, . . . ,Υν] ； O' ={Ql Q) ?1 ，Q= [qi(UUT)，q2(UUT)，· · ·，qP-1(UUT)]，qi(UUT)，i = 1，2, . . .，p_l表示矩陣1]1^第1大的特征值所對應(yīng)的特征向量。
[0027] 步驟二:采用優(yōu)化后的約束條件，尋找滿足初始條件的數(shù)據(jù)。
[0028] 1)在數(shù)據(jù)矩陣f中得到p個列向量：
[0029]
[0030] 式中，ΛΑ E鈀〃1~，向量選擇滿足：由這p個列向量作為頂點構(gòu)成的單形體體積比 其他任意P個列向量作為頂點構(gòu)成的單形體體積大的條件。
[0031] 2)在數(shù)據(jù)集f中尋找數(shù)據(jù)點集，數(shù)據(jù)集中任意一點滿足不在由M〇中的列向量作為 頂點構(gòu)成的單形體中的條件，最后得到滿足要求的數(shù)據(jù)點集r。
[0032]具體尋找方法如下，首先令：
[0033]
[0034]
[0035]
[0036] 其中，B e 股(;，-1)x(p-11，天 e ，si為s中的第i行元素。
[0037] 對所有兄e f, 1 ?? / ?? ;V (N為像素點個數(shù)）求取其對應(yīng)的向量i，若向量i不滿足 0_.·< ? <1的條件，則將另加入數(shù)據(jù)集合.fWwA_.中。假設(shè)最終有m(顯然p < m < N，且一般情況下m <<N)個不滿足上述條件的免，則最后得到的數(shù)據(jù)集合
[0038] 3)優(yōu)化的線性不等式約束條件為式(6)與式(7):
[0039]
[0040]
[0041]
[0042] 其中，qERlx(P-D，e W，： g eR("-11x1，Η與g為要優(yōu)化的變量矩陣。式(6) 與式(7)將不等式約束個數(shù)降為了 (p_l)X(m+l)個，不等式約束的個數(shù)被大大減小，從而大 大提高了求解的速度。
[0043] 本步驟通過求解線性規(guī)劃問題的方法，求解出一個滿足線性不等式約束式(6)與 式(7)關(guān)于變量矩陣Η與g的一組可行解矩陣Ho與go,作為后續(xù)步驟的初始值。
[0044] 本步驟中的<、<、2分別表示所有元素都<、<、2。
[0045] 步驟三:對線性規(guī)劃問題的目標(biāo)函數(shù)結(jié)合優(yōu)化后的約束條件進(jìn)行求解，計算中間 變量矩陣Hne3W、g n?。迭代計算中，變量矩陣Η與g的階段性最優(yōu)變量值存儲在矩陣Hne3W、gne3W中。
[0046] 構(gòu)造優(yōu)化函數(shù)，通過代數(shù)余子式的展開，將非線性的優(yōu)化問題改造成為線性規(guī)劃 子問題，通過周期循環(huán)迭代的方式，更新得到新的矩陣Hne3W、gne3W，直到滿足停止迭代的條件。 [0047] 1)構(gòu)造非線性目標(biāo)函數(shù)：
[0048]
[0049] 2)將非線性目標(biāo)函數(shù)轉(zhuǎn)化為線性的目標(biāo)函數(shù)：
[0050] 首先將Η的行列式det(H)通過代數(shù)余子式的方式展開：
[0051]
[0052] 式中是矩陣Η的第i行第j列的元素，是將矩陣Η移除第i行元素與第j列元素 后得到的子陣。進(jìn)而目標(biāo)函數(shù)轉(zhuǎn)化為：
[0053]
.
[0054]由于目標(biāo)優(yōu)化函數(shù)| det(H) |仍然是非線性的，因此將其分解為下面兩個線性目標(biāo) 函數(shù)：
[0055]
[0056]
[0057] 3)將線性目標(biāo)函數(shù)與優(yōu)化后的線性約束條件結(jié)合，構(gòu)造如下雙線性規(guī)劃子問題：
[0058]
[0059]
[0060]
[0061] .V/./7； r；i[/]-g, >0,'^.1',[/]￡ ΚΛ,,ι；ι,1；; ,// -1；,
[0062] 其中，lu為矩陣Η第i行的行向量，gl為向量g的第i個元素。
[0063] 4)求解上述兩個線性規(guī)劃子問題，根據(jù)f與cf的值，選擇性地更新對應(yīng)的矩陣元素 h^gl，直到將矩陣Η與g中所有的元素都更新一遍。
[0064] 具體更新hggl的規(guī)則為:若求解線性規(guī)劃問題式(14)得到的解為ξ、瓦，求解線 性規(guī)劃問題式（15)得到的解為，若|p1 > |q$|，則將hggl更新為[與瓦，否則，將lu與 gi更新為t與gi。
[0065] 在這一步驟中，逐行更新Η與g的值:令i = l, ...，p-l，依次得到hi與gi的值，最終得 到矩陣Η與g中所有元素都更新后的新矩陣!1_與gne3W。
[0066] 步驟四：根據(jù)變化率檢測終止條件判斷是否終止迭代計算，若不滿足終止條件，則 返回步驟三，繼續(xù)更新中間變量矩陣11_，g ne3W。
[0067] 1)設(shè)置變化率檢測終止條件如下：
[0068] { |det(Hnew) |-|det(Ho) | }/|det(H) | <ε (16)。
[0069] 式中，Ho代表在步驟三開始迭代之前，矩陣Η的初值。
[0070] 2)判斷中間變量矩陣Hnew是否滿足條件式（16)，如果滿足終止條件，則令!f = Hnew， g* = gnew，Η*與g*作為矩陣Η與g的最終更新結(jié)果。否貝lj，令H=Hnew，Ho = Hnew，g = gnew，轉(zhuǎn)入步驟 三的第4)步。
[0071 ]步驟五:還原端元矩陣，計算豐度系數(shù)，完成圖像的解混。
[0072] 根據(jù)f與g'還原端元矩陣，并計算豐度系數(shù)，完成圖像的解混，端元矩陣應(yīng)滿足非 負(fù)性的要求。
[0073] 1)令αρ = (Η*)-4*，[?丨,…丨]丨十（/廣）丨，α = [