基于NSST与改进稀疏表示的医学图像融合方法

朱 宏 伟

(吉林农业科技学院 网络信息中心, 吉林 吉林 132101)

目前, 图像融合技术在医学辅助诊断中的作用越来越重要, 使用不同成像机制检测不同的病变区域, 可为患者提供更全面的检查. 基于变换域的多尺度变换(MST)在多模态医学图像融合中表现优异. 传统的MST方法包括Laplace金字塔变换[1]、 小波变换[2]、 提升小波变换[3]和双树复数小波变换(DTCWT)[4]等. 作为多尺度变换的代表, 非下采样剪切波变换(NSST)使用剪切波在多尺度和多方向上进行图像分解, 确保分解过程中具有平移不变性[5]的优点. 陈广秋等[6]提出了一种多传感器图像融合方案, 该方案设计了用于低通子带和高通子带融合的多级方向引导滤波器; Bhatnagar等[7]提出了一种多模态医学图像的融合框架, 该融合框架采用两种基于相位一致性和方向对比的不同融合规则分别处理低频系数和高频系数. 虽然基于MST的方法性能较好, 但如何设计融合规则对于医学图像融合至关重要. 随着机器学习的发展, 稀疏表示具有良好的去噪能力和空间一致性等特征, 已成为图像融合的主流方法[8]. Nejati等[9]提出了一种使用字典学习的图像融合方法, 该方法从源图像的局部块中进行字典学习, 同时, 使用稀疏表示和Markov随机场生成用于图像融合的决策图. 本文提出一种基于非下采样剪切波变换与改进稀疏表示(ISR)的多模态医学图像融合方法. 实验结果表明, 该方法在主观视觉性能和客观评价上均具有优势.

1 稀疏表示

(1)

其中:x∈n表示原始信号;α表示稀疏系数;ε表示容错率.

2 改进稀疏表示

低频子带表示图像的轮廓信息, 可近似为源图像的平滑版本. 传统方法是从源图像中减去均值图像处理稀疏问题, 导致低频子带包含大量细节信息, 影响了低频子带的融合. 为解决该问题, 本文提出一种改进的稀疏表示算法处理低频子带. 首先, 建立一个过完备的字典训练集, 该训练集去除了低频子带中的细节信息, 从而提高了融合效率. 创建基于低频子带过完备字典的过程如下.

1) 用Sobel算子提取图像8个方向的信息, 公式如下：

I(s)=I-I⊗O(8),

(2)

其中:I(s)为删除细节信息后的图像；I为原始图像；O(8)为在8个不同方向上的Sobel算子; ⊗表示卷积运算.O(8)表示如下：

(3)

2) 通过引导滤波器[10]计算引导图像:

Fi=akIi+bk, ∀i∈ωk,

(4)

其中:I为引导图像；F为引导图像I的线性变换；ωk为窗口函数; 线性系数ak和bk是常数,

(5)

3) 均值滤波器用于对引导后的图像进行处理, 公式如下：

(6)

其中:g(x,y)为原始图像像素;m,n为窗口大小, 取值为3. 图1为创建过完备字典训练集的过程.

图1 过完备字典训练集

经过上述处理过程, 共处理40幅自然图片(自然图片数据库源自http://decsai.ugr.es/cvg/CG/base.html). 从40幅自然图像中选择10 000对8×8图像块作为字典学习的训练集, 使用K-SVD算法[11]对字典进行学习.

3 融合策略

本文算法融合框架如图2所示. 为表达方便, 以一对计算机断层扫描(CT)和核磁共振(MR)脑图像为例, 融合过程如下.

步骤1) NSST分解： 用NSST算法分解医学图像, 用aA和bB表示低频子带系数,cAs,l和cBs,l表示在第s层、 第l方向上的高频子带系数;

步骤2) 低频子带融合： 按照低频子带融合规则对低频子带进行融合;

步骤3) 高频子带融合： 按照高频子带融合规则对高频子带进行融合;

步骤4) NSST重建： 通过逆NSST算法获得融合图像F.

图2 算法流程

