本發(fā)明涉及通用圖像數(shù)據(jù)處理技術(shù)領(lǐng)域,具體涉及一種基于方向和尺度描述子的多圖譜圖像分割方法。

背景技術(shù)：

基于多圖譜配準(zhǔn)結(jié)合標(biāo)號融合技術(shù)的圖像分割方法，主要有以下分割方法：

一：直接將圖像全局相似性作為圖譜選擇依據(jù)。這種方法的理論依據(jù)可理解為：全局相似性越高，說明該圖譜灰度圖像與待分割圖像越相似，各組織結(jié)構(gòu)分布越相近。這種方法操作簡單，但對小器官分割時的精度較低。全局相似性主要由大器官決定，小器官局部相似性與全局相似性一致性差。Tong發(fā)現(xiàn)基于全局相似性圖譜選擇分割算法分割腹部MRI圖像，大器官如肝臟的分割結(jié)果明顯優(yōu)于小器官如脾臟等的分割結(jié)果。因此基于全局相似性選擇圖譜不適合用于分割小器官。

二：利用局部相似性作為標(biāo)號融合權(quán)重的分割算法，計算局部相似性的方法多種多樣，Aljabar等人用互信息作為相似性測度，Tong等利用殘差平方和作為相似性測度，這些算法一般都是直接用像素的局部灰度小塊為像素特征計算局部相似性，雖然充分利用了圖譜灰度圖像及待分割圖像的灰度信息，以及圖譜空間結(jié)構(gòu)信息，但是沒有考慮圖譜間亮度非均勻性。

三：構(gòu)建字典，利用字典元素的線性組合近似表達目標(biāo)像素，將線性系數(shù)作為標(biāo)號融合的權(quán)重。不同分割方法構(gòu)建字典的方式各異，Roy利用圖譜訓(xùn)練出可以近似表達所有像素的共有字典，雖然字典訓(xùn)練過程相對簡單，但分割精度低于為不同像素訓(xùn)練出特有字典的分割方法。

技術(shù)實現(xiàn)要素：

本發(fā)明所要解決問題為提供一種基于方向和尺度描述子的多圖譜圖像分割方法，能減小圖譜間亮度不均勻造成的分割誤差，提高分割精度。

本發(fā)明解決上述問題的技術(shù)解決方案是：

一種基于方向和尺度描述子的多圖譜圖像分割方法，包括以下步驟：

(1)讀入數(shù)據(jù)：讀入圖譜灰度圖像I_train和圖譜標(biāo)號圖像L，每一幅圖譜灰度圖像和相應(yīng)標(biāo)號圖像為一組圖譜。假設(shè)有K(本發(fā)明中K＝20)組圖譜，即有K幅灰度圖像I_train＝(I_train1,I_train2,I_train3,......,I_trainK)和K幅圖譜標(biāo)號圖像L＝(L₁,L₂,L₃,......,L_K)，再讀入待分割圖像I_target，假設(shè)圖譜灰度圖像大小為m×n×l。

(2)將圖譜灰度圖像配準(zhǔn)到待分割圖像，并利用形變后的標(biāo)號圖像的先驗信息確定感興趣區(qū)域,具體過程為：

以待分割圖像I_target為參考圖像，圖譜灰度圖像I_train為浮動圖像，將K幅圖譜灰度圖像I_train用DRAMMS配準(zhǔn)方法逐一配準(zhǔn)到待分割圖像I_target，得到使圖譜灰度圖像產(chǎn)生形變的形變場T＝(T₁,T₂,......,T_K)及形變后的圖譜灰度圖像I'_train＝(I'_train1,I'_train2,I'_train3,......,I'_trainK)，用形變場T＝(T₁,T₂,......,T_K)作用在對應(yīng)的圖譜標(biāo)號圖像上(使圖譜灰度圖像I_train1形變成I'_train1的形變場為T₁，將T₁作用到圖譜標(biāo)號圖像L₁上)，使L₁發(fā)生與I_train1相同的形變，L₁變形后的標(biāo)號圖像為L'₁，即L'₁＝T₁(L₁)，同理，將T₁,T₂,......,T_K分別作用在L₁,L₂,L₃,......,L_K上，得到變形后的圖譜標(biāo)號圖像L'＝(L'₁,L'₂,L'₃,......,L'_K)。利用形變后的圖譜標(biāo)號圖像的先驗信息確定感興趣區(qū)域大小為m_r×n_r×l_r，感興趣區(qū)域外的部分可以確定為非海馬體，感興趣區(qū)域即為將海馬體完全包圍的長方體框。

