唐艷慧，趙 鵬，王承琨

(東北林業(yè)大學(xué) 信息與計(jì)算機(jī)工程學(xué)院，黑龍江 哈爾濱 150040)

1 引 言

精準(zhǔn)快速地對(duì)木材進(jìn)行分類在工業(yè)生產(chǎn)中具有重大意義。使用計(jì)算機(jī)智能手段對(duì)木材進(jìn)行分類不僅可以減少工作人員的工作量，同時(shí)也可以彌補(bǔ)經(jīng)驗(yàn)和技術(shù)的不足減少木材分類出錯(cuò)率。目前使用計(jì)算機(jī)手段對(duì)木材分類的研究主要集中在數(shù)字圖像處理上。

目前木材樹種分類識(shí)別主要可以分為宏觀和微觀兩種分類方法，利用宏觀方法對(duì)木材樹種進(jìn)行分類主要是指使用木材表面的紋理結(jié)構(gòu)、木材的顏色等特征標(biāo)記木材種類。但是由于木材生長(zhǎng)環(huán)境復(fù)雜極易造成顏色、紋理等特征被篡改或破壞，使得木材錯(cuò)分率增加。另一種方法是使用木材的微觀結(jié)構(gòu)對(duì)木材樹種進(jìn)行識(shí)別分類，該類方法具有極高的分類正確率。近幾年來，許多學(xué)者利用木材的微觀紋理結(jié)構(gòu)對(duì)木材進(jìn)行分類并取得可喜的科研成果。印度學(xué)者Yadava[1-2]等人利用DWT變換和LBP(局部二值變換)對(duì)木材的微觀數(shù)字圖像進(jìn)行了識(shí)別，正確率均在90%以上。希臘學(xué)者Barmpoutis[3]將木材橫截面的二維數(shù)字圖像看成是一種多維信號(hào)的集合，多圖像進(jìn)行分解，利用多方向的紋理信息對(duì)木材進(jìn)行分類，正確率可接近91%。但是該類方法最大的問題在于無法做到無損檢測(cè)，而且往往需要制作切片，識(shí)別效率低下。本文使用的是小倍率放大，該方法可以不用破壞木材原有的形態(tài)結(jié)構(gòu)，同時(shí)也可以體現(xiàn)出木材的一定微觀結(jié)構(gòu)。

隨著圖像獲取設(shè)備的迭代更新，高光譜圖像成像儀成為了一種全新獲取圖像信息的手段。高光譜圖像具有“圖譜合一”的特性，在獲取物體圖像信息的同時(shí)也可以獲取該物體的光譜信息。目前便攜式小型高光譜成像儀走進(jìn)實(shí)驗(yàn)室，該類光譜儀小型輕便可以和數(shù)字放大鏡共同工作，方便采集木材微觀結(jié)構(gòu)的高光譜圖像。

本文借助于美國(guó)SOC公司生產(chǎn)的SOC710VP高光譜成像儀對(duì)8種不同生長(zhǎng)區(qū)域的木材進(jìn)行了樹種分類識(shí)別。首先使用波段選擇方法選擇了高光譜圖像中最有代表性的若干個(gè)波段，然后利用分形思想提出了這些波段圖像的多重分形曲線，最后融合了這些多重分形曲線找到了整個(gè)樣本的特征向量。實(shí)驗(yàn)證明本文方法具有良好的分類效果。

2 材料與方法

2.1 材料

本實(shí)驗(yàn)的所有數(shù)據(jù)采集使用了美國(guó)SOC公司生產(chǎn)的SOC710VP高光譜成像儀，該成像儀的光譜范圍為372.53～1 038.57 nm，共計(jì)128個(gè)波段，光譜分辨率為5.2 nm。拍攝到的顯微高光譜分辨率為200×200×128。

表1 樣本木材的詳細(xì)資料Tab.1 Details of the sample wood

續(xù) 表

本文以8種木材橫截面的高光譜圖像為研究對(duì)象，具體木材樣本資料見表1。每種木材采集120個(gè)樣本，這些樣本均來自不同的樹木，為方便采集，首先需要對(duì)樣本進(jìn)行處理，處理方法如下：

首先需要將木材切割成若干個(gè)長(zhǎng)與寬均為20 mm的條狀樣本，這些條狀樣本的側(cè)面是木材的弦切面，頂面和底面是木材的橫切面。然后使用圓盤鋸將這個(gè)條狀樣本切割成若干高度為30 mm的小長(zhǎng)方體，切割時(shí)注意其橫切面平滑無毛刺，最后制作的樣本長(zhǎng)寬高分別為20 mm×20 mm×30 mm的小長(zhǎng)方體。

