朱文維,李俊峰
(浙江理工大學 機械與自動控制學院,杭州 310018)
隨著醫學成像技術的不斷發展,醫學影像已成為臨床診斷和治療中不可或缺的一部分.而這就要求醫學影像能夠具備較高的空間分辨率,從而對高密度組織,軟組織有高質量的成像效果,這一要求僅僅依靠單模態醫學圖像無法實現.例如,計算機斷層掃描技術(CT)能夠對硬組織的病變位置進行準確定位,有利于顯現高密度的骨骼特征;磁共振血管成像(MRI)可以對低密度的軟組織,血管等清晰成像;單光子發射計算機斷層成像(SPECT)主要用來顯示人體器官,腫瘤組織血流和代謝功能[1].由于單一模態圖像成像機制不同會造成的圖像細節信息紋理缺失,有必要將不同模態的圖像進行信息融合,以實現優勢互補,使圖像信息更豐富,形成全面清晰的醫學影像.
融合后的圖像應該盡可能的保留源圖像的信息,不造成細節部分的缺失,同時也要展現紋理特征,更加有效地診斷醫學圖像.近來,基于變換域的圖像融合算法是研究熱點,常用的變換有小波變換[2],輪廓波變換,非下采樣輪廓波變換(NSCT)[3],剪切波變換等.這些變換或使圖像的細節邊緣部分表現能力不佳;或是不具備平移不變性,會引起圖像失真,無法得到清晰的邊緣;或是計算過于復雜,需耗費大量時間,使得實用性大大降低.為此,NSST被提出,作為剪切波變換的改進模型,它突出了圖像中的顯著特性,并且成功地應用于醫學和非醫學圖像處理應用中.NSST具有靈敏度高,方向性強,運算速度快等特點,優于小波,輪廓波等其它變換,具有處理多方向變化的能力,是一種非常強大的變換[1].
國內外學者對基于NSST的圖像融合進行了相關研究,并做出諸多貢獻.陳貞等人[4]采用基于NSST的醫學圖像融合算法,低頻采用區域能量加權,高頻采用區域能量和平均梯度加權的融合方法,取得較好的數據結果,但在一定程度上忽略了低頻子帶的邊緣信息;趙丹等人[5]提出區域特征值的融合方法,雖然可以綜合反映圖像的區域特征,但是沒有考慮到高低頻系數之間的影響;吳方等人[6]提出使用蜘蛛優化算法搜索函數的最優參數集合,進一步增強融合圖像的質量,但該算法的計算效率比較低;Mohammed等人[7]提出的一種在改進的NSCT中結合多模態圖像的方法,可以提供高質量的圖像,有助于臨床診斷的有效性,但是SR和PCNN計算比較費時;Padmavathi等人[8]提出的基于DTCWT-PCA和DTCWT-PSO的多模態醫學圖像融合方法,DTCWT-PCA與DTCWT-PSO算法具有較好的性能,特征向量尺寸較小,計算效率較高,但是實際融合效果不佳;文獻[9]提出利用空間域主成分分析(PCA)法從大數據集中提取相關信息,從小波變換中提取圖像的方向性細節,此方法雖然可以增強圖像的方向性特征,但是在一定程度上存在冗余細節,噪聲狀偽影,失真等問題;郭楚等人[10]通過多個坐標系的轉換關系將超聲坐標系配準到CT圖像坐標系里,最終將實時超聲圖像統一到CT圖像中,并在軟件中測量融合誤差,融合誤差得到了減小,但是實時超聲圖像與CT圖像的融合精度還是不夠高;鄧志華等人[11]提出在融合低頻成分時,采用傳統的絕對值取大的融合策略,以保留源圖像更多的亮度信息,但是卻沒有考慮到圖像整體的輪廓清晰程度;Yang等人[12]提出一種基于NSST和壓縮感知的融合算法,視覺效果較佳,但是會造成圖像信息的部分缺失.
為了提高圖像整體的區域強度,突出圖像的邊緣信息,提高圖像的對比度,本文提出一種新的基于NSST的融合算法.對于低頻子帶,合成局域平均能量與局域標準差進行子帶選擇;對于高頻子帶,利用改進的拉普拉斯能量和的規則進行融合.主觀視覺效果以及客觀評價均表明,融合后的圖像包含更豐富的紋理,輪廓清晰,對比度較高,視覺效果更佳.
剪切波是在合成小波理論的前提下,通過經典仿射將幾何與多分辨分析結合而成的系統理論.當維數n=2時,具有合成膨脹的仿射系統為[13,14]:

式中,Ψ ∈L2R2,A和(B均)為二維可逆方針,且|d etB|=1.若對于任意的f∈L2R2,AAB(Ψ)都構成一個Parseval框架(也叫緊支撐框架),即:

則稱AAB(Ψ)中的元素為合成小波.其中,矩陣Aj與尺度變換相關聯,矩陣Bl與區域保持的幾何變換相關聯,如旋轉與剪切變換.一般情況下,這種結構不僅允許像小波一樣構造基元素在各個尺度和位置上的Parseval框架,還可以構造基元素在各個方向上的Parseval框架.

