高光譜圖像分類

《機(jī)器學(xué)習(xí)》課程項(xiàng)目報(bào)告高光譜圖像分類——基于CNN和ELM學(xué)院信息工程學(xué)院專業(yè)電子與通信工程學(xué)號(hào)2111603035學(xué)生姓名曹發(fā)賢同組學(xué)生陳惠明、陳濤碩士導(dǎo)師楊志景2016年11月項(xiàng)目意義與價(jià)值高光譜遙感技術(shù)起源于20世紀(jì)80年代初，是在多光譜遙感技術(shù)基礎(chǔ)之上發(fā)展起來的[1]。高光譜遙感能夠通過成像光譜儀在可見光、近紅外、短波紅外、中紅外等電磁波譜范圍獲取近似連續(xù)的光譜曲線，將表征地物幾何位置關(guān)系的空間信息與表征地物屬性特征的光譜信息有機(jī)地融合在了一起，使得提取地物的細(xì)節(jié)信息成為可能。隨著新型成像光譜儀的光譜分辨率的提高，人們對(duì)相關(guān)地物的光譜屬性特征的了解也不斷深入，許多隱藏在狹窄光譜范圍內(nèi)的地物特性逐漸被人們所發(fā)現(xiàn)，這些因素大大加速了遙感技術(shù)的發(fā)展，使高光譜遙感成為21世紀(jì)遙感技術(shù)領(lǐng)域重要的研究方向之一。在將高光譜數(shù)據(jù)應(yīng)用于各領(lǐng)域之前，必須進(jìn)行必要的數(shù)據(jù)處理。常用的數(shù)據(jù)處理技術(shù)方法包括：數(shù)據(jù)降維、目標(biāo)檢測、變化檢測等。其中，分類是遙感數(shù)據(jù)處理中比較重要的環(huán)節(jié)，分類結(jié)果不但直接提取了影像數(shù)據(jù)有效信息，可以直接運(yùn)用于實(shí)際需求中，同時(shí)也是實(shí)現(xiàn)各種應(yīng)用的前提，為后續(xù)應(yīng)用提供有用的數(shù)據(jù)信息和技術(shù)支持，如為目標(biāo)檢測提供先驗(yàn)信息、為解混合提供端元信息等。相對(duì)于多光譜遙感而言，由于高光譜遙感的波譜覆蓋范圍較寬，因此我們可以根據(jù)需要選擇特定的波段來突顯地物特征，從而能夠精確地處理地物的光譜信[2]。目前，許多國家開展大量的科研項(xiàng)目對(duì)高光譜遙感進(jìn)行研究，研制出許多不同類型的成像光譜儀。高光譜遙感正逐步從地面遙感發(fā)展到航空遙感和航天遙感，并在地圖繪制、資源勘探、農(nóng)作物監(jiān)測、精細(xì)農(nóng)業(yè)、海洋環(huán)境監(jiān)測等領(lǐng)域發(fā)揮重要的作用。高光譜遙感技術(shù)雖然是遙感領(lǐng)域的新技術(shù)，但是高光譜圖像的分類一直制約著高光譜遙感的應(yīng)用[3,4]，因此對(duì)其進(jìn)行研究顯得尤為重要。高光譜遙感圖像較高的光譜分辨率給傳統(tǒng)的圖像分類識(shí)別算法提出嚴(yán)峻的挑戰(zhàn)。波段維數(shù)的增加不僅加重了數(shù)據(jù)的存儲(chǔ)與傳輸?shù)呢?fù)擔(dān)，同時(shí)也加劇了數(shù)據(jù)處理過程的復(fù)雜性，并且由于波段與波段間存在著大量的冗余信息，從而使得傳統(tǒng)圖像分類算法并不適用于高光譜遙感圖像的分類。傳統(tǒng)的分類方法往往需要很多數(shù)目的已知類別的訓(xùn)練樣本，從而導(dǎo)致計(jì)算量大，時(shí)間復(fù)一般地，C層（卷積層）為特征提取層，每個(gè)神經(jīng)元的輸入與前一層的局部感受野相連，并提取該局部的特征，一旦該局部特征被提取后，它與其他特征間的位置關(guān)系也隨之確定下來；S層（下采樣層）是特征映射層，網(wǎng)絡(luò)的每個(gè)計(jì)算層由多個(gè)特征映射組成，每個(gè)特征映射為一個(gè)平面，平面上所有神經(jīng)元的權(quán)值相等。特征映射結(jié)構(gòu)采用sigmoid函數(shù)等映射函數(shù)作為卷積網(wǎng)絡(luò)的激活函數(shù)，使得特征映射具有位移不變性。此外，由于一個(gè)映射面上的神經(jīng)元共享權(quán)值，因而減少了網(wǎng)絡(luò)自由參數(shù)的個(gè)數(shù)，降低了網(wǎng)絡(luò)參數(shù)選擇的復(fù)雜度。卷積神經(jīng)網(wǎng)絡(luò)中的每一個(gè)特征提取層（C-層）都緊跟著一個(gè)用來求局部平均與二次提取的計(jì)算層（S-層），這種特有的兩次特征提取結(jié)構(gòu)使網(wǎng)絡(luò)在識(shí)別時(shí)對(duì)輸入樣本有較高的畸變?nèi)萑棠芰Α＞矸e神經(jīng)網(wǎng)絡(luò)采用有監(jiān)督學(xué)習(xí)的方式進(jìn)行訓(xùn)練，即任何一個(gè)訓(xùn)練樣本的類別是已知的，訓(xùn)練樣本在空間中的分布不再是依據(jù)其自然分布傾向來劃分，而是根據(jù)同類樣本和不同類樣本中的空間關(guān)系進(jìn)行劃分，這需要不斷調(diào)整網(wǎng)絡(luò)模型的參數(shù)用以劃分樣本空間的分類邊界的位置，是一個(gè)耗時(shí)且復(fù)雜的學(xué)習(xí)訓(xùn)練過程[8]。神經(jīng)網(wǎng)絡(luò)在進(jìn)行訓(xùn)練時(shí)，所有的網(wǎng)絡(luò)權(quán)值都用一些不同的小隨機(jī)數(shù)進(jìn)行初始化，這些小的隨機(jī)數(shù)能偶保證網(wǎng)絡(luò)不會(huì)因?yàn)闄?quán)值過大而進(jìn)入飽和狀態(tài)，導(dǎo)致訓(xùn)練失敗。神經(jīng)網(wǎng)絡(luò)訓(xùn)練算法包括4個(gè)主要部分：（1）樣本集中取出樣本(X,yp)并將其輸入網(wǎng)絡(luò)，X代表圖像數(shù)組，（2）計(jì)算此次輸入相應(yīng)的實(shí)際輸出Op,這是一個(gè)前向傳播的過程（3）用一個(gè)指定的損失函數(shù)計(jì)算出實(shí)際輸出Op與理想輸出Y（4）按極小化誤差的方法反向傳播調(diào)整網(wǎng)絡(luò)權(quán)值。極限學(xué)習(xí)機(jī)極限學(xué)習(xí)機(jī)（extremelearningmachine）ELM是一種簡單易用、有效的單隱層前饋神經(jīng)網(wǎng)絡(luò)SLFNs學(xué)習(xí)算法。2004年由南洋理工大學(xué)黃廣斌副教授提出。傳統(tǒng)的神經(jīng)網(wǎng)絡(luò)學(xué)習(xí)算法（如BP算法）需要人為設(shè)置大量的網(wǎng)絡(luò)訓(xùn)練參數(shù)，并且很容易產(chǎn)生局部最優(yōu)解。極限學(xué)習(xí)機(jī)只需要設(shè)置網(wǎng)絡(luò)的隱層節(jié)點(diǎn)個(gè)數(shù)，在算法執(zhí)行過程中不需要調(diào)整網(wǎng)絡(luò)的輸入權(quán)值以及隱元的偏置，并且產(chǎn)生唯一的最優(yōu)解，因此具有學(xué)習(xí)速度快且泛化性能好的優(yōu)點(diǎn)。極限學(xué)習(xí)機(jī)的網(wǎng)絡(luò)訓(xùn)練模型采用前向單隱層結(jié)構(gòu)。設(shè)分別為網(wǎng)絡(luò)輸入層、隱含層和輸出層的節(jié)點(diǎn)數(shù)，是隱層神經(jīng)元的激活函數(shù)，為閾值。設(shè)有個(gè)不同樣本，，其中，則極限學(xué)習(xí)機(jī)的網(wǎng)絡(luò)訓(xùn)練模型如圖1所示。圖1極限學(xué)習(xí)機(jī)的網(wǎng)絡(luò)訓(xùn)練模型極限學(xué)習(xí)機(jī)的網(wǎng)絡(luò)模型可用數(shù)學(xué)表達(dá)式表示如下：式中，表示連接網(wǎng)絡(luò)輸入層節(jié)點(diǎn)與第i個(gè)隱層節(jié)點(diǎn)的輸入權(quán)值向量；表示連接第i個(gè)隱層節(jié)點(diǎn)與網(wǎng)絡(luò)輸出層節(jié)點(diǎn)的輸出權(quán)值向量；表示網(wǎng)絡(luò)輸出值。極限學(xué)習(xí)機(jī)的代價(jià)函數(shù)E可表示為式中，，包含了網(wǎng)絡(luò)輸入權(quán)值及隱層節(jié)點(diǎn)閾值。Huang等指出極限學(xué)習(xí)機(jī)的懸鏈目標(biāo)就是尋求最優(yōu)的S，β，使得網(wǎng)絡(luò)輸出值與對(duì)應(yīng)實(shí)際值誤差最小，即?？蛇M(jìn)一步寫為式中，H表示網(wǎng)絡(luò)關(guān)于樣本的隱層輸出矩陣，β表示輸出權(quán)值矩陣，T表示樣本集的目標(biāo)值矩陣，H，β，T分別定義如下：極限學(xué)習(xí)機(jī)的網(wǎng)絡(luò)訓(xùn)練過程可歸結(jié)為一個(gè)非線性優(yōu)化問題。當(dāng)網(wǎng)絡(luò)隱層節(jié)點(diǎn)的激活函數(shù)無限可微時(shí)，網(wǎng)絡(luò)的輸入權(quán)值和隱層節(jié)點(diǎn)閾值可隨機(jī)賦值，此時(shí)矩陣H為一常數(shù)矩陣，極限學(xué)習(xí)機(jī)的學(xué)習(xí)過程可等價(jià)為求取線性系統(tǒng)最小范數(shù)的最小二乘解，其計(jì)算式為式中時(shí)矩陣的MP廣義逆。具體實(shí)現(xiàn)及主要代碼訓(xùn)練的樣本及其樣本圖2、原圖及進(jìn)行分類后的圖和各個(gè)高光譜數(shù)據(jù)集每一類的分類精度TheresultofPaviaUniversityTheresultofPaviaCenterTheresultofSalinas3、主要代碼ThemaincodeofPaviaUniversityofMatlab%%ConvolutionNeuralNetworkAndExtremLearningMachine%%Test1forPaviaUniversity%%STEP1:InitializecnnConfig=config();[theta,meta]=cnnInitParams(cnnConfig);n_class=cnnConfig.layer{meta.numLayers}.dimension;loadPaviaU;loadPaviaU_gt;[images,labels,indexs,samimage]=loadtrain(paviaU,paviaU_gt);d=cnnConfig.layer{1}.dimension;images=reshape(images,d(1),d(2),d(3),[]);%%STEP2:Train(TheCnnAndElm)options.epochs=1;options.minibatch=30;options.alpha=0.05;newtheta=minFuncSGD(@(x,y,z)cnnCost(x,y,z,cnnConfig,meta),theta,images,labels,options);K=cnnExtract(newtheta,images,cnnConfig,meta);[TrainingTime,TrainingAccuracy,InputWeight,BiasofHiddenNeurons,OutputWeight,NumberofOutputNeurons]=elmtrain(K,labels',1,900);%%STEP3:Test[testImages,testLabels,testIndexs]=loadtest(paviaU,paviaU_gt);testImages=reshape(testImages,d(1),d(2),d(3),[]);[row,col]=size(paviaU_gt);testK=cnnExtract(newtheta,testImages,cnnConfig,meta);[TestingTime,TestingAccuracy,testoutputlabel,actualoutputs]=elmpredict(testK,testLabels',testIndexs,1,InputWeight,BiasofHiddenNeurons,OutputWeight,NumberofOutputNeurons,row,col);predimage1=zeros(row,col);predimage1(testIndexs)=testoutputlabel;figure,imagesc(predimage1);axisoff;[OA,kappa,AA,CA]=calcError(testLabels,predimage1(testIndexs),1:n_class);fprintf('cnn+elmOverallAccuracyis%f\n',OA);fprintf('cnn+elmAverageAccuracyis%f\n',AA);fprintf('cnn+elmCAAcuuyracyis%f\n',CA);fprintf('cnn+elmkappais%f\n',kappa);figure,imagesc(paviaU_gt);axisoff;ThemaincodeofPaviaCenterofMatlab%%ConvolutionNeuralNetworkAndExtremLearningMachine%%Test2forPaviacentercnnConfig=config3();[theta,meta]=cnnInitParams(cnnConfig);n_class=cnnConfig.layer{meta.numLayers}.dimension;loadPavia;loadPavia_gt;[images,labels,indexs]=loadtrain3(pavia,pavia_gt);d=cnnConfig.layer{1}.dimension;images=reshape(images,d(1),d(2),d(3),[]);%%STEP2:Train(TheCnnAndElm)options.epochs=1;options.minibatch=25;options.alpha=0.05;newtheta=minFuncSGD(@(x,y,z)cnnCost(x,y,z,cnnConfig,meta),theta,images,labels,options);K=cnnExtract(newtheta,images,cnnConfig,meta);[TrainingTime,TrainingAccuracy,InputWeight,BiasofHiddenNeurons,OutputWeight,NumberofOutputNeurons]=elmtrain(K,labels',1,850);%%STEP3:Test[testImages,testLabels,testIndexs]=loadtest3(pavia,pavia_gt);testImages=reshape(testImages,d(1),d(2),d(3),[]);[row,col]=size(pavia_gt);testK=cnnExtract(newtheta,testImages,cnnConfig,meta);[TestingTime,TestingAccuracy,testoutputlabel]=elmpredict(testK,testLabels',1,InputWeight,BiasofHiddenNeurons,OutputWeight,NumberofOutputNeurons,row,col);predimage1=zeros(row,col);predimage1(testIndexs)=testoutputlabel;figure,imagesc(predimage1);axisoff[OA,kappa,AA,CA]=calcError(testLabels,predimage1(testIndexs),1:n_class);fprintf('cnn+elmOverallAccuracyis%f\n',OA);fprintf('cnn+elmAverageAccuracyis%f\n',AA);fprintf('cnn+elmCAAcuuyracyis%f\n',CA);fprintf('cnn+elmkappais%f\n',kappa);figure,imagesc(pavia_gt);axisoffThemaincodeofSalinasofMatlab%%ConvolutionNeuralNetworkAndExtremLearningMachine%%Test1forSalinas%%STEP1:InitializecnnConfig=config();[theta,meta]=cnnInitParams(cnnConfig);n_class=cnnConfig.layer{meta.numLayers}.dimension;loadSalinas_corrected;loadSalinas_gt;[images,labels,indexs,samimage]=loadtrain(salinas_corrected,salinas_gt);d=cnnConfig.layer{1}.dimension;images=reshape(images,d(1),d(2),d(3),[]);%%STEP2:Train(TheCnnAndElm)options.epochs=1;options.minibatch=30;options.alpha=0.05;newtheta=minFuncSGD(@(x,y,z)cnnCost(x,y,z,cnnConfig,meta),theta,images,labels,options);K=cnnExtract(newtheta,images,cnnConfig,meta);[TrainingTime,TrainingAccuracy,InputWeight,BiasofHiddenNeurons,OutputWeight,NumberofOutputNeurons]=elmtrain(K,labels',1,2300);%%STEP3:Test[testImages,testLabels,testIndexs]=loadtest(salinas_corrected,salinas_gt);testImages=reshape(testImages,d(1),d(2),d(3),[]);[row,col]=size(salinas_gt);testK=cnnExtract(newtheta,testImages,cnnConfig,meta);[TestingTime,TestingAccuracy,testoutputlabel,actualoutputs]=elmpredict(testK,testLabels',testIndexs,1,InputWeight,BiasofHiddenNeurons,OutputWeight,NumberofOutputNeurons,row,col);predimage1=zeros(row,col);predimage1(testIndexs)=testoutputlabel;figure,imagesc(predimage1);axisoff;[OA,kappa,AA,CA]=calcError(testLabels,predimage1(testIndexs),1:n_class);fprintf('cnn+elmOverallAccuracyis%f\n',OA);fprintf('cnn+elmAverageAccuracyis%f\n',AA);fprintf('cnn+elmCAAcuuyracyis%f\n',CA);fprintf('cnn+elmkappais%f\n',kappa);figure,imagesc(salinas_gt);axisoff;參考文獻(xiàn)[1]童慶禧，張兵，鄭蘭芬.高光譜遙感原理.技術(shù)與應(yīng)用[M].北京：高等教育出版社,2006[2]S.Wang,C.I.Chang.Variable-NumberVariable-BandSelectionforFeatureCharacterizationinHyperspectralSignatures[J].IEEETransanctionsonGeoscienceandRemoteSensing,2007,45(9):2979-2992[3]Q.Du,C.I.Chang.Linearmixtureanalysis-basedcompressionforhyperspectralimageanalysis[J].IEEETransanctionsonGeoscienceandRemoteSensing,2004,42(4):875-891[4]H.Erives,N.B.Targhetta.