(3)提取圖譜圖像和待分割圖像感興趣區(qū)域內(nèi)目標(biāo)像素的像素特征。具體過程為：

本發(fā)明提出一種新的像素特征，即方向和尺度描述子S。S為像素與搜索域內(nèi)像素間的局部自相似性矩陣，像素間自相似性測度用e的冪函數(shù)衡量，用像素間局部塊之間的殘差平方和的負(fù)值為冪次方。在S中，各點的取值即為該像素與中心像素的局部自相似性。計算像素ν的方向和尺度描述子S過程為：首先確定S的大小p₂(本發(fā)明中當(dāng)S的大小為7×7×7時分割精度最佳，即p₂＝7)，搜索區(qū)域即為p₂×p₂×p₂，用殘差平方和(sum of squared difference,SSD)計算中心像素與p₂×p₂×p₂搜索域內(nèi)的所有像素之間的局部自相似性,以p₂＝7為例說明像素ν的方向和尺度描述子S的計算過程：首先，以像素ν為中心，7×7×7為搜索區(qū)域，提取搜索區(qū)域內(nèi)343(7×7×7＝343)個像素的局部小塊(局部小塊半徑為p₁)，若p₁＝3，即取出343個大小為3×3×3的局部小塊，計算所有343個局部小塊與中心像素局部小塊間的殘差平方和，得到中心像素與搜索區(qū)域內(nèi)所有像素之間的殘差平方和矩陣，記為D_p(I,x,p₂)；然后計算D_p(I,x,p₂)的均值V；最后用D_p(I,x,p₂)除以V，再取e的負(fù)冪次方即為該像素的方向和尺度描述子S。公式(Ⅰ)計算像素的方向和尺度描述子S：

為d_p(I,x,p₂)的均值。p₂為方向和尺度描述子大小，本發(fā)明中p₂＝7。D_p(I,x,p₂)表示搜索區(qū)域(以像素x為中心，搜索半徑為7的搜索區(qū)域)內(nèi)所有像素與中心像素之間的局部自相似性。分別以x與y為中心的兩局部小塊(本發(fā)明中取3×3×3)的殘差平方和，d_p(I,x,y，p₁)表示圖像I中，像素x與像素y之間的局部殘差和矩陣，d_p(I,x,y，p₁)的計算公式為

p₁為計算殘差平方和的局部小塊半徑。

(4)用圖譜圖像像素的方向和尺度描述子構(gòu)建待分割圖像目標(biāo)像素的方向和尺度描述子構(gòu)建特有字典；為每個待分割圖像的目標(biāo)像素構(gòu)建能近似表達像素特征的特有字典。具體過程為：

提取待分割圖像感興趣區(qū)域內(nèi)所有目標(biāo)像素的方向和尺度描述子S_test，S_test包含N(N＝m_r×n_r×l_r)個方向和尺度描述子，提取配準(zhǔn)后圖譜灰度圖像感興趣區(qū)域內(nèi)所有目標(biāo)像素的方向和尺度描述子，構(gòu)成訓(xùn)練圖譜方向和尺度描述子全集S_train，則S_train包含M(M＝N×K)個方向尺度描述子。從S_train中選出不同的像素特征即不同的方向尺度描述子，構(gòu)建待分割圖像目標(biāo)像素的特有字典，即從全集中按照一定原則選取出不同子集構(gòu)建成待分割目標(biāo)像素的字典，本發(fā)明中子集選取原則為取出所有圖譜圖像中與目標(biāo)像素空間位置一致的像素的鄰域內(nèi)像素的方向尺度描述子。構(gòu)建待分割圖像中目標(biāo)像素的方向和尺度描述子S_{test_i}的特有字典D_{_i}的具體過程為：首先確定待分割圖像中第i個目標(biāo)像素在待分割圖像I_target中的空間位置[pos_{_ix},pos_{_iy},pos_{_iz}],找出K幅圖譜圖像中空間位置為[pos_{_ix},pos_{_iy},pos_{_iz}]的像素提取以為中心，p₃×p₃×p₃鄰域內(nèi)所有像素的方向和尺度描述子構(gòu)建待分割圖像感興趣區(qū)域內(nèi)第i目標(biāo)像素的特有字典D_{_i}，將像素的方向和尺度描述子看成元素，則D_{_i}中包含的元素個數(shù)為D_num＝p₃×p₃×p₃×K，這D_num個方向和尺度描述子構(gòu)成的集合即為第i個目標(biāo)像素的字典，依次為待分割圖像感興趣區(qū)域內(nèi)N個目標(biāo)像素構(gòu)建出特有字典，即構(gòu)建N個特征字典。

