基于NSST域像素相关分析的医学图像融合

肖明尧 ,  李雄飞 ,  朱芮

吉林大学学报(工学版) ›› 2023, Vol. 53 ›› Issue (09) : 2640 -2648.

PDF (1648KB)
吉林大学学报(工学版) ›› 2023, Vol. 53 ›› Issue (09) : 2640 -2648. DOI: 10.13229/j.cnki.jdxbgxb.20200365
计算机科学与技术

基于NSST域像素相关分析的医学图像融合

作者信息 +

Medical image fusion based on pixel correlation analysis in NSST domain

Author information +
文章历史 +
PDF (1686K)

摘要

针对像素级多模态医学图像融合信息丢失的问题,提出了一种基于非下采样剪切波变换(NSST)的像素相关性分析(PCAS)的图像融合方法。首先,对源图像进行NSST分解,获得高低频子带。然后,利用提出的中心像素方差计算邻域像素与中心像素的强度相关因子,构建邻域像素相关系数矩阵,并提出将相关性拉普拉斯能量和作为高频方向子带的融合规则。再次,计算低频子带中心像素能量以及邻域像素能量梯度信息,得到低频融合决策图。最后,通过逆变换得到融合结果图像。磁共振图像(MRI)和计算机断层扫描(CT)、单光子发射计算机断层成像(PET)、正电子发射断层成像(SPECT)的脑部图像融合实验结果表明,本文融合方法可以很好地保留源图像的显著信息和纹理细节。

关键词

计算机应用 / 图像处理 / 图像融合 / 非下采样剪切波变换 / 像素相关性

Key words

computer application / image processing / image fusion / non-subsampled shearlet transform(NSST) / pixel correlation

引用本文

引用格式 ▾
肖明尧,李雄飞,朱芮. 基于NSST域像素相关分析的医学图像融合[J]. 吉林大学学报(工学版), 2023, 53(09): 2640-2648 DOI:10.13229/j.cnki.jdxbgxb.20200365

登录浏览全文

4963

注册一个新账户 忘记密码

0 引言

由于成像机制的多样性,不同模态的医学图像显示不同种类的信息,比如软组织、骨骼、细胞代谢等。多模态医学图像融合技术可以将不同模态的医学图像的显著信息进行整合,并通过单一图像进行表达,从而降低了信息的冗余性,提高了临床诊断的效率和准确性。

近年来,许多基于多尺度变换的融合方法被提出。其中,Yang1提出了一种基于离散小波变换的医学图像融合方法。Srivastava等2提出了基于曲波变换的计算图像局部能量的融合方法,取得了很好的效果。Koley等3提出了一种将轮廓波变换与模糊统计结合的多光谱磁共振图像(Magnetic resonance imaging,MRI)融合的方法,提高了脑肿瘤解剖病理信息的可视化效果。上述图像多尺度分解方法虽然取得了较好的效果,但是小波方法分解方向有限,不能有效地逼近图像边缘。轮廓波和曲波分解方法不具有平移不变性,融合图像中会产生伪吉布斯现象,影响图像的视觉效果。为了解决上述图像分解方法的局限性,Easley等4提出了具有平移不变性的非下采样剪切波变换(Nonsubsampled shearlet transform,NSST),实现了对图像的多尺度、多方向分解,能够以最稀疏的方式表示图像的轮廓信息。基于NSST分解的图像融合方法被广泛应用56,它可以根据图像的复杂度选择满足需求的分解方向的数量,能够更好地捕捉图像的纹理细节等重要信息,提升融合图像的质量。

传统融合频率子带的方法大多是根据像素本身或者整幅图像的灰度特征设计融合规则,例如,Zhu等7对高频和低频子图分别采用主成分分析和平均梯度法进行图像融合,该方法虽然有效,但效果并不理想。Devi等8利用平均法分别结合Gabor滤波器组和梯度算法融合低频系数和高频系数。刘哲等9采用绝对值最大规则融合高频子图。这些方法在获取融合决策图时并没有考虑待融合像素与相邻像素之间的相关关系,或者只是采取简单的加权平均和基于像素之间的距离关系进行计算,导致出现融合过程产生能量损失、结果图像对比度较低以及图像部分信息失真等现象。

