朱 宏 偉
(吉林農(nóng)業(yè)科技學(xué)院 網(wǎng)絡(luò)信息中心, 吉林 吉林 132101)
目前, 圖像融合技術(shù)在醫(yī)學(xué)輔助診斷中的作用越來越重要, 使用不同成像機(jī)制檢測不同的病變區(qū)域, 可為患者提供更全面的檢查. 基于變換域的多尺度變換(MST)在多模態(tài)醫(yī)學(xué)圖像融合中表現(xiàn)優(yōu)異. 傳統(tǒng)的MST方法包括Laplace金字塔變換[1]、 小波變換[2]、 提升小波變換[3]和雙樹復(fù)數(shù)小波變換(DTCWT)[4]等. 作為多尺度變換的代表, 非下采樣剪切波變換(NSST)使用剪切波在多尺度和多方向上進(jìn)行圖像分解, 確保分解過程中具有平移不變性[5]的優(yōu)點(diǎn). 陳廣秋等[6]提出了一種多傳感器圖像融合方案, 該方案設(shè)計了用于低通子帶和高通子帶融合的多級方向引導(dǎo)濾波器; Bhatnagar等[7]提出了一種多模態(tài)醫(yī)學(xué)圖像的融合框架, 該融合框架采用兩種基于相位一致性和方向?qū)Ρ鹊牟煌诤弦?guī)則分別處理低頻系數(shù)和高頻系數(shù). 雖然基于MST的方法性能較好, 但如何設(shè)計融合規(guī)則對于醫(yī)學(xué)圖像融合至關(guān)重要. 隨著機(jī)器學(xué)習(xí)的發(fā)展, 稀疏表示具有良好的去噪能力和空間一致性等特征, 已成為圖像融合的主流方法[8]. Nejati等[9]提出了一種使用字典學(xué)習(xí)的圖像融合方法, 該方法從源圖像的局部塊中進(jìn)行字典學(xué)習(xí), 同時, 使用稀疏表示和Markov隨機(jī)場生成用于圖像融合的決策圖. 本文提出一種基于非下采樣剪切波變換與改進(jìn)稀疏表示(ISR)的多模態(tài)醫(yī)學(xué)圖像融合方法. 實(shí)驗(yàn)結(jié)果表明, 該方法在主觀視覺性能和客觀評價上均具有優(yōu)勢.
(1)
其中:x∈n表示原始信號;α表示稀疏系數(shù);ε表示容錯率.
低頻子帶表示圖像的輪廓信息, 可近似為源圖像的平滑版本. 傳統(tǒng)方法是從源圖像中減去均值圖像處理稀疏問題, 導(dǎo)致低頻子帶包含大量細(xì)節(jié)信息, 影響了低頻子帶的融合. 為解決該問題, 本文提出一種改進(jìn)的稀疏表示算法處理低頻子帶. 首先, 建立一個過完備的字典訓(xùn)練集, 該訓(xùn)練集去除了低頻子帶中的細(xì)節(jié)信息, 從而提高了融合效率. 創(chuàng)建基于低頻子帶過完備字典的過程如下.
1) 用Sobel算子提取圖像8個方向的信息, 公式如下:
I(s)=I-I?O(8),
(2)
其中:I(s)為刪除細(xì)節(jié)信息后的圖像;I為原始圖像;O(8)為在8個不同方向上的Sobel算子; ?表示卷積運(yùn)算.O(8)表示如下:
(3)
2) 通過引導(dǎo)濾波器[10]計算引導(dǎo)圖像:
Fi=akIi+bk, ?i∈ωk,
(4)
其中:I為引導(dǎo)圖像;F為引導(dǎo)圖像I的線性變換;ωk為窗口函數(shù); 線性系數(shù)ak和bk是常數(shù),
(5)
3) 均值濾波器用于對引導(dǎo)后的圖像進(jìn)行處理, 公式如下:
(6)
其中:g(x,y)為原始圖像像素;m,n為窗口大小, 取值為3. 圖1為創(chuàng)建過完備字典訓(xùn)練集的過程.
圖1 過完備字典訓(xùn)練集
經(jīng)過上述處理過程, 共處理40幅自然圖片(自然圖片數(shù)據(jù)庫源自http://decsai.ugr.es/cvg/CG/base.html). 從40幅自然圖像中選擇10 000對8×8圖像塊作為字典學(xué)習(xí)的訓(xùn)練集, 使用K-SVD算法[11]對字典進(jìn)行學(xué)習(xí).
本文算法融合框架如圖2所示. 為表達(dá)方便, 以一對計算機(jī)斷層掃描(CT)和核磁共振(MR)腦圖像為例, 融合過程如下.
步驟1) NSST分解: 用NSST算法分解醫(yī)學(xué)圖像, 用aA和bB表示低頻子帶系數(shù),cAs,l和cBs,l表示在第s層、 第l方向上的高頻子帶系數(shù);
步驟2) 低頻子帶融合: 按照低頻子帶融合規(guī)則對低頻子帶進(jìn)行融合;
步驟3) 高頻子帶融合: 按照高頻子帶融合規(guī)則對高頻子帶進(jìn)行融合;
步驟4) NSST重建: 通過逆NSST算法獲得融合圖像F.
圖2 算法流程
利用NSST算法對源圖像進(jìn)行變換后, 可獲得相應(yīng)的高頻和低頻子帶. 低頻子帶包含大部分能量, 是原始圖像的近似表示. 但低頻子帶的系數(shù)不是稀疏的, 為了獲得更多的稀疏系數(shù), 需采用低頻子帶融合規(guī)則進(jìn)行融合. 低頻子帶融合規(guī)則如下:
(7)
(8)
3) 采用L1最大方法(max-L1)獲取融合后的稀疏系數(shù)矩陣為
(9)
(10)
5) 重復(fù)上述過程直到所有的圖像塊都被融合, 最后將所有融合后的圖像塊整合到最終的低頻子帶LF中.
假設(shè)cAs,l和cBs,l是在第s層、 第l方向上的高頻子帶系數(shù), 則高頻子帶融合規(guī)則為
(11)
本文實(shí)驗(yàn)采用MATLAB R2016a, CPU為Intel(R) Core (TM) i7-7700, 內(nèi)存為8 GB. 為了驗(yàn)證本文算法的有效性, 使用兩種不同模態(tài)的醫(yī)學(xué)圖像, 其中包括20對MR-T1與MR-T2圖像, 20對MR與CT圖像. 所有圖像均來自美國哈佛大學(xué)創(chuàng)建的The Whole Brain Atlas醫(yī)學(xué)圖像數(shù)據(jù)庫(http://www.med.harvard.edu/AANLIB/), 所有源圖像的大小均為256×256, 且實(shí)驗(yàn)所用圖像都經(jīng)過了配準(zhǔn). 為了驗(yàn)證算法的有效性, 將本文算法與基于邊界尋找方法(BF)[13]、 基于密集尺度不變特征變換(DSITF)[14]、 基于卷積稀疏表示方法(CSR)[15]、 基于卷積稀疏性的形態(tài)成分分析方法(CSMCA)[16]和基于四叉樹的方法(QB)[17]進(jìn)行比較. 由于融合質(zhì)量評估尚未形成統(tǒng)一的評價指標(biāo), 在實(shí)際應(yīng)用中只能根據(jù)經(jīng)驗(yàn)進(jìn)行選擇, 這將導(dǎo)致評估結(jié)果缺乏說服力. 因此, 本文采用4個不同方面的評價指標(biāo)進(jìn)行驗(yàn)證, 并對融合效果進(jìn)行綜合評估, 其中包括局部重要性質(zhì)量指數(shù)(Q0)[18]、 Piella的結(jié)構(gòu)相似性度量(Qw)、 熵(EN)和標(biāo)準(zhǔn)差(SD). 圖3為MR-T1與MR-T2圖像的融合結(jié)果.
圖3 MR-T1與MR-T2圖像融合結(jié)果
由圖3可見, 利用BF算法得到的融合圖像只存在MR-T1圖像的信息, 缺少M(fèi)R-T2圖像的信息; 利用DSIFT算法得到的融合圖像出現(xiàn)多處塊效應(yīng); 利用CSR與QB算法得到的融合圖像也存在塊效應(yīng); 利用CSMCA算法得到的融合圖像雖然沒有塊效應(yīng), 但亮度較低; 利用本文算法得到的融合圖像既不存在塊效應(yīng), 也未損失源圖像的重要信息.
表1列出了不同算法處理MR-T1與MR-T2圖像融合的評價指標(biāo), 所有值均為20對醫(yī)學(xué)圖像的均值. 由表1可見, 本文算法4種評價指標(biāo)的值最高, 表明本文算法在能量存儲和細(xì)節(jié)處理上都優(yōu)于其他5種對比算法. 圖4為CT與MR圖像融合的結(jié)果. 由圖4可見, 利用BF算法得到的融合圖像只包含CT圖像的信息, 缺少M(fèi)R圖像的信息, 融合效果較差, 不能為醫(yī)生提供輔助診斷作用; 利用DSIFT算法得到的融合圖像在眼球周圍存在明顯的塊效應(yīng); 利用CSR算法得到的融合圖像存在偽影; 利用CSMCA算法得到的融合圖像較暗, 未突顯CT圖像中的骨骼信息; 利用QB算法得到的融合圖像缺少CT圖像的重要信息, 錯誤融合的區(qū)域較多; 利用本文算法得到的融合圖像亮度大, 對比度高, 可為醫(yī)生診斷提供重要的參考. 表2為CT與MR圖像的融合結(jié)果. 由表2可見, 由于BF算法缺少M(fèi)R圖像的內(nèi)容, 所以Qw的值最低; CSMCA算法的亮度最低, 所以SD的值最低; 本文算法與其他算法相比4種評價指標(biāo)均最高.
表1 不同算法處理MR-T1和MR-T2圖像融合結(jié)果比較
表2 不同算法處理CT與MR圖像融合結(jié)果比較
圖4 CT與MR圖像融合結(jié)果
綜上可見, 本文算法無論在視覺感官上, 還是在客觀評價指標(biāo)上都優(yōu)于其他對比算法, 表明本文算法融合效率最高, 可為醫(yī)生判斷病情提供輔助診斷作用.