劉 震，朱耀琴

（南京理工大學(xué) 計(jì)算機(jī)科學(xué)與工程學(xué)院，江蘇 南京210094）

一種基于Spark的高光譜遙感圖像分類并行化方法

劉 震，朱耀琴

（南京理工大學(xué) 計(jì)算機(jī)科學(xué)與工程學(xué)院，江蘇 南京210094）

為了實(shí)現(xiàn)大數(shù)據(jù)量遙感圖像的分類，提出了一種Spark平臺(tái)下高光譜遙感圖像稀疏表示分類并行化方法PSCSRC，首先設(shè)計(jì)五元組形式的中間數(shù)據(jù)存儲(chǔ)結(jié)構(gòu)，其次在每次迭代中只是將稀疏矩陣中與中間矩陣各分塊對(duì)應(yīng)的子矩陣分發(fā)到各節(jié)點(diǎn)，這就減少了中間數(shù)據(jù)Shuffle和數(shù)據(jù)傳遞的時(shí)間消耗。實(shí)驗(yàn)結(jié)果表明，PSCSRC分類方法與SCSRC分類方法相比，在保證分類精度的同時(shí)，執(zhí)行速度有了明顯的提升。

高光譜；云計(jì)算；Spark平臺(tái)；分類；稀疏表示

目前，全球有超過1000顆衛(wèi)星實(shí)現(xiàn)對(duì)地觀測(cè)的任務(wù)[1]。伴隨著遙感技術(shù)的快速發(fā)展，遙感圖像的空間和光譜分辨率也越來越高，可以獲取的遙感數(shù)據(jù)量也迅猛增長。遙感數(shù)據(jù)量的爆炸性增長給傳統(tǒng)分類方法帶來了巨大的挑戰(zhàn)。為了解決該問題，目前常用的方法有基于GPU加速的高性能計(jì)算系統(tǒng)和基于云計(jì)算平臺(tái)的分布式集群。

GPU在浮點(diǎn)運(yùn)算和并行計(jì)算方面相對(duì)于CPU有著巨大的優(yōu)勢(shì)，因此有學(xué)者利用GPU技術(shù)進(jìn)行并行計(jì)算。王啟聰[2]提出了基于GPU的空譜聯(lián)合核稀疏表示高光譜分類并行優(yōu)化算法。Santos L[3]等人基于高并行GPU框架實(shí)現(xiàn)了機(jī)載有損高光譜圖像壓縮算法，得到了顯著的加速效果。不過，GPU技術(shù)對(duì)于數(shù)據(jù)密集型的計(jì)算效率卻不高，而且代碼編寫調(diào)試復(fù)雜。

遙感圖像分類算法在云計(jì)算平臺(tái)上的應(yīng)用多處在研究階段。Li Y[4]等人設(shè)計(jì)了Spark[5-6]平臺(tái)下基于主成份分析的高光譜解混方法，在保證精度的同時(shí)，取得了很好的性能提升。Mohamed H.Almeer[7]在Hadoop[8]云平臺(tái)下實(shí)現(xiàn)了高容量遙感數(shù)字圖像的索貝爾濾波、銳化、對(duì)比度增強(qiáng)等算法。常生鵬[9]等人利用Hadoop技術(shù)設(shè)計(jì)改進(jìn)了Meanshift圖像邊緣分割算法，實(shí)驗(yàn)結(jié)果表明，在Hadoop環(huán)境下的高分辨率衛(wèi)星圖像數(shù)據(jù)處理速度有了明顯的改善。

文中將云計(jì)算引入到高光譜遙感圖像分類中，在Spark平臺(tái)上實(shí)現(xiàn)了基于空間相關(guān)性正則化的稀疏表示分類并行化方法 （Parallel Spatial Correlationregularized Sparse Representation Classification，PSCSRC），進(jìn)行遙感圖像的分類，并通過實(shí)驗(yàn)驗(yàn)證了Spark處理高光譜遙感圖像的可行性和高效性。

1 PSCSRC算法設(shè)計(jì)

1.1 Spark云計(jì)算平臺(tái)介紹