本文分别提出了新的高频子带、低频子带融合规则。其中,针对高频子带,提出了一种基于中心像素方差的计算方法,通过衡量邻域像素与中心像素之间的强度相关性得到强度相关因子,并结合Sigmoid函数构建了邻域像素融合权重矩阵。同时,提出采用相关性拉普拉斯能量和作为子带融合规则。针对低频子带,依据中心像素能量和邻域像素能量梯度信息计算融合权重,减少了融合过程中的能量损失,保持了融合图像的对比度,提高了其视觉质量。

1 非下采样剪切波变换

NSST对图像的分解操作包含两个步骤:①利用拉普拉斯金字塔对图像进行多尺度分解(Non-subsampled laplacian pyramid,NSLP)。图像经过一次NSLP分解,会得到一个与源图像尺度相同的低频子带LF和一个高频子带HF,以此方式在低频子带上迭代k次,会得到1个低频子带和k+1个高频子带。②利用局部剪切波滤波器(Shearlet filter,SF)对分解得到的高频子带系数进行卷积操作实现方向局部化,可以得到NSST l级分解第d个方向子带 H F l , d

NSST对图像进行2次分解的过程如图1所示。从图1可以看出,2级NSST分解图像可以得到2个方向子带和1个低频子带。仔细观察获得的方向子带,方向数目可以根据图像进行自由选择,能够更灵活地捕获源图像的边缘信息。对医学图像进行NSST分解的结果如图2所示。为了更好地观察方向子带系数,图2中选择1级分解,且方向数为8。其中,图2(a)为计算机断层扫描(CT)源图像;图2(b)~(i)依次为1次NSST分解获得的8个方向上的方向子带。由图2可以看出,NSST分解方法可以有效地提取图像的边缘细节。所以本文采用NSST方法分解源图像A、B,分别获得高频子带 { H F A l , d , H F B l , d }和低频子带 { L F A , L F B }

2 图像融合算法

有效的图像分解算法应尽可能将源图像的低频信息与高频信息分离,使得针对不同层设计的融合规则可以保留该层的显著信息。其中,低频子带保留了源图像最多的能量信息,是源图像的近似图像。图3为基于像素相关分析的医学图像融合框架。

2.1 高频融合

2.1.1 像素相关性分析

为了保留边缘特征的同时减少融合过程中可能发生的能量损失,保持图像的对比度,本文提出了基于中心像素方差 σ '的概念,其定义如下:

σ ( x , y ) ' 2 = 1 M ( m , n ) Ω I ( m , n ) - I ( x , y ) 2

式中: ( x , y )为中心像素坐标; Ω为中心像素的邻域; ( m , n )为邻域内任一点;M为邻域内像素的个数。

图4为中心像素与领域像素空间关系示意图。图4中,中心像素为 I ( x , y ),其邻域 Ω大小为 N × N(本文选择窗口大小为 3 × 3),邻域内像素为黑色。因为一幅图像的灰度值方差可以反映整幅图像的能量分布,所以类似地,提出的中心像素方差 σ ( x , y ) '可以衡量邻域像素相对中心像素的灰度变化程度, σ ( x , y ) '值越大,说明该中心像素的邻域灰度值跳跃性越大; σ ( x , y ) '值越小,说明邻域相对越平滑。基于此,本文提出衡量邻域像素与中心像素的强度相关因子 α,其定义如下:

1 α = I ( m , n ) - I ( x , y ) 2 σ ( x , y ) ' 2 + ε

式中: ε为修正常量,为了防止计算过程中出现除数为0的情况。

在同一中心像素的邻域内, σ ( x , y ) '值固定, I ( m , n ) - I ( x , y ) 2值越小,说明位于 ( m , n )的邻域像素的强度与中心像素越相近,二者强度相关性越高,在计算融合权重时贡献越大;反之, I ( m , n ) - I ( x , y ) 2值越大,二者强度相关性越低,在计算融合权重时贡献越小。所以,利用强度相关因子 α计算邻域像素融合系数 ω