(5)用圖譜圖像感興趣區(qū)域內(nèi)目標(biāo)像素特征的線性組合近似表達待分割圖像感興趣區(qū)域內(nèi)目標(biāo)像素的像素特征。具體過程為：

用LAE法實現(xiàn)用特有字典D_{_i}中元素的線性組合近似表達待分割圖像的目標(biāo)像素的方向和尺度描述子S_{test_i}。即將特有字典D_{_i}中的元素看成基向量，將待分割圖像的目標(biāo)像素的方向尺度描述子S_{test_i}看成基向量組成的空間中的向量，線性組合的原理即找出S_{test_i}在各基向量上的投影長度，完成線性表達，投影長度即為線性系數(shù)或編碼系數(shù)。根據(jù)公式(Ⅲ)得到稀疏編碼系數(shù)C，即完成稀疏編碼過程。

S_{test_i}＝c₁S_i1+c₂S₂+c₃S_i3+L+c_jS_ij+L+c_MS_iM+ε (Ⅲ)

ε＝(c₁S_i1+c₂S₂+c₃S_i3+L L+c_MS_iM)-S_{test_i}

ε為誤差。將ε作為目標(biāo)函數(shù)，使得ε最小的C即為所求。

(6)加權(quán)融合圖譜標(biāo)號，估計待分割圖像感興趣區(qū)域內(nèi)目標(biāo)像素的標(biāo)號。

將C作為標(biāo)號融合的權(quán)重，融合圖譜標(biāo)號。即將(5)中得到的線性系數(shù)作為標(biāo)號加權(quán)融合時的權(quán)重。特有字典D_{_i}中像素j的方向和尺度描述子S_{train_ij}在線性表達S_{test_i}時所占比重即為線性系數(shù)c_j，而j的標(biāo)號為L'_j，則第j個像素對待分割像素i的標(biāo)號所作貢獻即為c_j×L'_j，D_{_i}中元素在重建待分割圖像中的目標(biāo)像素i的方向和尺度描述子S_{test_i}的貢獻即為編碼系數(shù)，編碼系數(shù)與D_{_i}中對應(yīng)元素的標(biāo)號加權(quán)融合的結(jié)果即為待分割像素i的標(biāo)號估計值。根據(jù)公式(Ⅳ)估計待分割圖像目標(biāo)像素標(biāo)號L_estimate。

L_estimate＝c₁L'₁+c₂L'₂+L L+c_jL'_j+L L+c_ML'_M (Ⅳ)

(7)閾值處理待分割圖像感興趣區(qū)域內(nèi)目標(biāo)像素的標(biāo)號估計值，完成分割。

設(shè)定閾值Th，若標(biāo)號估計值L_estimate大于閾值Th，則將待分割圖像的目標(biāo)像素i標(biāo)記為1，即認(rèn)為該像素是海馬體，同理，若標(biāo)號估計值L_estimate小于閾值Th，則將待分割圖像的目標(biāo)像素i標(biāo)記為0，即認(rèn)為該像素為非海馬體。根據(jù)公式(Ⅴ)判定L_estimate得到目標(biāo)像素標(biāo)號L_test，完成分割。

上述方案中計算殘差平方和矩陣的局部塊大小為p₁×p₁×p₁，方向和尺度描述子搜索域為p₂×p₂×p₂，字典鄰域大小為p₃×p₃×p₃，LAE迭代次數(shù)為10，稀疏編碼過程中非零系數(shù)的個數(shù)為q，閾值為Th。

本發(fā)明的有益效果為：能減小圖譜間亮度不均勻造成的分割誤差，提高分割精度。

附圖說明