Spark是一個(gè)快速、易用的分布式大規(guī)模數(shù)據(jù)計(jì)算框架，其核心是彈性分布式數(shù)據(jù)集（Resilient Distributed Datasets，RDD），RDD 是對(duì)集群上并行處理數(shù)據(jù)的分布式內(nèi)存抽象，通過在集群中的多臺(tái)機(jī)器上進(jìn)行數(shù)據(jù)分區(qū)實(shí)現(xiàn)并行計(jì)算。RDD有轉(zhuǎn)換算子和動(dòng)作算子兩類操作，通過轉(zhuǎn)換算子生成新RDD，形成階段（Stage），最終通過動(dòng)作（Action）提交 Spark程序。

1.2 稀疏表示分類方法簡介

分類是高光譜圖像處理中最重要的內(nèi)容之一。隨著壓縮感知和稀疏表示理論的發(fā)展，越來越多的學(xué)者提出采用稀疏表示進(jìn)行分類[10]，其假設(shè)條件是測(cè)試樣本能夠由一些訓(xùn)練樣本所組成的結(jié)構(gòu)化字典進(jìn)行線性組合稀疏表示。該方法首先通過求解稀疏性約束的最優(yōu)化重構(gòu)模型，得到字典原子權(quán)值組成的稀疏向量，然后根據(jù)恢復(fù)稀疏向量的特征，確定測(cè)試樣本的類別。

1.3 基于空間相關(guān)性正則化的稀疏表示分類方法

基于空間相關(guān)性正則化的稀疏表示分類方法[11]（Spatial Correlation regularized Sparse Representation Classification，SCSRC）是一種基于空間相關(guān)性正則化稀疏表示的迭代型分類算法，該方法首先建立高光譜圖像分類的稀疏表示模型，通過添加像元間的光滑性正則項(xiàng)和訓(xùn)練樣本的空間位置信息，用以提高分類精度，用ADMM[12]求解模型，最終進(jìn)行分類。

1.4PSCSRC算法

算法分為數(shù)據(jù)準(zhǔn)備階段和計(jì)算階段，數(shù)據(jù)準(zhǔn)備階段的主要工作是讀入高光譜圖像數(shù)據(jù)，再分塊存儲(chǔ)在HDFS中。計(jì)算階段的主要工作是從HDFS讀取矩陣文件，迭代計(jì)算稀疏矩陣和若干中間矩陣，若達(dá)到指定迭代次數(shù)，則迭代結(jié)束，使用稀疏矩陣對(duì)測(cè)試樣本進(jìn)行分類，否則，使用稀疏矩陣更新中間矩陣，繼續(xù)迭代。

1.4.1 數(shù)據(jù)準(zhǔn)備階段

數(shù)據(jù)準(zhǔn)備階段就是一個(gè)數(shù)據(jù)預(yù)處理的過程。首先將高光譜遙感圖像讀入內(nèi)存以三維數(shù)組形式存放，為了處理方便將其轉(zhuǎn)化為二維數(shù)組X并歸一化，其中 X=[X1，X2，…，XN]，Xi=[x1，x2，…xL]T，Xi表示一個(gè)像元在各個(gè)波段上的光譜響應(yīng)，N為總樣本數(shù)量，L為波段數(shù)。從X中隨機(jī)取得若干列向量作為訓(xùn)練樣本存入 A 中，A=[A1，A2，…，AM]，其中 Ai=[a1，a2，…aL]T，M為訓(xùn)練樣本數(shù)量。利用A和X進(jìn)一步計(jì)算，得到F和B，預(yù)處理過程的偽代碼如下所示：

預(yù)處理之后得到了矩陣F和矩陣B。B=[B1，B2，…，BM]T，其中 Bi=[b1，b2，…bN]，將矩陣 B 使用 HDFS的SequenceFile以鍵值對(duì)形式分塊存儲(chǔ)在集群中，其中 i是矩陣行號(hào)，Bi是行向量。 F=[F1，F(xiàn)2，…，F(xiàn)M]，其中 Fi=[f1，f2，…fM]T，直接使用文本文件保存在本地即可。

1.4.2 計(jì)算階段

根據(jù)RDD提供的算子以及自定義的函數(shù)設(shè)計(jì)PSCSRC算法，其具體流程如下所示。