ω = 1 1 + e α

因为Sigmoid函数 Y = 1 1 + e - X定义域为实数域,值域为 ( 0,1 ),适合作为权重系数的基函数。

2.1.2 相关性拉普拉斯能量和

传统的拉普拉斯能量和(Sum of modified laplacian,SML)只计算水平方向和垂直方向的邻域像素信息,忽略了对角邻域像素的灰度相关性。Yin等10增加了对角邻域像素信息的计算,但是其权重是根据像素之间的欧式距离计算的,并没有考虑像素之间能量的相关性。为此,本文提出了相关性拉普拉斯能量和(Correlation-sum of modified laplacian,C-SML)的概念:

C - S M L ( x , y ) = i = - N N   j = - N N C M L ( x + i , y + j )
C M L ( x , y ) = ( ω 1 + ω 2 ) I ( x , y ) - ω 1 I ( x - 1 , y ) - ω 2 I ( x + 1 , y ) + ( ω 3 + ω 4 ) I ( x , y ) - ω 3 I ( x , y - 1 ) - ω 4 I ( x , y + 1 ) + ( ω 5 + ω 6 ) I ( x , y ) - ω 5 I ( x - 1 , y - 1 ) - ω 6 I ( x + 1 , y + 1 ) + ( ω 7 + ω 8 ) I ( x , y ) - ω 7 I ( x - 1 , y + 1 ) - ω 8 I ( x + 1 , y - 1 )

式中: ω ii=1,2,…,8)为水平方向、垂直方向和对角方向的权重系数。

将灰度能量纳入权重计算可以更好地保留图像的细节信息以及保持融合图像的对比度。

2.1.3 高频子带融合规则

根据式(6),比较待融合图像中相同位置像素的C-SML值,可以获得图像A高频子带的融合决策图 m a p A l , d,通过式(7)可以获得融合后的高频子带 H F F l , d

m a p A l , d = 1 , C - S M L A l , d C - S M L B l , d 0 , 其他
H F F l , d = H F A l , d × m a p A l , d + H F B l , d × ( 1 - m a p A l , d )

2.2 低频融合

正如前面所说,低频信息保留的多少能够影响融合图像整体的视觉效果。如果能量在融合的过程中发生损失,则融合图像会产生斑驳或者失真的现象。为了更好地保留像素点携带的能量,同时考虑邻域像素的能量变化,本文提出了基于邻域像素能量梯度的低频融合规则。

2.2.1 中心像素能量

为了保持图像的对比度和视觉质量,首先计算中心像素能量 E ( x , y )

E ( x , y ) = i = - N   N j = - N N W L × L F ( x + i , y + j ) 2

因为低频子带图像接近平滑,方差较小,可以直接利用像素间的距离关系计算能量衰减,所以式(8)中的 W L定义如下所示。

以位于 ( x , y )的像素点为中心,其邻域像素 ( m , n )对应的权重 W L , ( m , n )为:

W L , ( m , n ) = 2 2 × N - ( x - m + y - n ) ( a , b ) Ω 2 2 × N - ( x - a + y - b )

式中:邻域窗口 Ω的尺寸为 N × N

2.2.2 邻域像素能量梯度

在低频图像融合过程中,只考虑源图像单一像素能量会导致融合图像灰度信息的不连续,所以本文引入邻域像素能量梯度信息,进一步修正低频子带的融合权重。与C-SML类似,这里选择 3 × 3的邻域窗口计算邻域像素能量梯度(Neighborhood energy gradient,NEG):

N E G ( x , y ) = i = - N   N j = - N N L E G ( x + i , y + j ) 2
L E G ( x , y ) 2 = ( m , n ) Ω L F ( x , y ) - L F ( m , n ) 2

式中:局部能量梯度(Local energy gradient,LEG)为邻域像素与中心像素的能量梯度差值。

2.2.3 低频子带融合规则