3.1 低频子带融合规则

利用NSST算法对源图像进行变换后, 可获得相应的高频和低频子带. 低频子带包含大部分能量, 是原始图像的近似表示. 但低频子带的系数不是稀疏的, 为了获得更多的稀疏系数, 需采用低频子带融合规则进行融合. 低频子带融合规则如下：

(7)

(8)

3) 采用L1最大方法(max-L1)获取融合后的稀疏系数矩阵为

(9)

(10)

5) 重复上述过程直到所有的图像块都被融合, 最后将所有融合后的图像块整合到最终的低频子带LF中.

3.2 高频子带融合规则

假设cAs,l和cBs,l是在第s层、 第l方向上的高频子带系数, 则高频子带融合规则为

(11)

4 实验结果与分析

本文实验采用MATLAB R2016a, CPU为Intel(R) Core (TM) i7-7700, 内存为8 GB. 为了验证本文算法的有效性, 使用两种不同模态的医学图像, 其中包括20对MR-T1与MR-T2图像, 20对MR与CT图像. 所有图像均来自美国哈佛大学创建的The Whole Brain Atlas医学图像数据库(http://www.med.harvard.edu/AANLIB/), 所有源图像的大小均为256×256, 且实验所用图像都经过了配准. 为了验证算法的有效性, 将本文算法与基于边界寻找方法(BF)[13]、 基于密集尺度不变特征变换(DSITF)[14]、 基于卷积稀疏表示方法(CSR)[15]、 基于卷积稀疏性的形态成分分析方法(CSMCA)[16]和基于四叉树的方法(QB)[17]进行比较. 由于融合质量评估尚未形成统一的评价指标, 在实际应用中只能根据经验进行选择, 这将导致评估结果缺乏说服力. 因此, 本文采用4个不同方面的评价指标进行验证, 并对融合效果进行综合评估, 其中包括局部重要性质量指数(Q0)[18]、 Piella的结构相似性度量(Qw)、 熵(EN)和标准差(SD). 图3为MR-T1与MR-T2图像的融合结果.

图3 MR-T1与MR-T2图像融合结果

由图3可见, 利用BF算法得到的融合图像只存在MR-T1图像的信息, 缺少MR-T2图像的信息; 利用DSIFT算法得到的融合图像出现多处块效应; 利用CSR与QB算法得到的融合图像也存在块效应; 利用CSMCA算法得到的融合图像虽然没有块效应, 但亮度较低; 利用本文算法得到的融合图像既不存在块效应, 也未损失源图像的重要信息.

表1列出了不同算法处理MR-T1与MR-T2图像融合的评价指标, 所有值均为20对医学图像的均值. 由表1可见, 本文算法4种评价指标的值最高, 表明本文算法在能量存储和细节处理上都优于其他5种对比算法. 图4为CT与MR图像融合的结果. 由图4可见, 利用BF算法得到的融合图像只包含CT图像的信息, 缺少MR图像的信息, 融合效果较差, 不能为医生提供辅助诊断作用; 利用DSIFT算法得到的融合图像在眼球周围存在明显的块效应; 利用CSR算法得到的融合图像存在伪影; 利用CSMCA算法得到的融合图像较暗, 未突显CT图像中的骨骼信息; 利用QB算法得到的融合图像缺少CT图像的重要信息, 错误融合的区域较多; 利用本文算法得到的融合图像亮度大, 对比度高, 可为医生诊断提供重要的参考. 表2为CT与MR图像的融合结果. 由表2可见, 由于BF算法缺少MR图像的内容, 所以Qw的值最低; CSMCA算法的亮度最低, 所以SD的值最低; 本文算法与其他算法相比4种评价指标均最高.

表1 不同算法处理MR-T1和MR-T2图像融合结果比较

表2 不同算法处理CT与MR图像融合结果比较

图4 CT与MR图像融合结果

综上可见, 本文算法无论在视觉感官上, 还是在客观评价指标上都优于其他对比算法, 表明本文算法融合效率最高, 可为医生判断病情提供辅助诊断作用.