圖1為本發(fā)明所述的基于方向和尺度描述子的多圖譜圖像分割方法的流程圖。

圖2為本發(fā)明所述的基于方向和尺度描述子的多圖譜圖像分割方法的流程圖示意圖。

圖3為本發(fā)明所述的基于方向和尺度描述子的多圖譜圖像分割方法的為圖譜灰度圖像(a)和為圖譜標(biāo)號圖像(b)。

圖4為本發(fā)明所述的基于方向和尺度描述子的多圖譜圖像分割方法配準(zhǔn)過程及配準(zhǔn)后的圖譜灰度圖像和標(biāo)號圖像。

圖5為本發(fā)明所述的基于方向和尺度描述子的多圖譜圖像分割方法的目標(biāo)像素確定過程示意圖，曲線包圍的所有像素即為目標(biāo)像素。

圖6為本發(fā)明所述的基于方向和尺度描述子的多圖譜圖像分割方法的分割左海馬時的感興趣區(qū)域圖，矩形框內(nèi)為分割左海馬時的感興趣區(qū)域。

圖7為本發(fā)明所述的基于方向和尺度描述子的多圖譜圖像分割方法的構(gòu)建特有字典及稀疏編碼過程圖。

圖8為本發(fā)明所述的基于方向和尺度描述子的多圖譜圖像分割方法的分割結(jié)果二維對比，左側(cè)線條代表本發(fā)明分割結(jié)果，右側(cè)線條代表專家勾畫分割結(jié)果。

圖9本發(fā)明所述的基于方向和尺度描述子的多圖譜圖像分割方法的分割結(jié)果三維視圖。

圖10本發(fā)明所述的基于方向和尺度描述子的多圖譜圖像分割方法的分割結(jié)果與專家勾畫分割的差異三維視圖。

具體實施方式

實施例1

如圖1-圖10所示，以分割頭部MRI圖像中海馬體為例描述本發(fā)明所述分割方法。具體分割過程如下所述：

步驟1，讀入圖譜灰度圖像I_train和圖譜圖譜標(biāo)號圖像L，每一幅圖譜灰度圖像和對應(yīng)的標(biāo)號圖像為一組圖譜。本實驗中有20組圖譜，即有20幅灰度圖像I_train＝(I_train1,I_train2,I_train3,......,I_traK)和20幅對應(yīng)的圖譜標(biāo)號圖像L＝(L₁,L₂,L₃,......,L_K)，再讀入待分割圖像I_target，圖譜大小為256×256×277。

步驟2，將待分割圖像作為參考圖像，圖譜灰度圖像作為浮動圖像，將所有圖譜灰度圖像I_train用DRAMMS配準(zhǔn)方法逐一配準(zhǔn)到待分割圖像，得到配準(zhǔn)后的圖譜灰度圖像I'_train＝(I'_train1,I'_train2,I'_train3,......,I'_trainK)，然后將配準(zhǔn)得到的形變場作用在相應(yīng)圖譜標(biāo)號圖像上，得到變形后的圖譜標(biāo)號圖像L'＝(L'₁,L'₂,L'₃,......,L'_K)。

步驟3，為了減少計算量，提取出包含所有目標(biāo)像素的感興趣區(qū)域，只處理感興趣區(qū)域內(nèi)的像素，提取感興趣區(qū)域的過程為：將形變后的L'疊加，得到海馬體解剖位置概率圖，即將每一幅形變后的標(biāo)號圖像看成專家，對所有像素判斷，若20幅圖像在某像素位置的標(biāo)號值都為1，即所有專家將該像素判定為海馬體，同理，若20幅圖像在某像素位置的標(biāo)號值都為0，即所有專家將該像素判定為非海馬體，若K_l∈(0,20)幅圖像在某像素位置的標(biāo)號值為1，該像素是海馬體的概率為將所有概率大于0的像素確定為目標(biāo)像素，本實驗中若某像素為海馬體的概率大于0.8，直接將該像素判定為海馬體，同理，將概率小于0.2的像素直接判定為非海馬體。需進一步分割的目標(biāo)像素為概率介于0.2～0.8之間的像素，即為待分割的目標(biāo)像素，假設(shè)目標(biāo)像素為N個。感興趣區(qū)域為包含所有目標(biāo)像素的區(qū)域，如圖5中，紅色曲線包圍的區(qū)域為目標(biāo)像素，藍色矩形框為包含所有目標(biāo)像素的感興趣區(qū)域。