其中,φa,s,t(x)為剪切波,式(3)為剪切波系統.可以看出,剪切波是一個涵蓋位置,尺寸和方向的函數集合.
NSST的離散化過程主要分為兩步:多尺度分解和方向局部化.其中,多尺度分解是通過非下采樣金字塔濾波器組來實現的,圖像經K級非下采樣金字塔濾波器分解成1個低頻和K個高頻子帶圖像,而且每個子帶圖像均與原圖像大小相同;方向局部化是通過剪切濾波器實現的[13].標準的剪切波變換中使用的剪切濾波器是在偽極化坐標系中通過窗函數的平移操作實現的,此過程中包含下采樣操作,因而不具有平移不變性.而NSST將標準的剪切濾波器從偽極化網絡系統映射回到笛卡爾左邊系統,從而有效地摒棄了下采樣操作實現其平移不變性[15].如圖1為非下采樣剪切波三層分解示意圖,圖2為NSST梯形支撐區間.
本文提出的醫學圖像融合算法流程圖如圖3所示.首先,將從標準庫取出的已配準源圖像A,B進行NSST分解,分別得到圖像的低頻子帶系數和,高頻子帶系數和(這里,k為方向級數,l為分解層數);接著,結合低頻子帶系數和各層高頻方向子帶系數的特點,通過合成局域平均能量與局域標準差的方法進行低頻子帶的選擇,并利用改進的拉普拉斯能量和的策略進行高頻子帶的融合;最后,對獲得的圖像融合系數進行逆NSST變換,得到融合后的圖像.

圖1 NSST三層分解示意圖

圖2 NSST梯形支撐區間

圖3 醫學圖像融合算法流程圖
低頻子帶反映圖像的輪廓與基本信息.融合后的圖像要想要最大程度地保留源圖像的基礎信息,重點在于恰當地處理低頻系數.常用的低頻區域融合方法有:加權平均法,絕對值取大法,標準差選擇法等[16].加權平均法容易造成圖像整體亮度的下降,視覺效果上,對比度略顯不足,導致圖像邊緣紋理和內部細節信息的模糊;絕對值取大法的優勢在于,可以較大程度地保留源圖像的邊緣紋理細節信息,但是低亮度區域的紋理細節弱化較多;標準差選擇法可以優先選擇邊緣細節豐富的區域信息,但是對噪聲現象抵御力較弱[5,14].為了解決這一問題,本文采取局域平均能量與局域標準差相結合的方法進行低頻部分的融合.具體步驟如下:
1)計算低頻系數的均值CL,區域能量EL以及局域標準差DL.其中:


窗口大小為M×N,窗口函數為
2)根據區域能量EL及局域標準差DL構造圖像區域特征H1I,H2I(I=A,B),及RP(P=A,B).其中,

式(7)中的“1”用于防止DL為零時融合系數失效.R的值越大,對應圖像包含的信息就越豐富.
3)設計低頻系數融合方案如下:

其中,q1=RA/RA+RB+0.01,q2=1?q1為權值.
高頻子帶圖像含有邊緣,紋理等細節信息,對準確地診斷病情十分重要[17].改進的拉普拉斯能量和是圖像清晰度的有效度量,與方差,梯度能量,拉普拉斯能量等清晰度評價指標相比,性能更為優越[16,18],融合后圖像的邊緣和細節保留更加豐富.改進的拉普拉斯(Modified Laplacian,ML)定義為:


f(i,j)為帶通子帶圖像在(i,j)處的系數,“1”為步長.
最終,設計高頻系數融合方案如下:

其中,NS ML1是源圖像A的改進拉普拉斯能量和值,NS ML2是源圖像B的改進拉普拉斯能量和值,另外o1=NS ML1/(NS ML1+NS ML2+0.01),o2=1?o1為權值.
為了驗證本文算法的有效性,選用3種已有算法的數據作為實驗對比數據,并選擇信息熵,空間頻率,標準差和平均梯度對融合后的圖像進行評價.圖4(c)和圖5(c)是文獻[4]提出的基于NSST的方法分解圖像,低頻子帶根據區域能量加權融合系數進行融合,高頻子帶利用平均梯度和區域能量加權的融合策略:圖4(d)和圖5(d)是文獻[20]所提出的低頻子帶利用區域能量加權,高頻子帶利用區域能量取大的融合方法[20];圖4(e)和圖5(e)是文獻[21]提出的低頻區域采用信息熵加權,高頻子帶圖像梯度能量加權的融合方法[21];圖4(f)和圖5(f)是應用本文算法所融合得到的融合圖像.表1和表2是應用不同融合算法產生的醫學融合圖像的質量評價數據.

圖4 正常腦部CT/MRI圖像融合結果圖