1）PSCSRC算法首先將矩陣F作為廣播變量F_broadcast發(fā)送到各個(gè)Worker節(jié)點(diǎn)，再利用RDD提供的SequenceFile函數(shù)讀取矩陣B并轉(zhuǎn)換成RDD數(shù)據(jù)模型B_RDD，同時(shí)使用partition參數(shù)設(shè)置分塊數(shù)目，記分塊數(shù)量為partitions。接著filterWith算子將B_RDD每個(gè)分塊轉(zhuǎn)換成新的RDD，再分別調(diào)用count函數(shù)得到行向量的數(shù)量，由于矩陣B在HDFS上按照行號(hào)順序存儲(chǔ)，因此可以得到B_RDD中每個(gè)分塊中行向量的起始行號(hào)和終止行號(hào)

2）接下來對(duì)B_RDD使用mapValues算子創(chuàng)建五元祖group_RDD，分塊中元素鍵值對(duì)具體形式為>，.其中 P=[P1，P2，…，PM]T，Pi=[p1，p2，…pN]，Q=[Q1，Q2，…，QM]T，Qi=[q1，q2，…qN]，R=[R1，R2，…，RM]T，Ri=[r1，r2，…rN]，T=[T1，T2，…，TM]T，Ti=[t1，t2，…tN]，使用這種五元組的形式存儲(chǔ)中間矩陣不僅便于接下來對(duì)其各分塊的更新，而且減少了集群中數(shù)據(jù)的流動(dòng)。將group_RDD緩存到內(nèi)存中以備下輪迭代使用。

3）group_RDD使用mapValues算子將分區(qū)中的五元組相加得到plus_RDD。接著，plus_RDD使用map Partitions算子，通過外積法[13]與廣播變量F_broadcast相乘得到中間矩陣，然后使用reduceByKey算子將中間矩陣按照key累加得到S_RDD。

4）使用collect函數(shù)將S_RDD匯總到Driver端并稀疏化得到SArray，如果已達(dá)到指定迭代次數(shù)，轉(zhuǎn)5，否則，由于P、Q、R、T各分塊與B的各分塊結(jié)構(gòu)相同，所以根據(jù)P、Q、R、T中的各個(gè)分塊的起始行號(hào)和終止行號(hào)，將SArray中對(duì)應(yīng)的子矩陣發(fā)送到各分塊所在的Worker節(jié)點(diǎn)更新分塊，其中對(duì)P的分塊使用Gauss-Seidel[14，15，16]方法更新，對(duì) Q 的分塊使用軟閾值方法更新，對(duì)R、T的分塊使用線性方法更新，更新完成后轉(zhuǎn)3繼續(xù)迭代。這種首先計(jì)算各中間矩陣每個(gè)分塊對(duì)應(yīng)的子矩陣，然后再把子矩陣分發(fā)到各分塊所在節(jié)點(diǎn)的方法，相比于直接將稀疏矩陣分發(fā)到各節(jié)點(diǎn)之后再提取出子矩陣更新分塊的方法來說，將集群中移動(dòng)的數(shù)據(jù)量降低為1/partitions。

5）復(fù)用SCSRC算法中根據(jù)稀疏矩陣SArray計(jì)算OA和Kappa系數(shù)的部分，總精度（OA）反映了一個(gè)隨機(jī)樣本的分類結(jié)果與真實(shí)標(biāo)記類型相一致的概率，其可以表述為如下形式：

為了更能體現(xiàn)圖像分類的整體誤差，使用Kappa系數(shù)k來評(píng)價(jià)分類精度。Kappa系數(shù)計(jì)算公式表述為如下形式：

這里C表示類別數(shù)量，M為一個(gè)C×C的混淆矩陣，N表示測(cè)試樣本數(shù)量。

2 實(shí)驗(yàn)結(jié)果及分析

本實(shí)驗(yàn)使用機(jī)載可見紅外成像光譜儀采集的美國印第安納州Indian Pines實(shí)驗(yàn)區(qū)的高光譜遙感圖像。該圖像共包含220個(gè)波段，空間分辨率為20米，圖像大小為145×145（數(shù)據(jù)大小為6 MB）。并使用ENVI軟件的Mosaicking功能將此數(shù)據(jù)上下鑲嵌為29 000×145 （數(shù)據(jù)大小為 1.21 GB）、116 000×145（4.83 GB）的數(shù)據(jù)。

實(shí)驗(yàn)運(yùn)行在由9臺(tái)服務(wù)器組成的Spark集群上，各節(jié)點(diǎn)配置為Intel Xeon E7-4807，6核，內(nèi)存24 GB，Spark版本為1.4.1，操作系統(tǒng)為Ubuntu12.04。