每種木材樣本的前96個(gè)小長(zhǎng)方體作為訓(xùn)練集，從剩下的24個(gè)小長(zhǎng)方體作為測(cè)試集。將這些樣本在SOC710VP高光譜成像儀下完成高光譜圖像采集,放大倍數(shù)為45×，圖1為上述8種樣本所采集的橫截面圖像。

圖1 8種木材樣本的橫截面Fig.1 Cross-section of eight wood samples

2.2 高光譜圖像的波段選擇

高光譜圖像可以看作是一個(gè)三維數(shù)據(jù)立方體。該三維立方體可以看成是由多張圖片累積而成，每張圖片的長(zhǎng)與寬對(duì)應(yīng)高光譜圖像的長(zhǎng)M與寬N，圖片的數(shù)量多少對(duì)應(yīng)高光譜圖像的高L，圖2中給出了高光譜圖像的示意圖。

圖2 高光譜三維立方體示意圖Fig.2 Hyperspectral three-dimensional cube diagram

高光譜圖像往往具有幾十至上百個(gè)波段，如果直接對(duì)這些波段進(jìn)行處理，不僅處理速度慢也會(huì)影響分類精度的大小。所以在對(duì)高光譜圖像進(jìn)行處理時(shí)往往需要降維處理。目前高光譜圖像的降維主要包含兩類方法，一類是將高光譜圖像做某種變換進(jìn)行特征提取的方法，即尋找一個(gè)高光譜圖像的低維投影，這類方法一般會(huì)改變高光譜圖像的整體結(jié)構(gòu)，降維后的高光譜圖像與任意波段均不相同。代表方法有主成分分析(principal component analysis,PCA)、線性判別分析(linear discriminant analysis,LDA)、最小噪聲分?jǐn)?shù)(minimum noise fraction rotation,MNF Rotation)等。另一類方法是利用某種手段從現(xiàn)有的高光譜圖像中挑選出若干波段的波段選擇方法，這類算法挑選的結(jié)果是整幅高光譜圖像的某一子集。代表方法有自適應(yīng)波段選擇(adaptive band selection ABS)、基于K-L的波段選擇(Kullback-Leibler divergence,K-L)。本文主要側(cè)重高光譜圖像紋理分析與研究，故采用波段選擇方法對(duì)高光譜圖像進(jìn)行降維。

2.2.1 自適應(yīng)波段選擇(Adaptive Band Selection)

自適應(yīng)波段選擇主要是通過波段信息量和波段與波段之間的相關(guān)性進(jìn)行波段選擇，具體描述如下[4]：

首先我們定義高光譜圖像為A={b1,b2,b3…bi,…b128}，其中bi代表高光譜圖像的第i個(gè)波段。則第i個(gè)波段與第j個(gè)波段的相關(guān)系數(shù)定義為：

(1)

(2)

定義第i個(gè)波段的索引指數(shù)Ii為：

(3)

其中：上式中的Ri,j-1與Ri,j+1分別代表第i個(gè)波段與其前后兩個(gè)波段的相關(guān)系數(shù)。值得注意的是這里的i,j∈[2,127]且i,j∈N*,最后選擇的波段不包含第一個(gè)波段與最后一個(gè)波段。

Ii充分考慮了每幅圖像的信息富集程度與相鄰波段的相似性，Ii越大，相應(yīng)波段的信息量就越大，也就是說該波段在整幅圖像中越具有代表性。假設(shè)選擇出的波段數(shù)量為n,則選擇Ii較大的前n個(gè)波段即可完成波段選擇[5-6]。

2.2.2 K-L散度(Kullback-Leibler divergence)

K-L散度是一種在信息論中得到廣泛應(yīng)用的信號(hào)相似度的度量方法，其定義為如下：設(shè)兩個(gè)離散的概率分布函數(shù)分別為P=[P1,…,Pn,…PN]T和Q=[Q1,…,Qn,…QN]T，則Q相對(duì)于P的K-L散度為[5-6]：

(4)

上式的物理意義是使用Q中的元素來表示P的所有元素所需要的額外信息量的大小，其單位與log的底數(shù)有關(guān),當(dāng)?shù)讛?shù)為自然常數(shù)時(shí)，其單位為nat；當(dāng)?shù)讛?shù)為2時(shí)，其單位為bit。K-L散度越大,說明使用Q來表示P越困難[7]。

(5)

以此計(jì)算高光譜中任意兩個(gè)波段的散度，得到矩陣對(duì)稱矩陣D，設(shè)D中的第i行為ai=[D(i,1)…D(i,L)]，ai中元素的平均值代表將第i個(gè)波段從數(shù)據(jù)集中刪除后，所帶來的損失。找到元素平均值最小的ai，去除該波段后重復(fù)上述算法流程再次去除波段，就可以得到含有信息量最大的n個(gè)波段[8]。