圖5 多發性腦梗塞MR-T1/MR-T2圖像融合結果圖
為了實現本文算法的可比性,所有實驗均采用Window 10操作系統,使用Matlab R2014a平臺進行仿真.實驗選取3組經過配準的源圖像,大小均為256×256,顏色深度為8 bit灰度[5].如圖所示,圖4(a)和4(b)是正常腦部CT/MRI圖像,圖5(a)和圖5(b)是多發性腦梗塞MR-T1/MR-T2圖像.
從主觀角度分析,圖4(d)和圖4(e)的組織信息亮度稍顯暗淡,圖4(d)中輪廓較模糊,而圖4(c)在輪廓和對比度上要優于4(d)和4(f),但是紋理邊緣不夠清晰.,圖5(c)細節結構不夠清晰,輪廓信息不完整,對比度稍優于圖5(d)和圖5(e),圖5(d)對比度明顯不夠高,邊緣輪廓也不夠清晰,圖5(e)對比度優于5(d),但是細節信息有所缺失.而從本文算法融合的圖4(f)和圖5(f)可以看出融合質量較高,不僅很大程度上彌補了其他圖像的不足,而且細節信息全面,輪廓細致,紋理清晰,視覺效果良好.

表1 圖4中應用不同融合算法產生的醫學融合圖像的質量評價

表2 圖5中應用不同融合算法產生的醫學融合圖像的質量評價
從客觀角度分析,觀察表1和表2,圖4(d)和圖5(d)以及圖4(e)和圖5(e)的各項評價指標數值普遍較低,尤其是圖5(d)的標準差和平均梯度指數過低,說明對應算法得到的融合圖像紋理特征不夠清晰,對比度不高;圖5(e)的信息熵稍稍高于前面兩種,但是標準差指數較低,細節信息丟失過多.相比其他3種算法,本文算法在客觀評價指數上大部分都優于其余幾種算法,使融合后的圖像信息豐富,輪廓完整,更具實用性.
灰度圖像過后,為了進一步驗證本文算法的通用性,將本文算法運用于彩色圖像與灰度圖像之間的融合.由于彩色圖像是三維的,所以需要先將24位彩色圖像從RGB空間轉換到HIS空間,然后利用NSST分解出I分量圖像與灰度圖像,并用本文算法融合,將融合后得到的灰度圖像代替原來的I分量,最后通過HIS逆變換得到最終的RGB彩色醫學圖像.融合步驟如圖6所示,本文選取兩組SPECT/MRI彩色圖像進行仿真實驗,如圖7(a)和圖7(b)是轉移性肺MRI/SPECT圖像,而圖8(a)和圖8(b)是正常腦組織MRI/SPECT圖像.

圖6 MRI/SPECT圖像融合流程圖

圖7 轉移性肺癌MRI/SPECT圖像融合結果圖

圖8 正常腦組織MRI/SPECT圖像融合結果圖
從主觀角度分析,圖7(c),圖8(c)與圖7(e),圖8(e)整體輪廓不夠清晰,內部信息丟失過多,紋理邊緣不夠理想,圖7(d),圖8(d)整體性效果稍稍好于圖7(c),圖8(c)與圖7(e),圖8(e),但是該圖像的對比度不佳,降低了臨床診斷的實用性.而本文算法融合后的圖像,輪廓清晰,紋理清楚,信息全面.
從客觀數據分析,觀察表3和表4,本文算法融合的圖像空間頻率,平均梯度等評價指數都是較高的,表明圖像像素活躍程度高,圖像清晰.表3中(d),(e)的標準差較低,說明圖像細節有所缺失,表4的(d)中信息熵與空間頻率較低,說明圖像中包含的信息量較少,空間活躍度較低.本文算法得出的評價指標大部分都優于其他算法,所以本文算法也適用于彩色圖像融合,且融合質量較高.

表3 圖7中應用不同融合算法產生的醫學融合圖像的質量評價
綜上分析,本文提出的基于NSST域的醫學圖像融合算法,在有效保留源圖像概貌的同時,融合圖像的灰度分布均勻,細節信息豐富,與源圖像之間相似度高,融合質量優,效果好.是一種效果良好,具有一定有效性與通用性地醫學影像融合方法.

表4 圖8中應用不同融合算法產生的醫學融合圖像的質量評價
醫學圖像融合要求融合后的圖像盡可能保真,準確豐富,全面的體現源圖像含有的信息.基于此,本文提出一種基于NSST的醫學圖像融合算法,低頻子帶系數采用局域平均能量與局域標準差相結合的方法進行融合,這樣圖像的基本信息可以較為完善的保留;高頻子帶系數采用改進的拉普拉斯能量和的策略進行融合,有利于提高圖像抓取細節的能力.實驗結果表明,本文融合算法得到的圖像信息全面,輪廓完整,紋理清晰.在人眼視覺體驗和客觀評價數據上都優于其他3種算法.但是也存在一定的局限性,本文算法是針對醫學圖像展開研究的,將此算法拓寬到其他領域如遙感,紅外與可見光等是下一階段研究的重點.