一种基于区域对比的医学图像融合方法
【专利摘要】本发明公开了一种基于区域对比的医学图像融合方法，属于医学图像处理领域。该方法包括以下步骤：a、获取已配准的多模医学图像，根据多模医学图像的大小制定分块策略，按照从上到下、从左到右的原则对图像进行分块操作；b、计算每一像素在不同分块下的对比度和分块图像的直方图矩阵，并以直方图矩阵作为系数，生成每一像素的区域对比度矩阵；c、应用典型相关分析方法计算每一像素区域对比度矩阵的相关系数，求解像素区域对比度矩阵最相关时的系数矩阵；d、以求解得到的系数矩阵作为权值，生成融合图像。与现有技术相比，本发明的基于区域对比的医学图像融合方法具有融合效果好、执行效率高等特点，具有很好的推广应用价值。
【专利说明】
-种基于区域对比的医学图像融合方法
技术领域
[0001] 本发明设及医学图像处理领域，具体地说是一种基于区域对比的医学图像融合方 法。
【背景技术】
[0002] 随着成像技术的发展，医学图像已成为现代医学诊断和治疗的重要手段。目前，临 床中常用的医学图像主要分为解剖图像和功能影像两种。解剖图像(如CT、MR等)能清晰的 提供器官的解剖结构、形态特征，并且具有较高的分辨率;功能影像(如PETJMRI等)可W准 确提供器官的功能和新陈代谢信息，但是其分辨率较低。单独一种模态的图像只能提供患 者某些方面的特征，不能反映其全部病理信息。因此，在临床实践中，往往需要有机整合多 种影像中的互补信息，融合多种模态的医学图像，为精确诊断和精准治疗提供依据。
[0003] 现有技术中，多模医学图像融合方法按照处理策略可W分为像素级、特征级、决策 级=个层次，像素级融合方法主要有空域融合、变换域融合、智能融合等。空域方法通过对 像素值本身的加权平均或形态学处理进行融合，其计算简单，融合效果较差。变换域方法首 先对医学图像进行金字塔、小波或多尺度几何变换，提取变换域子带系数，针对子带系数制 定不同的融合规则，最后通过反变换生成融合图像，其融合效果较好，但计算复杂度高，实 时性差。智能融合方法通过神经网络、统计学习、模糊集等算法迭代计算融合系数，生成融 合图像，其融合效果好坏往往依赖于学习算法的设计，且同样存在计算复杂度高，实时性差 等问题。总的来讲，现有的像素级融合方法都是基于多模医学图像的结构、形态、颜色、时频 信息等特征制定融合策略，缺乏对像素空间位置关系及对比度等因素的分析；同时，现有方 法一般采取加权平均或绝对值取大作为融合规则，易使融合结果陷入全局最优，导致局部 细节丢失或模糊，难W满足实际应用需求。
【发明内容】
[0004] 本发明的技术任务是针对上述现有技术的不足，提供一种应用像素空间位置关系 及对比度制定融合规则的基于区域对比的医学图像融合方法。
[0005] 本发明的技术任务是按W下方式实现的：一种基于区域对比的医学图像融合方 法，其特点是包括W下步骤：
[0006] a、获取已配准的多模医学图像，根据多模医学图像的大小制定分块策略，按照从 上到下、从左到右的原则对图像进行分块操作；
[0007] b、计算每一像素在不同分块下的对比度和分块图像的直方图矩阵，并W直方图矩 阵作为系数，生成每一像素的区域对比度矩阵；
[000引C、应用典型相关分析方法计算每一像素区域对比度矩阵的相关系数，求解像素区 域对比度矩阵最相关时的系数矩阵；
[0009] d、W求解得到的系数矩阵作为权值，生成融合图像。
[0010] 作为优选，步骤a中优选采用W下方式对图像进行分块操作：图像分块操作分t次进 行，每次分块尺寸分别为多模医学图像大小的1,1/4,1/8…l/2t倍。t为大于等于2的正整数，定 义配准后多模医学图像A和B大小均为mXn(单位:像素），则t的最大取值为 良P f'巨[2, 。
[0011] 步骤b中采用W下方法计算每一像素区域对比度矩阵：
[0012] 定义多模医学图像A和B的分块分别表示为BAi和BBi, i表示第i次分块，i的取值为 大于等于2且小于等于t的正整数，分块内像素点(x，y)的对比度CU, y)i为：
[0013]
[0014]
[001引分块BAi和BBi的直方图矩阵HAi、HBi，贝晦一像素区域对比度矩阵CU, y)为：
[0016] CA(x,y) = [HA(x,y)iXCA(x,y)r--HA(x,y)iXCA(x,y)i---HA(x,y)tXCA(x,y)t]^
[0017] CB(x,y) = [皿（x,y)lXCB(x,y)l???皿（x,y)iXCB(x,y)i???皿（x,y)tXCB(x,y)t]'。
[0018] 步骤C利用典型相关分析求解区域对比度矩阵的系数矩阵的方法包括：
[0019] 首先，生成待融合图像像素点的对比度矩阵向量Z:
[0020] Z=[CA(x,y)CB(x,y)]
[0021] 计算其协方差矩阵S :
[0022]
[0023] 其中，乙ii = Cov(CA(x,y) ,CA(x,y))
[0024] E12= E2i = Cov(CA(x,y) ,CB(x,y))
[0025] E22 = Cov(CB(x,y) ,CB(x,y))
[00%]其次，令U=AiCAU, y)，V=BiCB(x，y)作为两组变量的线性组合，则有，
[0027] Cov(U,V)=AiEi2Bi'
[002引其中，Ai、Bi为相关系数矩阵，即需要求解的系数矩阵，Bi'为Bi的转置矩阵。
[0029] 步骤d的融合结果为：
[0030] C(x,y) =wiA(x,y)+W2B(x,y)
[003^1]其中，wi、W2 为 Ai、Bi 前 r 维系数，:r 为 CA(x,y)和 CB(x,y)的秩。
[0032] 与现有技术相比，本发明的基于区域对比的医学图像融合方法具有W下突出的有 益效果：
[0033] ( - )结合了区域相似性、像素空间位置关系等因素作为融合策略，利用相关系数 制定融合规则，使每一像素点融合结果与原多模图像的相关性最高；
[0034] (二)采用特有的分块方式，兼顾像素的全局及局部对比度，避免了融合结果的局 部最优。
【附图说明】
[0035] 附图1是本发明基于区域对比的医学图像融合方法的流程图；
[0036] 附图2是本发明提供的一例烦脑CT图像实例；
[0037] 附图3是本发明提供的一例烦脑MR图像实例；
[0038] 附图4是基于PCA方法对图2和图捕虫合结果实例；
[0039] 附图5是基于DCT变换方法对图2和图捕虫合结果实例；
[0040] 附图6是基于NSCT方法对图2和图捕虫合结果实例；
[0041] 附图7是基于LPSR方法对图2和图捕虫合结果实例；
[0042] 附图8是基于GFF方法对图2和图捕虫合结果实例；
[0043] 附图9是基于本发明方法对图2和图捕虫合结果实例。
【具体实施方式】
[0044] 参照说明书附图W具体实施例对本发明的基于区域对比的医学图像融合方法作 W下详细地说明。
[0045] 【实施例】
[0046] 如附图1所示，本发明的基于区域对比的医学图像融合方法的具体实施步骤为：
[0047] (1)获取已配准多模态医学图像图2和图3,图2为烦脑CT图像，图3为MR图像，大小 均为512X512,分别将其表示为A和B。按照从左到右、从上到下的原则，对图像进行t次不同 大小的分块操作，/ e [2丄l0g_，(min(512,512))j]，即t G [2，10]。t值可W根据实际需要自由选 取，为了便于计算，本实施例中t = 5。
[004引设置分块尺寸分别为原始图像大小的1，1/4,1/8...1/2*倍。将A和B的每一分块表 示为BAi和BBi, i表示第i次分块操作，即大小为原始图像的1/21倍时的分块表示。
[0049] (2)计算分块BAi和BBi的直方图矩阵HAi、皿i，直方图矩阵反映了每一像素在不同 的分快内的概率分布。之后，根据不同分块内的像素均值MAi、MBi，计算分块内像素点（X，y) 的对比度，
[(K)加]
[0化1 ]
[0052] 按照分块大小进行排列，W直方图矩阵HAi、皿1作为系数，生成每一像素的区域对 比度矩阵。
[0053] CA(x,y) = [HA(x,y)i XCA(x,y)r--HA(x,y)iXCA(x,y)i---HA(x,y)tXCA(x,y)t]^
[0054] CB(x,y)=[皿（x,y)i XCB(x,y)i...皿（x,y)iXCB(x,y)i...皿（x,y)tXCB(x,y)t]'
[0055] (3)利用典型相关分析方法计算每一像素点对比度矩阵的相关系数。
[0056] 首先，生成待融合图像像素点的对比度矩阵向量
[0057] Z=[CA(x,y)CB(x,y)]
[0058] 计算其协方差矩阵
[0化9]
[0060]其中，乙ii = Cov(CA(x,y) ,CA(x,y))
[0061] Ei2= E2i = Cov(CA(x,y) ,CB(x,y))
[0062] E22 = Cov(CB(x,y) ,CB(x,y))
[00创其次，令1] = 4沁4山7)，￥ = 8側山7)作为两组变量的线性组合，41、81为相关系数 矩阵。则有，
[0064] Cov(U,V)=AiEi2Bi'
[0065] 最后，根据Cov(U，V)，求解得到Ai、Bi。运里求解的方法有很多种，但结果总是唯一， 本实施例应用Lagrange乘数法进行求解。
[0066] (4) W相关系数矩阵Ai、Bi作为权值，生成融合图像。
[0067] C(x,y) =wiA(x,y)+W2B(x,y)
[006引其中，wi、W2为Ai、Bi前r维系数，:T为CA(x,y)和CB(x,y)的秩，本实施例中r=l。
[0069] 【实验及结果】
[0070] 下面结合实施例对本发明融合效果进行说明。图2和图3为已完成配准的烦脑CT和 MR影像，图像大小为512X512。图4、图5、图6、图7、图8分别为利用主成分分析(PCA)、离散余 弦变换(DCT)、非下采样轮廓波变换(NSCT)、稀疏拉普拉斯变换化PSR)、边缘保持滤波(GF) 等方法融合后生成的实例图像，图9为利用本发明方法生成的融合图像。从实施例中可W看 出，图4、图6、图8融合结果各丢失了部分烦内细节信息，图5、图7烦内软组织细节模糊、对比 度低，图9具有较高的对比度，烦内软组织细节清晰。
[0071] 表1比较了各组融合方法的性能，选取弧度/角度相似性(Qabif)、联合赌化ntropy) 及运行时间（Time)作为参数，实验平台为普通电脑，酷眷i5四核处理器，8G内存，Win? 64位 操作系统，MATLAB R2014A运行环境。
[0072] _
[0073] 表1各组融合方法性能参数
[0074] 从表1中可W看出，与上述多种算法相比，本发明方法联合赌、弧度/角度相似性明 显提高，具有较好的融合效果；同时，计算时间较小，可满足实时性需求，具有较高的执行效 率。
[0075] 综上所述，上述实施方式仅为本发明较佳的【具体实施方式】，但本发明的保护范围 并不局限于此，任何熟悉本技术领域的技术人员在本发明披露的技术范围内，进行的无创 新性修改和替换，都应涵盖在本发明的保护范围之内。
【主权项】
1. 一种基于区域对比的医学图像融合方法，其特征在于包括以下步骤： a、 获取已配准的多模医学图像，根据多模医学图像的大小制定分块策略，按照从上到 下、从左到右的原则对图像进行分块操作； b、 计算每一像素在不同分块下的对比度和分块图像的直方图矩阵，并以直方图矩阵作 为系数，生成每一像素的区域对比度矩阵； c、 应用典型相关分析方法计算每一像素区域对比度矩阵的相关系数，求解像素区域对 比度矩阵最相关时的系数矩阵； d、 以求解得到的系数矩阵作为权值，生成融合图像。2. 根据权利要求1所述的基于区域对比的医学图像融合方法，其特征在于步骤a中图像分 ±夬操作分t次进行，每次分块尺寸分别为多模医学图像大小的1，1/4，1/8…l/2t倍，t为大于等于 2的正整数，定义配准后多模医学图像A和B大小均为mXn，则t的最大取值为[_l〇g2( minU，；7))」。3. 根据权利要求1或2所述的基于区域对比的医学图像融合方法，其特征在于步骤b中 采用以下方法计算每一像素区域对比度矩阵： 定义多模医学图像A和B的分块分别表示为MjPBB1, i表示第i次分块，i的取值为大于 等于2且小于等于t的正整数，分块内像素点(x，y)的对比度C(^y)1S:分块BAi和BBi的直方图矩阵，则每一像素区域对比度矩阵C(x，y)为： CA(x，y) = [HA(x，y)iXCA(x，y)i …HA(x，y)iXCA(x，y)i …HA(Xj)tXCA(Xj)t]' CB(x，y) = [HB(x，y)iXCB(x，y)i …HB(x，y)iXCB(x，y)i …HB(Xj)tXCB(Xj)t]^4. 根据权利要求3所述的基于区域对比的医学图像融合方法，其特征在于步骤c利用典 型相关分析求解区域对比度矩阵的系数矩阵的方法为： 首先，生成待融合图像像素点的对比度矩阵向量Z: Z=[CA(x,y) CB(x,y)] 计筧其协方差矩阵Σ :其次，令U=A1CA (X，y )，V = B1CB (X，y)作为两组变量的线性组合，则有， Cov(U1V)=A^i2Bi' 其中，A^B1为需要求解的系数矩阵，B1'为B1的转置矩阵。5. 根据权利要求4所述的基于区域对比的医学图像融合方法，其特征在于步骤d的融合 结果为： C(x,y)=wiA(x,y)+W2B(x,y) 其中，wi、W2为Ai、Bi前r维系数，r为CA(x,y)和CB(x,y)的秩。
【文档编号】G06T5/00GK106023126SQ201610318413
【公开日】2016年10月12日
【申请日】2016年5月12日
【发明人】曹强, 李宝生, 毛羽, 李振江
【申请人】山东省肿瘤医院