相同的執(zhí)行環(huán)境下基于CPU的SCSRC程序運(yùn)行結(jié)果（分片數(shù)為1）和基于Spark的PSCSRC并行程序運(yùn)行結(jié)果如表1所示，結(jié)果為實(shí)驗(yàn)重復(fù)執(zhí)行10次后的平均值。圖1給出了兩種方法的視覺比較，其中（a）是各類地物的真實(shí)分布，（b）是SCSRC分類結(jié)果，（c）是 PSCSRC 分類結(jié)果。

表1 Spark下不同分片數(shù)處理結(jié)果

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

由表1可知，PSCSRC算法在保證分類精度的同時(shí)，隨著分塊數(shù)量的增加執(zhí)行時(shí)間不斷降低，只是加速比效果不太理想。并行算法下能保證分類精度是因?yàn)镾CSRC算法分類時(shí)使用極大似然估計(jì)的思想使用稀疏矩陣來分類，雖然Spark平臺(tái)下得到的稀疏矩陣與串行平臺(tái)下有微小的差別，不過并不影響分類效果。加速比效果不太理想是因?yàn)镾park會(huì)為每個(gè)分塊啟動(dòng)一個(gè)線程來處理該分塊的數(shù)據(jù)，而本實(shí)驗(yàn)中每個(gè)分塊數(shù)據(jù)量太小，創(chuàng)建線程的系統(tǒng)開銷相對(duì)分片數(shù)據(jù)的計(jì)算時(shí)間占較大比重，因此加速效果不明顯。

由圖1的b和c可以看出，PSCSRC算法的分類結(jié)果與SCSRC的分類結(jié)果完全相同。

分別使用 29 000×145、116 000×145 的數(shù)據(jù)作為實(shí)驗(yàn)數(shù)據(jù)集重復(fù)以上實(shí)驗(yàn)，實(shí)驗(yàn)結(jié)果如圖2所示。

圖2 加速比性能測(cè)試結(jié)果

從圖2可知，當(dāng)數(shù)據(jù)量變大時(shí)，加速效果越來越明顯，數(shù)據(jù)量越大，實(shí)際的加速比越接近線性加速比，并且分塊數(shù)量越大，加速效果越好。原因有兩個(gè)：1）算法中盡量避免了會(huì)導(dǎo)致shuffle操作的RDD轉(zhuǎn)換算子，大大減少了集群中數(shù)據(jù)的跨節(jié)點(diǎn)流動(dòng)；2）每次迭代過程中使用稀疏矩陣SArray更新P、Q、R、T時(shí)，只是將各分塊對(duì)應(yīng)的子矩陣分發(fā)到各個(gè)節(jié)點(diǎn)，而不是整個(gè)稀疏矩陣，同樣是減少了數(shù)據(jù)在集群中的流動(dòng)。

3 結(jié) 論

SCSRC中涉及到多個(gè)矩陣相加再相乘的部分可能由于單機(jī)尺度的限制而無法運(yùn)行，因此本文提出了基于Spark的PSCSRC算法，將多個(gè)矩陣以行向量的形式存儲(chǔ)至五元組，再分發(fā)到集群中各個(gè)Worker節(jié)點(diǎn)。結(jié)合高速緩沖存儲(chǔ)器Cache的特性，采用利于分布式下矩陣并行計(jì)算的外積法來提升矩陣相乘效率，其中對(duì)于大矩陣與小矩陣相乘的情況，通過采用廣播變量將小矩陣分發(fā)到各Worker節(jié)點(diǎn)進(jìn)行快速計(jì)算。利用Spark中數(shù)據(jù)本地性的任務(wù)分配策略，將迭代中常用的RDD進(jìn)行緩存，利用Spark內(nèi)存計(jì)算的特點(diǎn)進(jìn)行快速迭代計(jì)算。最后，通過分塊中行向量的起始行號(hào)和終止行號(hào)，使用大矩陣的對(duì)應(yīng)子矩陣更新4個(gè)中間矩陣，而不是將整個(gè)大矩陣分發(fā)到各Worker節(jié)點(diǎn)，大大減少了節(jié)點(diǎn)間的流量傳輸，加快程序的執(zhí)行速度。