低频子带是源图像的近似图像,包含少量细节信息,同时考虑图像像素强度和强度变化梯度,可以更好地捕获低频子带的显著信息。所以,低频子带融合决策图定义为:

m a p A = 1 , E A · N E G A E B · N E G B 0 , 其他

式中: P · Q表示矩阵对应位置元素相乘。

通过式(13)可以获得融合后的低频子带 L F F

L F F = L F A × m a p A + L F B × ( 1 - m a p A )

2.3 图像融合

经过融合低频子带和高频子带,分别得到了对应子带的融合子图。然后通过逆NSST变换,可以获得最终的融合图像。本文提出的融合方法可以总结为算法1。

算法1 基于NSST域像素相关分析的融合算法(NSST_PCAS)

输入:待融合源图像A、B。

输出:融合图像F。

主要步骤:

(1)源图像 I ( I = A B )经过NSST分解,得到高频子带 H F I l , d和低频子带 L F I

(2)根据式(1)(2)(3)计算高频子带 H F I l , d邻域像素融合系数ω

(3)根据式(4)(5)计算高频子带 H F I l , d相关性拉普拉斯能量和(C-SML)。

(4)根据式(6)计算高频子带的融合决策图 m a p I l , d

(5)根据式(7)获得融合后的高频子带 H F F l , d

(6)根据式(9)计算邻域像素权重矩阵 W L,并代入式(8)计算像素能量。

(7)根据式(10)(11)计算低频子带 L F I的邻域像素能量梯度。

(8)将式(8)(10)的计算结果代入式(12),计算低频子带的融合决策图 m a p I

(9)根据式(13)获得融合后的低频子带 L F F

(10)通过逆NSST变换获得融合图像F。

3 实验结果

3.1 实验设置

为了全面评估本文提出的融合方法,从主观和客观两个方面对融合图像进行评价。仿真实验环境为Intel(R)Core(TM)i7-9700,3.00 GHz CPU,8.00 GB RAM,MatlabR2016b。