步驟4，提取像素特征描述子——方向和尺度描述子。根據(jù)公式(Ⅰ)計算圖譜灰度圖像感興趣區(qū)域內(nèi)所有像素的方向和尺度描述子S。

p₂為方向和尺度描述子半徑，本發(fā)明中p₂＝5，d_p(I,x,p₂)為以x為中心，5×5×5內(nèi)所有像素與中心像素x之間的殘差平方和矩陣，根據(jù)公式(Ⅱ)得到兩像素間殘差平方和矩陣

d_p(I,x,y)＝∑(I(x+p₁)-I(y+p₁))² (Ⅱ)

式(Ι)中d_p(I,x,y)即表示在圖像I中，像素x與像素y之間的殘差平方和矩陣，p₁為計算像素間殘差平方和的局部小塊大小，本發(fā)明中p₁＝3，I(x+p₁)與I(y+p₁)分別表示以x與y為中心，3×3×3為大小的灰度小塊，d_p(I,x,y)為兩灰度小塊內(nèi)對應(yīng)元素之差的平方和，為d_p(I,x,y)的均值。

步驟5，提取待分割圖像感興趣區(qū)域內(nèi)目標(biāo)像素的方向和尺度描述子S_test，S_test包含N個方向和尺度描述子，提取圖譜圖像所有像素的S_train，從S_train中選擇相應(yīng)的S_train構(gòu)建待分割圖像中感興趣區(qū)域內(nèi)目標(biāo)像素的特有字典。如構(gòu)建待分割圖像第i個目標(biāo)像素S_{test_i}的特有字典D_{_i}的過程為：首先確定第i個目標(biāo)像素在目標(biāo)圖像I_target中的空間位置x_i：[pos_{_ix},pos_{_iy},pos_{_iz}],找出20幅圖譜灰度圖像中空間位置為[pos_{_ix},pos_{_iy},pos_{_iz}]的像素y_i＝(y_i1,y_i2,y_i3,......,y_iK)，提取以y_i為中心，搜索區(qū)域為3×3×3內(nèi)的所有像素的方向和尺度描述子，構(gòu)建待分割圖像第i個目標(biāo)像素的特有字典D_{_i}，D_{_i}中包含的元素個數(shù)為M＝3×3×3×20＝540。將S_{test_i}展成列向量，S_{test_i}為125×1的列向量，將字典中所有S_train展成列向量，字典D_{_i}的大小為125×540。

即

步驟6，用LAE法實現(xiàn)用特有字典D_{_i}中元素的線性組合近似表達待分割圖像中目標(biāo)像素的方向和尺度描述子S_{tes_i}，根據(jù)公式(Ⅲ)得到稀疏編碼系數(shù)C，即完成稀疏編碼過程。

S_{test_i}＝c₁S_i1+c₂S_i2+c₃S_i3+L L+c_MS_iM+ε (Ⅲ)

ε為誤差。即使得ε最小的C即為所求。即

步驟7，將C作為標(biāo)號融合的權(quán)重，加權(quán)融合估計待分割圖像感興趣區(qū)域內(nèi)目標(biāo)像素的標(biāo)號。根據(jù)公式(Ⅳ)估計待分割圖像感興趣區(qū)域內(nèi)目標(biāo)像素標(biāo)號值Li_estimate。

Li_estimate＝c_i1L'_i1+c_i2L'_i2+c_i3L'_i3+L+c_ijL'_ij+L+c_i540L'_i540 (Ⅳ)

步驟8,重復(fù)步驟(5)、(6)和步驟(7)估計出待分割圖像感興趣區(qū)域內(nèi)N目標(biāo)像素的標(biāo)號值L_estimate。

步驟9，設(shè)定閾值0.5，根據(jù)公式(Ⅴ)判定L_estimate得到目標(biāo)像素標(biāo)號L_test。

步驟10，輸出標(biāo)號圖像，完成分割。

以上所述，僅為本發(fā)明的具體實施方式，但本發(fā)明的保護范圍并不局限于此，任何不經(jīng)過創(chuàng)造性勞動想到的變化或替換，都應(yīng)涵蓋在本發(fā)明的保護范圍之內(nèi)。