2.3 基于多重分形的高光譜圖像紋理分類方法

分形和不規(guī)則的形狀有關(guān)，是對(duì)圖像自相似性的描述，由法國(guó)數(shù)學(xué)家Mandelbrot[9]于1975年提出。通過對(duì)自然景物紋理圖像的研究進(jìn)展，人們發(fā)現(xiàn)分形維數(shù)與圖像粗糙度具有十分緊密的聯(lián)系，作為一種非線性理論，分形理論可以很好的應(yīng)用于研究不同木材的紋理特征。

其中計(jì)算分形維數(shù)的方法有很多，常用的有盒分法(box-counting)、隨機(jī)游走法(Fractional Brownian motion)以及頻域法[10-13]。針對(duì)不同的研究對(duì)象選擇相應(yīng)的分形維數(shù)作為度量，本文使用盒分法來計(jì)算分形維數(shù)其具體定義如下：

設(shè)F是Rn上任意非空的有界子集，Nδ(F)是直徑最大為δ，可以覆蓋F集的最少個(gè)數(shù)，則F的盒維數(shù)記為：

(6)

上述的盒分法可以計(jì)算一幅數(shù)字圖像的分形維數(shù)。但是對(duì)于一幅數(shù)字圖像來說，如果只使用一個(gè)分形維數(shù)對(duì)其進(jìn)行描述，顯然不能完全刻畫其結(jié)構(gòu)和特點(diǎn)。為此我們將一幅圖像分成若干個(gè)小尺度，分別計(jì)算其分形維數(shù)，這樣就可以獲得多個(gè)分形維數(shù)，為更好的實(shí)現(xiàn)這一想法，我們引入圖像密度函數(shù)。

圖像密度函數(shù)表示了在測(cè)量點(diǎn)附近區(qū)域的一種不均勻性，其具體表達(dá)式如下：

(7)

其中：X(x,y)∈R2代表二維數(shù)字圖像,B(X(x,y),r)代表一個(gè)中心為X(x,y)半徑為r的正方形區(qū)域，f(B(X(x,y),r))代表對(duì)這個(gè)區(qū)域的一種運(yùn)算，具體包含了以下3種方法[14]：

(1)基于圖像亮度(強(qiáng)度)的函數(shù)運(yùn)算：

f(B(X(x,y),r))=
?B(X(x,y),r)(Gr*I)dX(x,y).

(8)

(2)基于圖像梯度的函數(shù)運(yùn)算,其中fk代表橫向梯度與縱向梯度：

(9)

(3)基于拉普拉斯算子的函數(shù)運(yùn)算：

f(B(X(x,y),r))=

?B(X(x,y),r)|2(Gr*I)|dX(x,y).

(10)

以上3個(gè)式子中的“*”運(yùn)算代表二維卷積運(yùn)算，Gr代表高斯平滑核，其具體表達(dá)式如下所示,其作用是對(duì)圖像起平滑作用。

(11)

圖3是美國(guó)紅橡的密度圖像(r=4,δ=2)其中圖3(a)利用圖像梯度函數(shù)得到的密度圖像，圖3(b)是利用亮度函數(shù)得到的密度圖像，圖3(c)是利用拉普拉斯算子函數(shù)得到的密度圖像。

圖3 密度利用圖像梯度函數(shù)(a)、亮度函數(shù)(b)和拉普拉斯算子函數(shù)得到的(c)圖像函數(shù)。Fig.3 Density images by using image gradient function(a),luminance function(b),and Laplace operator functicn(c),respectively.

根據(jù)圖3可知密度函數(shù)圖像的像素取值范圍主要集中在0～4之間。設(shè)步長(zhǎng)為s=0.1,則整個(gè)密度函數(shù)圖像的像素取值可被分為m個(gè)小區(qū)間如：(0,0.1],(0.1,0.2]….(0.1*(i-1),0.1*i]…(0.1*(m-1),0.1*m],整個(gè)圖像可按照上述區(qū)間范圍被分解成m個(gè)二值圖像，其具體方法按式(12)計(jì)算：

(12)

我們使用計(jì)盒維數(shù)來估計(jì)圖像的分形維數(shù)，計(jì)盒維數(shù)的算法如下:

首先設(shè)一個(gè)長(zhǎng)度為r(r∈{2,4,8,…,2n}∧2n

分別對(duì)這m個(gè)二值圖像使用盒分法求解計(jì)盒維數(shù)，以小區(qū)間(0.1*(i-1),0.1*i]的右端點(diǎn)為橫坐標(biāo)α以該區(qū)間所對(duì)應(yīng)二值圖像的計(jì)盒維數(shù)f(α)=k作為縱坐標(biāo)可以建立多重分形譜曲線。圖4是根據(jù)圖3采用不同密度圖像所得到的多重分形譜曲線。

圖4 3種不同密度圖像所對(duì)應(yīng)的多重分形譜曲線Fig.4 Multifractal spectral curves corresponding to three different density images

利用該方法可以對(duì)普通的數(shù)字圖像進(jìn)行分類，但通過實(shí)驗(yàn)發(fā)現(xiàn)，有多處重合并不能在木材微觀結(jié)構(gòu)的數(shù)字圖像上得到理想的分類效果，這一點(diǎn)會(huì)在第三節(jié)中給出結(jié)論。

高光譜圖像往往擁有百余個(gè)波段，并不能像正常上述的灰度圖像那樣求解其密度圖像，為了能夠提高效率并盡可能多地使用高光譜圖像中的有用信息，必須使用降維方法對(duì)高光譜圖像的多個(gè)波段進(jìn)行選擇，對(duì)所選波段求解多重分形譜。

假設(shè)通過2.3節(jié)的波段選擇方法得到的波段為{x1,x2,…xn}，使用上述方法求解每一個(gè)波段xi所對(duì)應(yīng)的多重分形曲線yi(i∈[1,n])，這樣就可以得到n個(gè)多重分形曲線，圖5(a)中給出了美國(guó)紅橡木材的一個(gè)樣本在n=10時(shí)使用梯度密度函數(shù)圖像所得到的10條多重分形曲線。

圖5 高光譜圖像10個(gè)波段的多重分形譜Fig.5 Multifractal spectra of 10 bands in hyperspectral images

3 結(jié)果與討論

在本節(jié)我們將與傳統(tǒng)常用的紋理分類方法的準(zhǔn)確率和運(yùn)行時(shí)間作對(duì)比，本文與應(yīng)用較多的傳統(tǒng)方法灰度共生矩陣(GLCM)做對(duì)比實(shí)驗(yàn)，該方法的基本思想是使用PCA對(duì)高光譜圖像進(jìn)行降維，提取高光譜圖像的第一主成分，使用第一主成分對(duì)應(yīng)的圖像利用灰度共生矩陣(GLCM)提取特征值，最后送入分類器進(jìn)行分類。與多重分形做紋理對(duì)比試驗(yàn)，分類精度提升顯著。

3.1 不同算法的時(shí)間對(duì)比

本文實(shí)驗(yàn)所使用的計(jì)算機(jī)配置為CPU：Intel I7-6700，內(nèi)存：8 G，顯卡：AMD Radeon R7 200 硬盤容量為1 TB。使用上述配置得到的各個(gè)算法的運(yùn)行時(shí)間如表2所示。

表2 不同算法在各個(gè)階段的運(yùn)行時(shí)間Tab.2 Running time of different algorithms in each stage

從表2中可以看出主要影響運(yùn)行時(shí)間的因素包含降維方法與函數(shù)密度計(jì)算，顯然ABS降維方法的效率要遠(yuǎn)遠(yuǎn)高于K-L散度。函數(shù)密度中的亮度運(yùn)行速度最快，拉普拉斯最慢。下面我們來看一下使用各種算法處理一個(gè)實(shí)驗(yàn)樣本所需的時(shí)間。表3中給出了提取的波段數(shù)量為10時(shí)使用不同方法處理一個(gè)實(shí)驗(yàn)樣本所需的時(shí)間。

表3 處理一個(gè)樣本所需時(shí)間Tab.3 The time required to process a sample

從表3中可以看出，使用ABS作為降維方法和亮度密度函數(shù)求解分形曲線，速度最快，使用K-L散度作為降維方法和拉普拉斯密度函數(shù)求解分形曲線速度最慢。

3.2 使用PCA算法求解高光譜圖像的紋理分類結(jié)果

我們使用PCA對(duì)高光譜圖像進(jìn)行降維，并將其第一主成分作為圖像紋理分類對(duì)象。表4中給出了在SVM分類器下的混合分類正確率。

表4 使用PCA降維方法的分類正確率
Tab.4 Classification accuracy using PCA dimensionality reduction method

方法訓(xùn)練集交叉正確率/%測(cè)試集正確率/%亮度密度函數(shù)69.0179.16梯度密度函數(shù)69.0180.20拉普拉斯密度函數(shù)73.3078.12GLCM67.6072.91

從表4中可以看出上述4種紋理分析法在PCA降維方法下的分類正確率略低，一方面原因在于使用PCA降維方法對(duì)高光譜圖像進(jìn)行降維時(shí)，高光譜的大量信息被壓縮，這就導(dǎo)致正確率的下降。另一方面該方法對(duì)圖像的尺寸大小具有嚴(yán)格要求，不同木材的橫截面微觀結(jié)構(gòu)也具有一定的相似度。顯然使用PCA提取高光譜的第一主成分這種方法在木材高光譜圖像中并不十分適合。

3.3 基于波段選擇的多重分形譜紋理分類結(jié)果

圖6是美國(guó)紅橡樹種使用亮度密度函數(shù)提取的波段選擇數(shù)量n=15，8，2時(shí)的多重分形曲線。

圖6 不同波段選擇數(shù)量下的多重分形曲線。(a)n=15;(b)n=8;(c)n=2。Fig.6 Multifractal curves with different number of bands selected.(a)n=15;(b)n=8;(c)n=2.

通過圖6可以看出，波段選擇數(shù)量越大分形曲線越平滑，每條多重分形曲線的類內(nèi)距離更近，這樣在分類時(shí)就能得到更高的分類精度，這一結(jié)論可以在后面的分類正確率上得到證實(shí)。

為了能夠更加準(zhǔn)確地說明分類正確率，本文使用SVM與BP神經(jīng)網(wǎng)絡(luò)兩種分類器對(duì)其進(jìn)行分類。由于BP神經(jīng)網(wǎng)絡(luò)每一次的分類正確率并不穩(wěn)定，所以本文的分類精度是在訓(xùn)練25次后所得到的平均分類精度。

圖7、圖8中給出了使用K-L散度和ABS作為波段選擇方法的SVM分類精度與BP混合分類精度，可以看出使用梯度作為圖像密度的K-L降維方法分類精度最高，在SVM下最高分類精度可達(dá)97.91%，在BP下最高分類精度可達(dá)90.83%。另外亮度和拉普拉斯方法均滿足正確率先上升后下降的趨勢(shì)，梯度方法隨著波段選擇數(shù)量增加而增加。造成這一因素的主要原因在于過多的波段選擇數(shù)量會(huì)將一部分噪聲融入平均分形曲線中從而造成正確率先升后降的現(xiàn)象。從整體上看在SVM分類器下梯度方法正確率高于其他兩種方法，在BP分類器下亮度方法正確率高于其它兩種方法。

圖7 基于K-L降維方法的分類精度。(a)SVM;(b)BP。Fig.7 Classification accuracy based on k-l dimension reduction method.(a)SVM;(b)BP.

圖8 基于ABS降維方法的分類精度。(a)SVM;(b)BP。Fig.8 Classification accuracy based on ABS dimension reduction method.(a)SVM;(b)BP.

表5中給出了不同方法下的最佳分類正確率所對(duì)應(yīng)的波段選擇數(shù)量。將表5和表4進(jìn)行對(duì)比可發(fā)現(xiàn)，使用K-L波段選則方法和多重分形曲線對(duì)對(duì)木材的紋理進(jìn)行提取可擁有更高的分類正確率。分類準(zhǔn)確率最高可達(dá)到97.91%。

表5 不同方法下的最佳分類正確率與波段選擇數(shù)量Tab.5 Optimal classification correctness and band selection number under different methods

4 結(jié) 論

本文重點(diǎn)研究高光譜圖像中的紋理信息，對(duì)比分析了傳統(tǒng)紋理識(shí)別方法GLCM和本文應(yīng)用的多重分形提取紋理的方法做對(duì)比。結(jié)果表明，從運(yùn)行速度角度分析，ABS降維方法優(yōu)于K-L降維方法；從運(yùn)行準(zhǔn)確率角度分析，K-L降維方法高于ABS降維方法。其中在分類選擇上我們都選擇有監(jiān)督分類，SVM分類效果高于BP，其中基于梯度密度函數(shù)的SVM準(zhǔn)確率最高為97.91%。

本文采用紋理圖像密度函數(shù)中的圖像亮度函數(shù)、圖像梯度函數(shù)、圖像拉普拉斯算子函數(shù)運(yùn)算不同盒維數(shù)值具有顯著相關(guān)性?？梢宰鳛檠芯繉?duì)象不同屬性和指標(biāo)進(jìn)行研究。多重分形值可以作為基于紋理識(shí)別木材的重要參數(shù)，多重分形值與其它指標(biāo)相結(jié)合進(jìn)行木材識(shí)別與提取研究將成為下一步研究側(cè)重點(diǎn)。