本文使用的数据来自哈佛医学院全脑图库(http://www.med.harvard.edu/aanlib/home.html)。实验采用了10组MRI和CT图像、4组MRI和功能图像进行实验。

将本文提出的医学图像融合算法与以下9种方法进行对比:双边滤波法(CBF)11、稀疏表示和形态学分析法(CSMCA)12、曲波变换法(CVT)13、离散小波变换法(DWT)、引导滤波法(GFF)14、多尺度变换稀疏表示法(MSTSR)15、非下采样轮廓波变换法(NSCT)、局部拉普拉斯能量法(NSCT_LLE)16和自适应脉冲耦合神经网络法(NSST_PCNN)10。这些方法的实验参数都是作者给出的默认值,没有经过其他改动。

3.2 主观评价

鉴于医学图像应用的特殊性,其主观视觉效果至关重要。图5图6分别为2组MRI和CT图像在不同算法下的融合结果。

图5图6中可以看出:CVT、DWT和NSCT算法的融合图像整体对比度较低,丢失了源图像的部分能量信息。CBF、GFF方法保留源图像的边缘信息较好,但是在CT图像的轮廓部分和MRI图像的中心部分能量损失较多;而MSTSR保留信息较CBF和GFF方法多一些,但是视觉效果仍然较差。对比以上方法的融合图像,CSMCA、NSCT_LLE、NSST_PCNN和本文提出的方法在保留源图像的信息方面表现得更好一些。其中,CSMCA算法的融合图像对比度弱一些,而NSST_PCNN方法的融合结果中,颅骨轮廓部分残留过多MRI的纹理信息,丢失了对应位置CT图像的低频能量。

3.3 客观评价指标

除了主观视觉效果,图像评价指标可以客观地评价融合图像的质量。本文选择了6个评价指标:熵(Entropy,EN)、融合质量指数(Q 0)、结构相似度(Q EQ W)、结构相似性(Structural similarity, SSIM)和视觉信息保真度融合(Visual information fidelity fusion,VIFF)。指标值越大,说明融合图像包含的源图像的信息越多,与源图像的相似性越高,融合质量越好。

表1表2分别为上述两组MRI-CT融合图像的客观评价结果,粗体表示最高分。从表1表2可以看出,本文方法在Q 0Q EQ W、SSIM和VIFF指标上都得到了最高分。说明本文算法得到的融合图像总体质量最好、保留源图像的纹理细节最好、结构相似性最高,融合效果最好。综合主观评价结果分析得出,本文方法优于视觉效果相当的NSCT_LLE融合算法。此外,表1表2的客观评价结果表明,NSST_PCNN方法在熵值取得了最高分,而本文方法得分位居第二,但是两个评分差距并不大,在可接受范围内,说明本文融合算法可以较好地保留源图像的信息和能量,具有较高的性能。

3.4 MRI图像与功能图像融合

除了解剖图像以外,功能图像在临床诊断中也应用广泛,其中包括单光子发射计算机断层成像(Single-photon emission computed tomography,SPECT)和正电子发射断层成像(Positron emission tomography,PET),统称为发射型计算机断层成像(Emission computed tomography,ECT)。由于不同的成像机制,功能图像可以提供身体组织部分代谢信息,所以MRI图像与ECT图像融合也有很大的应用价值。

为了验证本文方法的有效性和可扩展性,做以下实验:①将ECT图像通过颜色空间转换从RGB颜色空间映射到YCrCb颜色空间,并提取光照分量作为待融合图像;②利用本文算法融合MRI图像与步骤①提取的光照分量图像;③将融合后的图像作为新的光照分量,再次转换回RGB颜色空间,获得最终融合图像。

图7图8分别为不同算法下MRI图像与SPECT图像、PET图像的融合结果。其中,GFF图像虽然很好地保留了MRI图像的纹理,但是丢失了大部分的彩色信息,效果不理想。CSMCA和MSTSR图像出现局部性图像失真,说明基于稀疏表示的方法在融合过程中会丢失部分MRI图像的能量。

CVT、DWT和NSCT的结果均偏暗,但是纹理信息保持得相对完整,说明基于变换域的方法对细节信息敏感,但是灰度信息会有些许损失。在图7中,CBF的右下角有部分MRI图像冗余的纹理信息。NSCT_LLE、NSST_PCNN和本文方法在主观感受上都能够保留MRI图像的细节部分,以及ECT图像的伪彩信息。表3表4为MRI-ECT融合图像的客观评价结果,粗体表示最高分。从表3可以看出:GFF方法在Q E指标上表现最好,本文方法排在第2位,说明GFF在保留纹理信息方面表现最优,但是结合其他客观指标和主观感受,本文方法优于GFF方法。在表4中,虽然在Q w指标上,本文方法得分低于NSST_PCNN和NSCT_LLE方法,但是差距并不大,属于可接受范围内。综合分析以上实验结果,说明本文方法在MRI图像与ECT图像融合方面也具有可行性。

4 结束语

本文提出了在NSST域内,基于像素相关性分析(Pixel correlation analysis,PCAS)的医学图像融合算法。首先,利用NSST分解待融合的医学图像,获得高频和低频子带。在高频子带,利用邻域像素相关性计算邻域像素相关系数矩阵,设计了融合系数矩阵,并提出采用相关性拉普拉斯能量和作为融合规则。在低频子带,计算像素能量以及邻域像素能量梯度信息,构建能够保持能量信息的融合决策图。最终,通过逆NSST变换获得融合医学图像。主观评价和客观指标结果表明,本文提出的图像融合算法具有有效性。此外,通过RGB-YCrCb颜色空间变换,验证了本文算法在融合MRI图像与功能图像方面也有很好的稳定性,能够保留MRI图像的纹理信息和功能图像的代谢信息。

参考文献

[1]

Yang Y. Multimodal medical image fusion through a new DWT based technique[C]//2010 4th International Conference on Bioinformatics and Biomedical Engineering, Chengdu, China, 2010: 1-4.

[2]

Srivastava R, Prakash O, Khare A. Local energy-based multimodal medical image fusion in curvelet domain[J]. IET Comput Vision, 2016, 10(6): 513-527.

[3]

Koley S, Galande A, Kelkar B, et al. Multispectral MRI image fusion for enhanced visualization of meningioma brain tumors and edema using contourlet transform and fuzzy statistics[J]. Journal of Medical and Biological Engineering, 2016, 36(4): 470-484.

[4]

Easley G, Labate D, Lim W Q. Sparse directional image representations using the discrete shearlet transform[J]. Applied and Computational Harmonic Analysis, 2008, 25(1): 25-46.

[5]

高印寒, 陈广秋, 刘妍妍. 基于图像质量评价参数的非下采样剪切波域自适应图像融合[J]. 吉林大学学报:工学版, 2014, 44(1): 225-234.

[6]

Gao Yin-han, Chen Guang-qiu, Liu Yan-yan. Adaptive image fusion based on image quality assessment parameter in NSST system[J]. Journal of Jilin University(Engineering and Technology Edition), 2014, 44(1): 225-234.

[7]

Vishwakarma A, Bhuyan M K. Image fusion using adjustable non-subsampled shearlet transform[J]. IEEE Transactions on Instrumentation and Measurement, 2018, 68(9): 3367-3378.

[8]

Zhu Y L, Zhou X Y, Li X W, et al. Algorithm of medical image fusion based on Laplasse pyramid and PCA[C]//IOP Conference Series: Materials Science and Engineering, Shanghai, China, 2019: 490(4): No.042030.

[9]

Devi M S, Balamurugan P. Local energy match based non-sub sampled contourlet Transform for multi modal medical image fusion[J]. International Journal of Engineering and Technology, 2018, 7(2): 165-169.

[10]

刘哲, 徐涛, 宋余庆, 基于NSCT变换和相似信息鲁棒主成分分析模型的图像融合技术[J]. 吉林大学学报:工学版, 2018,48(5): 1614-1620.

[11]

Liu Zhe, Xu Tao, Song Yu-qing, et al. Image fusion technology based on NSCT transform and robust principal component analysis model of similarity information[J]. Journal of Jilin University (Engineering and Technology Edition), 2018,48(5):1614-1620.

[12]

Yin Ming, Liu Xiao-ning, Liu Yu, et al. Medical image fusion with parameter-adaptive pulse coupled neural network in nonsubsampled shearlet transform domain[J]. IEEE Transactions on Instrumentation and Measurement, 2018, 68(1): 49-64.

[13]

Kumar B K S. Image fusion based on pixel significance using cross bilateral filter[J]. Signal, Image and Video Processing, 2015, 9(5): 1193-1204.

[14]

Liu Y, Chen X, Ward R K, et al. Medical image fusion via convolutional sparsity based morphological component analysis[J]. IEEE Signal Processing Letters, 2019, 26(3): 485-489.

[15]

Nencini F, Garzelli A, Baronti S, et al. Remote sensing image fusion using the curvelet transform[J]. Information Fusion, 2007, 8(2): 143-156.

[16]

Li S, Kang X, Hu J. Image fusion with guided filtering[J]. IEEE Transactions on Image Processing, 2013, 22(7): 2864-2875.

[17]

Liu Yu, Liu Shu-ping, Wang Zeng-fu. A general framework for image fusion based on multi-scale transform and sparse representation[J]. Information Fusion, 2015, 24: 147-164.

[18]

Zhu Zhi-qin, Zheng Ming-yao, Qi Guan-qiu, et al. A phase congruency and local Laplacian energy based multi-modality medical image fusion method in NSCT domain[J]. IEEE Access, 2019, 7: 20811-20824.

基金资助

国家自然科学基金项目(61801190)

吉林省自然科学基金项目(20180101055JC)

国家博士后科研基金项目(2017M611323)

吉林省教育厅科学研究项目(JJKH20230920KJ)

AI Summary AI Mindmap
PDF (1648KB)

292

访问

0

被引

详细

导航
相关文章

AI思维导图

/