[1]姚禹，向晶.全球在軌衛(wèi)星數(shù)量突破1000顆大關(guān)[J].中國無線電，2012（11）:77-77.

[2]王啟聰，吳澤彬，劉建軍，等.基于GPU的空譜聯(lián)合核稀疏表示高光譜分類并行優(yōu)化[J].計(jì)算機(jī)工程與科學(xué)， 2014，36（12）:2321-2330.

[3]Santos L， Magli E， Vitulli R，et al.Highlyparallel GPU architecture for lossy hyperspectral image compression[J].IEEE Journal of Selected Topics in Applied Earth Observations&Remote Sensing， 2013，6（2）:670-681.

[4]Li Y， Wu Z， Wei J， et al.Fast principal componentanalysis for hyperspectralimaging based on cloud computing[C]//IEEE International Geoscience and Remote Sensing Symposium.IEEE，2015:513-516.

[5]HOLDEN KARAU.Spark快速數(shù)據(jù)處理[M].余璜，譯.北京：機(jī)械工業(yè)出版社，2014.

[6]夏俊鸞.Spark大數(shù)據(jù)處理技術(shù)[M].北京：電子工業(yè)出版社，2015.

[7]Almeer M H.Cloud hadoop map reduce for remote sensing image analysis[J].Journal of Emerging Trends in Computing& Information Sciences，2012，3（4）:637-644.

[8]蔡斌.Hadoop技術(shù)內(nèi)幕[M].北京：機(jī)械工業(yè)出版社，2013.

[9]常生鵬，馬億旿，蔡立軍，等.一種基于Hadoop的高分辨率遙感圖像處理方法[J].計(jì)算機(jī)工程與應(yīng)用，2015，51（11）:167-171.

[10]Chen Y，Nasrabadi N M，Tran T D.Hyperspectral image classification using dictionary-based sparse representation[J].IEEE Transactions on Geoscience&Remote Sensing，2011，49（10）:3973-3985.

[11]劉建軍.高光譜圖像非負(fù)矩陣分解解混與稀疏表示分類研究[D].南京：南京理工大學(xué)，2014.

[12]Yang J，Zhang Y.Alternating direction algorithms for l1-problems in compressive sensing[J].Siam Journal on Scientific Computing， 2011，33 （1）:250-278.

[13]孫遠(yuǎn)帥，陳垚，官新均，等.基于Hadoop的大矩陣乘法處理方法 [J].計(jì)算機(jī)應(yīng)用，2013，33（12）:3339-3344.

[14]Li S， Zhang B， Li A， et al.Hyperspectral imagery clustering with neighborhood constraints[J].IEEE Geoscience&Remote Sensing Letters，2013， 10（3）:588-592.

[15]Rakotomamonjy A.Surveying and comparing simultaneous sparse approximation （or group-lasso） algorithms[J].Signal Processing， 2011， 91（7）:1505-1526.

[16]Iordache M D，Bioucas-Dias J M，Plaza A.Total variation spatial regularization for sparse hyperspectralunmixing [J].IEEE Transactions on Geoscience&Remote Sensing， 2012， 50 （11）:4484-4502.

A kind of hyperspectral image classification parallel method based on Spark

LIU Zhen，ZHU Yao-qin
（College of Computer Science and Engineering，Nanjing University of Science and Technology，Nanjing 210094，China）

In order to achieve the classification of large amounts of hyperspectral image，proposed a hyperspectral sparse representation classification parallel method on Spark，designing a 5-tuple intermediate data storage structure and then distributing sub-matrix of spare matrix corresponded to each block of intermediate matrices to each node in each iteration，which reduce time consumption caused by intermediate data Shuffle and data transmission.Experiment results show that execution speed of sparse representation classification parallel method under Spark platform has been significantly improved.

hyperspectral； cloud computing； spark platform； classification； sparse representation

TN919

A

1674－6236（2017）12-0019-04

2016-10-18稿件編號(hào)：201610087

國家重點(diǎn)實(shí)驗(yàn)室開放研究基金，復(fù)雜產(chǎn)品智能制造系統(tǒng)技術(shù)國家重點(diǎn)實(shí)驗(yàn)室開放基金資助（QYYE1603）

劉 震（1990—），男，河南睢縣人，碩士研究生。研究方向：云計(jì)算。