CN115578361A - 一种基于卷积稀疏表示的脑影像融合方法 - Google Patents

一种基于卷积稀疏表示的脑影像融合方法 Download PDF

Info

Publication number
CN115578361A
CN115578361A CN202211305327.0A CN202211305327A CN115578361A CN 115578361 A CN115578361 A CN 115578361A CN 202211305327 A CN202211305327 A CN 202211305327A CN 115578361 A CN115578361 A CN 115578361A
Authority
CN
China
Prior art keywords
sparse
dictionary
convolution
follows
brain image
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202211305327.0A
Other languages
English (en)
Inventor
胡秋
蔡卫明
胡绍海
白杨
邹欢清
金婧
马新莉
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhejiang University of Science and Technology ZUST
Original Assignee
Zhejiang University of Science and Technology ZUST
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Zhejiang University of Science and Technology ZUST filed Critical Zhejiang University of Science and Technology ZUST
Priority to CN202211305327.0A priority Critical patent/CN115578361A/zh
Publication of CN115578361A publication Critical patent/CN115578361A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/0002—Inspection of images, e.g. flaw detection
    • G06T7/0012—Biomedical image inspection
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20048—Transform domain processing
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20081—Training; Learning
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20212—Image combination
    • G06T2207/20221—Image fusion; Image merging
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/30—Subject of image; Context of image processing
    • G06T2207/30004—Biomedical image processing
    • G06T2207/30016—Brain

Landscapes

  • Engineering & Computer Science (AREA)
  • Health & Medical Sciences (AREA)
  • General Health & Medical Sciences (AREA)
  • Medical Informatics (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Radiology & Medical Imaging (AREA)
  • Quality & Reliability (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Image Analysis (AREA)

Abstract

本发明公开了一种基于卷积稀疏表示的脑影像融合方法，该方法将自适应卷积稀疏算法作为融合框架，将脑影像既作为训练样本得到预训练字典，又作为测试样本参与融合过程；采用非下采样轮廓波变换，即NSCT，结合子带关联性构造新的活性测度；最后，通过二尺度重构得到融合的脑影像；通过以上方法可有效改善局部失配和纹理结构丢失对卷积稀疏脑影像融合的影响，并可显著降低算法耗时。

Description

一种基于卷积稀疏表示的脑影像融合方法
技术领域
本发明涉及脑功能成像技术领域，具体而言，涉及一种基于卷积稀疏表示的脑影像融合方法。
背景技术
包括人类在内的动物行为涉及动机、感知、决策、运动控制、反馈等一系列复杂的神经环路功能，如何更好地理解动物的神经环路活动与行为之间的联系，是系统神经科学的重点目标。以鱼类动物的捕食行为为例，现有的对鱼类动物的全脑神经活动的记录都是在头部固定有检测装置的活体鱼上完成的，对脑功能结构的研究有助于提高对鱼类全脑神经活动的理解。脑功能成像处理技术中的脑影像融合，是通过将多模态的脑功能成像关于同一场景的互补信息融合到单一影像中，能够有效提升对脑功能结构信息的立体感知，从而对理解鱼类全脑神经活动与行为之间的联系提供了必要的研究基础。因多模态脑影像成像机制不同，脑影像间往往存在明显的亮暗对比、纹理差异和相似度低等特点。脑影像自身的成像特点使得已有融合方法效果并不理想，因此研究适合于多模态脑影像成像特点的融合方法非常必要。
近年来，研究者对多模态脑影像融合进行了广泛的研究。根据所采用的图像变换策略不同，此类研究主要集中四个方面：基于多尺度变换的方法、基于稀疏表示的方法、基于深度学习的方法以及三者的混合方案。不同于基于基函数假设的多尺度变换，也不同于深度特征选择依赖于人工设计规则的深度学习，稀疏表示通过对不同样本集的训练来学习获得完备字典并用于影像融合，能够实现对源影像更加稳定和灵活的表示，并且稀疏表示方式更符合人类视觉系统的生理机制。然而，传统的稀疏表示因采用分块处理方式进行影像重构，会存在多值性的问题，未能实现影像的全局优化。而卷积稀疏表示(CSR)中字典数组和稀疏系数通过卷积能够实现影像全局优化的目的，所得重构影像具有单值性的特点，既有利于融合结果中的纹理信息，又提高了对影像适配的鲁棒性。
基于CSR的脑影像融合性能主要由过完备字典和活性测度表征两个关键因素决定。其中，过完备字典的获取方式通常可分为两种，一种是通过离散余弦变换和小波变换等预先构造的固定字典，另一种是由高质量的自然图像训练得到的学习字典。已有研究表明由源影像训练得到的自适应字典在脑影像融合性能改善方面要优于固定字典。基于此，研究者提出了多种高效的字典学习算法，如联合稀疏模型、同步正交匹配追踪算法、全局字典学习算法和基于梯度信息的分类自适应字典学习算法等。除了过完备字典因素外，活性测度主要用于识别融合过程中源影像的显著特征。不合理的活性测度策略会降低稀疏系数融合的权重精度，如已有基于CSR的脑影像融合中采用L1范数的局部窗平均法构造活性测度，因L1范数采用相同权重使得稀疏显著特征不能有效突出，易造成融合结果中有效信息丢失，此外局部窗平均减弱了CSR图像变换的平移不变性，易造成融合结果在失配区域产生振铃伪影即伪吉布斯现象，并且人为控制的窗口大小也很难获得对配准理想的鲁棒性效果。基于此，结合多模态脑影像成像特点，设计一种字典冗余度低、对配准鲁棒性高的卷积稀疏脑影像融合方法是非常必要的。
发明内容
本发明解决的问题是提供一种基于卷积稀疏表示的脑影像融合方法，字典冗余度低、对配准鲁棒性高，有效改善局部失配和纹理信息丢失对基于卷积稀疏表示的脑影像融合的影响，使得融合结果保留更完备的纹理结构信息，获得满足鱼类捕食行为分析的脑功能成像。
为解决上述问题，本发明提供一种基于卷积稀疏表示的脑影像融合方法,，包括以下步骤：
S1、输入源影像Ir,r∈{1,...,R}，R表示影像样本数，r表示第r个影像样本数；
S2、将源影像Ir作为训练样本通过卷积稀疏表示算法得到预训练字典数组dm,m∈{1,...,M}，M表示字典数组中元素个数，m表示字典数组中的第m个元素，并对源影像Ir进行二尺度分解得到基层Br和若干个细节层Dr；
S3、每一个细节层Dr均对预训练字典数组dm进行稀疏编码得到稀疏系数图
Figure BDA0003905246290000021
并通过NSCT结合子带关联性构造的活性测度得到融合的稀疏系数图
Figure BDA0003905246290000022
S4、基层Br通过算术平均法得到基层融合结果
Figure BDA0003905246290000023
每一个细节层Dr均对融合的稀疏系数图
Figure BDA0003905246290000024
进行稀疏重构得到细节层融合结果
Figure BDA0003905246290000025
S5、对基层融合结果
Figure BDA0003905246290000031
以及细节层融合结果
Figure BDA0003905246290000032
进行二尺度重构得到融合脑影像IF。
进一步地，步骤S2中的卷积稀疏表示算法包括以下步骤：
S201、定义初始字典数组dm，m∈{1,...,M}与带约束的辅助变量ym，并构建满足交替方向乘子法的稀疏迭代优化代价函数如下：
Figure BDA0003905246290000033
其中，
Figure BDA0003905246290000034
表示向量的内积操作，*表示卷积操作，λ表示正则化参数，
S202、对稀疏迭代优化代价函数中的辅助变量ym以及稀疏系数sm求解，并通过交替方向乘子法对辅助变量ym和稀疏系数sm进行迭代优化以及稀疏系数更新；
S203、定义带约束的辅助变量gm，构建满足交替方向乘子法的字典迭代优化代价函数如下：
Figure BDA0003905246290000035
其中，字典范数约束集CPN＝{s∈RN:(I-PPT)s＝0||s||2＝1}，可保证优化得到的字典数组结果在空间域具有适当的约束支持，P表示补零操作，I表示单位矩阵,约束集CPN的指标函数
Figure BDA0003905246290000036
定义为：
Figure BDA0003905246290000037
S204、对字典迭代优化代价函数中的初始字典数组dm和辅助变量gm求解，并通过交替方向乘子法对dm和辅助变量gm进行迭代优化以及字典数组更新；
S205、当辅助变量ym等于稀疏系数sm且初始字典数组dm等于辅助变量gm时，结束迭代优化更新，得到预训练字典数组dm。
进一步地，步骤S3中的稀疏编码的公式如下：
Figure BDA0003905246290000038
NSCT对每一稀疏系数图
Figure BDA0003905246290000039
进行多尺度多方向分解的公式为：
Figure BDA00039052462900000310
其中，
Figure BDA00039052462900000311
表示进行K层NSCT分解操作。
Figure BDA00039052462900000312
表示低频子带部分，
Figure BDA00039052462900000313
表示具有2K个分解方向的高频子带部分。
子带关联性测度定义为：
Figure BDA0003905246290000041
其中，
Figure BDA0003905246290000042
和
Figure BDA0003905246290000043
分别表示第k和k’方向的相邻子带；
Figure BDA0003905246290000044
和
Figure BDA0003905246290000045
分别表示
Figure BDA0003905246290000046
和
Figure BDA0003905246290000047
的均值；i和j表示多尺度变换子带系数图的大小。
进一步地，步骤S2中的源影像Ir通过Tikhonov正则化方法进行二尺度分解得到基层Br和细节层Dr，公式如下：
Figure BDA0003905246290000048
其中，Tx＝[-11]和Ty＝[-11]T分别表示沿x方向和y方向的梯度算子，得到基层Br后，最终得到相应的细节层为：
Dr＝Ir-Br。
进一步地，步骤S4中的算术平均法公式为：
Figure BDA0003905246290000049
得到细节层融合结果的公式为：
Figure BDA00039052462900000410
进一步地，步骤S5中的二尺度重构即
Figure BDA00039052462900000411
进一步地，步骤S202中的稀疏系数sm求解的具体过程如下，将步骤S202中的稀疏迭代优化代价函数中稀疏系数sm的迭代优化过程表示为：
Figure BDA00039052462900000412
并化简为：
Figure BDA00039052462900000413
其中，令zm＝ym-um，定义线性操作子Dm，则Dmsm＝dm*sm，此时上式通过DFT卷积定理可转化为：
Figure BDA00039052462900000414
其中，
Figure BDA00039052462900000415
和
Figure BDA00039052462900000416
分别为Dm、sm、I和zm的离散傅里叶变换，将该式转化为集合形式为：
Figure BDA00039052462900000417
其中，定义：
Figure BDA0003905246290000051
此时，可得稀疏系数解的线性形式：
Figure BDA0003905246290000052
其中，矩阵
Figure BDA0003905246290000053
由M个串联的R×R对角矩阵组成。
进一步地，步骤S202中的辅助变量ym通过收缩软阈值法求解，公式如下：
Figure BDA0003905246290000054
其中，收缩阈值算子满足
Figure BDA0003905246290000055
⊙表示矩阵元素乘法操作。
进一步地，步骤S203中的初始字典数组dm求解的具体过程如下，将步骤S203中的字典迭代优化代价函数中的初始字典数组dm迭代优化过程表示为：
Figure BDA0003905246290000056
并化简为：
Figure BDA0003905246290000057
此时上式通过DFT卷积定理可转化为：
Figure BDA0003905246290000058
其中，满足
Figure BDA0003905246290000059
Figure BDA00039052462900000510
和
Figure BDA00039052462900000511
分别为
Figure BDA00039052462900000512
Ir和zm的离散傅里叶变换；将上式转化为集合形式为：
Figure BDA00039052462900000513
其中，定义：
Figure BDA00039052462900000514
此时，可得字典数组解的线性形式：
Figure BDA00039052462900000515
进一步地，步骤204中的辅助变量gm通过prox函数求解，公式如下：
Figure BDA00039052462900000516
且需满足：
Figure BDA0003905246290000061
本发明带来的有益效果为，采用自适应卷积稀疏表示来降低算法耗时和融合性能的不确定性；采用NSCT结合子带关联性构造新的活性测度来改善融合结果对失配的敏感性和保留更完备的有效信息。通过以上措施可有效改善局部失配和纹理结构丢失对卷积稀疏脑影像融合的影响，并可显著降低算法耗时
附图说明
图1为基于NSCT和子带关联性的自适应卷积稀疏脑影像融合系统框架图；
图2为本发明方法的流程示意图；
图3为活性测度表征对融合性能的影响示意图；
图4为局部失配对融合性能的影响示意图；
图5为对CT/MRI脑影像对采用不同融合方法的性能对比示意图；
图6为对CT/MR-T1脑影像对采用不同融合方法的性能对比示意图；
图7为对CT/MR-T2脑影像对采用不同融合方法的性能对比示意图；
图8为不同融合方法的整体客观统计结果示意图；
图9为不同融合方法在三组脑影像对下的客观性能对比示意图；
图10为不同方法融合空间分辨率为256×256的两脑影像时的CPU耗时示意图。
具体实施方式
为使本发明的上述目的、特征和优点能够更为明显易懂，下面结合附图对本发明的具体实施例做详细的说明。
需要说明的是，以下实施例假设多模态脑影像间已完成严格的配准，该假设是此类融合方案的典型情况，非本发明的特殊要求；以下定义的第一辅助变量ym和第二辅助变量gm仅是为了区分两个辅助变量，并不是表示该两个变量有某种先后关系；以下实施例将交替方向乘子法简称为ADMM算法，离散傅里叶变换卷积定理简称为DFT卷积定理。
本实施例提供一种基于卷积稀疏表示的脑影像融合方法,包括以下步骤：
输入源影像Ir,r∈{1,...,R}；其中，R表示影像样本数，r表示第r个影像样本数；
将输入的源影像Ir作为训练样本通过卷积稀疏表示算法得到预训练字典数组dm,m∈{1,...,M}；
对源影像Ir通过Tikhonov正则化方法进行二尺度分解得到基层Br和细节层Dr；
每一所述细节层Dr均通过对预训练字典数组dm进行稀疏编码得到稀疏系数图
Figure BDA0003905246290000071
并通过NSCT结合子带关联性构造的活性测度得到融合的稀疏系数图
Figure BDA0003905246290000072
基层Br通过算术平均法得到基层融合结果
Figure BDA0003905246290000073
细节层Dr通过对稀疏系数图
Figure BDA0003905246290000074
进行稀疏重构得到细节层融合结果
Figure BDA0003905246290000075
对基层融合结果
Figure BDA0003905246290000076
以及细节层融合结果
Figure BDA0003905246290000077
进行二尺度重构得到融合脑影像IF。
其中，卷积稀疏表示算法包括以下步骤：
定义初始字典数组Dm，m∈{1,...,M}与带约束的第一辅助变量，M表示字典数组中元素的个数，m表示第m个元素，并构建满足所述交替方向乘子法的稀疏迭代优化代价函数；
对稀疏迭代优化代价函数中的第一辅助变量通过收缩软阈值法求解，对稀疏系数通过离散傅里叶变换卷积定理求解，并通过ADMM算法对第一辅助变量和稀疏系数进行迭代优化以及稀疏系数更新；
定义字典数组Dm中带约束的第二辅助变量，构建满足ADMM算法的字典迭代优化代价函数；
对字典迭代优化代价函数中的字典数组Dm通过离散傅里叶变换卷积定理求解，对第二辅助变量prox函数求解，并通过ADMN算法对Dm和第二辅助变量进行迭代优化以及字典数组更新；
当时第一辅助变量等于稀疏系数时，得到预训练字典数组dm,m∈{1,...,M}。
还需要说明的是，以下为简化描述，将初始字典数组Dm与预训练字典数组dm均用dm表示，上述的各算法的具体实施方式如下：
卷积稀疏表示算法：
卷积稀疏表示算法采用一组稀疏系数图
Figure BDA0003905246290000081
和相应的字典数组dm,m∈{1,...,M}的卷积和的形式对脑影像Ir,r∈{1,...,R}进行近似，公式如下：
Figure BDA0003905246290000082
其中，
Figure BDA0003905246290000083
表示向量的内积操作，*表示卷积操作，λ表示正则化参数，M表示字典数组中元素个数，m表示字典数组中的第m个元素，R表示影像样本数，r表示第r个影像样本数。
具有稀疏惩罚项||s||1形式的稀疏编码算法，即卷积基追踪去噪问题，定义为：
Figure BDA0003905246290000084
通过引入带约束的第一辅助变量ym以采用ADMM算法迭代求解上述优化问题，此时稀疏编码的代价函数可转化为：
Figure BDA0003905246290000085
相应的变量迭代优化过程可表示为：
Figure BDA0003905246290000086
Figure BDA0003905246290000087
Figure BDA0003905246290000088
上述迭代循环过程中，变量ym可采用收缩软阈值法来求解：
Figure BDA0003905246290000089
其中，收缩阈值算子满足
Figure BDA00039052462900000810
⊙表示矩阵元素乘法操作。
稀疏系数sm可采用DFT卷积定理来求解，此时(2)式可简化为：
Figure BDA00039052462900000811
其中，令zm＝ym-um。定义线性操作子Dm，则Dmsm＝dm*sm，此时(5)式通过DFT卷积定理可转化为：
Figure BDA00039052462900000812
其中，
Figure BDA00039052462900000813
和
Figure BDA00039052462900000814
分别为Dm、sm、I和zm的离散傅里叶变换。将(6)式转化为集合形式为：
Figure BDA00039052462900000815
其中，定义：
Figure BDA0003905246290000091
此时，通过(7)式可得稀疏系数解的线性形式：
Figure BDA0003905246290000092
其中，矩阵
Figure BDA0003905246290000093
由M个串联的R×R对角矩阵组成，M表示字典数组中的元素个数，R表示图像
Figure BDA0003905246290000094
中的样本数。I表示单位矩阵。
由于(8)式中的线性形式符合Sherman-Morrison公式求解线性系统的要求，因此采用Sherman-Morrison公式快速得到迭代的稀疏系数。最后，通过离散傅里叶逆变换得到迭代的稀疏系数。(2)-(4)式迭代循环实现稀疏系数更新，直至第一辅助变量ym等于稀疏系数sm，结束迭代循环。
稀疏编码得到稀疏系数后，相应的字典更新算法为最优方向法(MOD)的卷积形式，定义为：
Figure BDA0003905246290000095
其中，字典范数约束集CPN＝{s∈RN:(I-PPT)s＝0,||s||2＝1}可保证优化得到的字典数组结果在空间域具有适当的约束支持，P表示补零操作，I表示单位矩阵。
通过引入带约束的第二辅助变量gm以采用ADMM算法迭代求解上述优化问题，此时字典更新的代价函数可转化为：
Figure BDA0003905246290000096
其中，约束集CPN的指标函数
Figure BDA0003905246290000097
定义为：
Figure BDA0003905246290000098
此时，利用ADMM算法原理得到的变量迭代优化过程可表示为：
Figure BDA0003905246290000099
Figure BDA00039052462900000910
Figure BDA00039052462900000911
上述迭代循环过程中，变量gm的求解满足几何关系：
Figure BDA00039052462900000912
则可直接采用prox函数来求解该问题，满足：
Figure BDA00039052462900000913
字典数组dm的求解可通过DFT卷积定理来求解，此时(9)式可简化为：
Figure BDA0003905246290000101
通过DFT卷积定理，(12)式可转化为：
Figure BDA0003905246290000102
其中，满足
Figure BDA0003905246290000103
和
Figure BDA0003905246290000104
分别为
Figure BDA0003905246290000105
Ir和zm的离散傅里叶变换。将(13)式转化为集合形式为：
Figure BDA0003905246290000106
其中，定义：
Figure BDA0003905246290000107
此时，通过(14)式可得字典数组解的线性形式：
Figure BDA0003905246290000108
由于(15)式中的线性形式符合Sherman-Morrison公式求解线性系统的要求，因此采用Sherman-Morrison公式快速得到迭代的字典数组。利用上述DFT卷积定理和重复Sherman-Morrison公式循环，可得到(8)式中迭代的字典数组。
通过上述稀疏编码和字典更新的循环迭代，直至第一辅助变量ym等于稀疏系数sm，结束迭代循环，获得卷积稀疏表示下的预训练字典。
影像二尺度分解：
每一源影像Ir,r∈{1,...,R}通过如下代价函数分解得到基层Br和细节层Dr：
Figure BDA0003905246290000109
其中，Tx＝[-11]和Ty＝[-11]T分别表示沿x方向和y方向的梯度算子。这是一个Tikhonov正则化问题，可通过傅里叶变换得到有效解决。得到基层Br后，可得到相应的细节层为：
Dr＝Ir-Br。
细节层Dr融合：
对于每一细节层Dr，其对应的稀疏系数图
Figure BDA00039052462900001010
通过稀疏编码中的(1)式得到。由于多尺度变换通常利用低分辨率成分来分析图像的近似特征，利用高分辨率成分来分析图像的细节特征，即通过多尺度和多方向分解得到的含有不同分辨率的分解层包含了图像的不同特征。受此启发，可将多尺度变换理解为特征分类过程，并且同一分解层的相似特征的融合有助于提高融合结果的准确性。因此，首先采用NSCT对每一稀疏系数图进行多尺度多方向分解：
Figure BDA00039052462900001011
其中，
Figure BDA0003905246290000111
表示进行K层NSCT分解操作。
Figure BDA0003905246290000112
表示低频子带部分，
Figure BDA0003905246290000113
表示具有2K个分解方向的高频子带部分。
源影像对应稀疏系数图的低频子带部分
Figure BDA0003905246290000114
采用平均策略进行融合：
Figure BDA0003905246290000115
高频子带部分
Figure BDA0003905246290000116
采用子带关联性测度(即稀疏系数图的相同尺度下相邻方向子带的相关性)来表征活性测度，其中子带关联性测度定义为：
Figure BDA0003905246290000117
其中，
Figure BDA0003905246290000118
和
Figure BDA0003905246290000119
分别表示第k和k’方向的相邻子带。
Figure BDA00039052462900001110
和
Figure BDA00039052462900001111
分别表示
Figure BDA00039052462900001112
和
Figure BDA00039052462900001113
的均值。i和j表示多尺度变换子带系数图的大小。此时，对应稀疏系数图的NSCT高频子带
Figure BDA00039052462900001114
的融合规则为：
Figure BDA00039052462900001115
在得到融合的NSCT低频子带和高频子带后，通过NSCT逆变换得到融合的稀疏系数图
Figure BDA00039052462900001116
最后，通过稀疏重构得到细节层的融合结果为：
Figure BDA00039052462900001117
基层Br融合：
由于多模态脑影像间在同一位置的亮度和纹理差异普遍较大，传统的取大融合策略可能会由于局部像素的微小变化导致视觉不一致性，因此，采用平均策略得到基层的融合结果为：
Figure BDA00039052462900001118
二尺度重构：
得到融合的基层
Figure BDA00039052462900001119
和融合的细节层
Figure BDA00039052462900001120
后，通过二尺度重构可得最终的融合影像为：
Figure BDA00039052462900001121
为了验证本发明所提方法中活性测度对融合性能的有效性，实验设计以人脑影像为例，实验设置条件如下：(1)多模态脑影像来自哈佛医学院创建的全脑图谱脑影像数据库(WBAMI)，且影像的空间分辨率均为256×256；(2)在本发明的融合框架下，基于单变量控制的思想，对比活性测度包括：子带关联性测度、结合NSCT和L1范数的测度和本发明所提结合NSCT和子带关联性的测度；(3)实验均是在电脑配置为3.3GHz CPU和16.0GHz RAM，且MATLAB版本为R2017a和Win7 64位操作系统的实验环境下进行。活性测度表征对融合性能的影响实验如附图3所示。由附图3可见，对比本发明所提活性测度和子带关联性测度所得融合结果，可知本发明所提活性测度的融合结果能够保留更完备的结构信息，因此NSCT有助于促进结构信息的保留。对比本发明所提活性测度和结合NSCT和L1范数的测度所得融合结果，可知本发明所提活性测度的融合结果能够保留更完备的纹理信息，因此子带关联性相比于L1范数有助于纹理信息的保留。综上，本发明所提融合方法中的活性测度有助于保留更完备的纹理和结构信息。
为了验证本发明所提方法的融合结果对失配的鲁棒性和保证实验的公平性，实验设计采用了与对比方法在失配的鲁棒性验证时相同的多聚焦图像对，实验设置条件如下：(1)对比方法包括：采用与本发明相同融合框架的基于子带关联性测度的方法，有代表性的基于卷积稀疏的脑影像融合方法CSMCA法和CSR法；(2)多聚焦图像对的空间分辨率均为256×256；(3)实验均是在电脑配置为3.3GHz CPU和16.0GHz RAM，且MATLAB版本为R2017a和Win7 64位操作系统的实验环境下进行。局部失配对脑影像融合性能的影响实验如附图4所示。由附图4可见，CSMCA法和CSR法的融合结果在失配区域边缘均存在明显的振铃伪影，即伪吉布斯现象。这是因为CSMCA和CSR采用的局部窗平均法减弱了图像变换的移不变性，这使二者所得融合结果对失配的敏感性降低，并且人为控制的窗口大小也很难获得对配准理想的鲁棒性效果，而在基于子带关联测度的方法和本发明所提方法中，子带关联性采用全局子带保证了移不变性，可有效消除伪吉布斯现象。
为了验证本发明所提融合方法的合理性，实验设计以人脑影像为例，实验条件的设置包括：(1)从哈佛医学院创建的全脑图谱脑影像数据库的“Acute stroke”数据集、“Hypertensive encephalopathy”数据集和“Multiple embolic infarctions”数据集中随机选取了30组多模态脑影像，且影像的空间分辨率均为256×256；(2)对比方法包括：采用与本发明相同融合框架的基于子带关联性测度的方法，有代表性的基于卷积稀疏的脑影像融合方法CSMCA法和CSR法；(3)算法参数设置：每一字典数组的空间大小为8×8，字典数组设置为32，正则化参数λ设置为0.01，上述参数设置均与CSR法和CSMCA法的设置相同。(5)客观评价指标包括：熵EN、联合熵QU、互信息MI、结构相似度SSIM，上述客观指标的值越大表明融合性能越好；(6)实验均是在电脑配置为3.3GHz CPU和16.0GHz RAM，且MATLAB版本为R2017a和Win7 64位操作系统的实验环境下进行。此外，采用主观评价和客观评价相结合的方式进行融合性能分析。其中，主观评价的融合结果如附图5、附图6和附图7所示，客观评价的性能指标如附图9(加粗字体表示最优指标，下划线表示次优指标)、附图8和附图10所示。由附图5、附图6和附图7可见，CSMCA法和CSR方的融合结果在纹理清晰度和亮度信息方面保留更完备，而边缘纹理的连续性遭到破坏，但均存在局部细节丢失的情况。基于子带关联性测度的方法在细节纹理的清晰度和结构完整性方面都不占有优势，但在边缘结构纹理方面保留更完备。相比之下，本发明所提融合方法保留了更完备的结构和纹理信息。由附图9、图8和附图10可见，本发明所提方法在算法耗时偏低情况下，四个客观指标在三个脑影像数据集实验中得分均比较突出。其中，本发明所提方法的MI指标优势最明显，表明本发明所提方法既能够保留更完备的纹理结构信息，又能够保持较高的算法执行效率。
虽然本公开披露如上，但本公开的保护范围并非仅限于此。本领域技术人员，在不脱离本公开的精神和范围的前提下，可进行各种变更与修改，这些变更与修改均将落入本发明的保护范围。

Claims (10)

1.一种基于卷积稀疏表示的脑影像融合方法,其特征在于，包括以下步骤：
S1、输入源影像Ir,r∈{1,...,R}，R表示影像样本数，r表示第r个影像样本数；
S2、将所述源影像Ir作为训练样本通过卷积稀疏表示算法得到预训练字典数组dm,m∈{1,...,M}，M表示字典数组中元素个数，m表示字典数组中的第m个元素，并对源影像Ir进行二尺度分解得到基层Br和若干个细节层Dr；
S3、每一个所述细节层Dr均通过所述预训练字典数组dm进行稀疏编码得到稀疏系数图
Figure FDA0003905246280000011
并通过NSCT结合子带关联性构造的活性测度得到融合的稀疏系数图
Figure FDA0003905246280000012
S4、所述基层Br通过算术平均法得到基层融合结果
Figure FDA0003905246280000013
每一个所述细节层Dr均对所述融合的稀疏系数图
Figure FDA0003905246280000014
进行稀疏重构得到细节层融合结果
Figure FDA0003905246280000015
S5、对所述基层融合结果
Figure FDA0003905246280000016
以及所述细节层融合结果
Figure FDA0003905246280000017
进行二尺度重构得到融合脑影像IF。
2.根据权利要求1所述的一种基于卷积稀疏表示的脑影像融合方法,其特征在于，
所述步骤S2中的所述卷积稀疏表示算法包括以下步骤：
S201、定义初始字典数组dm，m∈{1,...,M}与带约束的辅助变量ym，并构建满足交替方向乘子法的稀疏迭代优化代价函数如下：
Figure FDA0003905246280000018
其中，
Figure FDA0003905246280000019
表示向量的内积操作，*表示卷积操作，λ表示正则化参数，I即源影像Ir，sm为稀疏系数；
S202、对所述稀疏迭代优化代价函数中的所述辅助变量ym以及稀疏系数sm求解，并通过所述交替方向乘子法对所述辅助变量ym和所述稀疏系数sm进行迭代优化以及稀疏系数更新；
S203、定义带约束的辅助变量gm，构建满足所述交替方向乘子法的字典迭代优化代价函数如下：
Figure FDA00039052462800000110
其中，字典范数约束集CPN＝{s∈RN:(I-PPT)s＝0,||s||2＝1}，可保证优化得到的字典数组结果在空间域具有适当的约束支持，P表示补零操作，I表示单位矩阵,s表示字典数组dm的范数，和所述辅助变量求解，并通过所述交替方向乘子法对和所述辅助变量进行迭代优约束集CPN的指标函数
Figure FDA0003905246280000021
定义为：
Figure FDA0003905246280000022
S204、对所述字典迭代优化代价函数中的所述初始字典数组dm和所述辅助变量gm求解，并通过所述交替方向乘子法对dm和所述辅助变量gm进行迭代优化以及字典数组更新；
S205、当所述辅助变量ym等于所述稀疏系数sm且所述初始字典数组dm等于所述辅助变量gm时，结束迭代优化更新，得到所述预训练字典数组dm。
3.根据权利要求1所述的一种基于卷积稀疏表示的脑影像融合方法,其特征在于，
所述步骤S3中的稀疏编码的公式如下：
Figure FDA0003905246280000023
NSCT对每一稀疏系数图
Figure FDA0003905246280000024
进行多尺度多方向分解的公式为：
Figure FDA0003905246280000025
其中，
Figure FDA0003905246280000026
表示进行K层NSCT分解操作，
Figure FDA0003905246280000027
表示低频子带部分，
Figure FDA0003905246280000028
表示具有2K个分解方向的高频子带部分；
子带关联性测度定义为：
Figure FDA0003905246280000029
其中
Figure FDA00039052462800000210
和
Figure FDA00039052462800000211
分别表示第k和k’方向的相邻子带；
Figure FDA00039052462800000212
和
Figure FDA00039052462800000213
分别表示
Figure FDA00039052462800000214
和
Figure FDA00039052462800000215
的均值；i和j表示多尺度变换子带系数图的大小。
4.根据权利要求1所述的一种基于卷积稀疏表示的脑影像融合方法,其特征在于，
所述步骤S2中的所述源影像Ir通过Tikhonov正则化方法进行二尺度分解得到基层Br和细节层Dr，公式如下：
Figure FDA00039052462800000216
其中，Tx＝[-1 1]和Ty＝[-1 1]T分别表示沿x方向和y方向的梯度算子，得到基层Br后，最终得到相应的细节层为：
Dr＝Ir-Br。
5.根据权利要求1所述的一种基于卷积稀疏表示的脑影像融合方法,其特征在于，
所述步骤S4中的算术平均法公式为：
Figure FDA0003905246280000031
得到细节层融合结果的公式为：
Figure FDA0003905246280000032
6.根据权利要求1所述的一种基于卷积稀疏表示的脑影像融合方法,其特征在于，
所述步骤S5中的二尺度重构即
Figure FDA0003905246280000033
7.根据权利要求2所述的一种基于卷积稀疏表示的脑影像融合方法,其特征在于，
所述步骤S202中的所述稀疏系数Sm求解的具体过程如下，将所述步骤S202中的稀疏迭代优化代价函数中所述稀疏系数sm的迭代优化过程表示为：
Figure FDA0003905246280000034
并化简为：
Figure FDA0003905246280000035
其中，令zm＝ym-um，定义线性操作子Dm，则Dmsm＝dm*sm，此时上式通过DFT卷积定理可转化为：
Figure FDA0003905246280000036
其中，
Figure FDA0003905246280000037
和
Figure FDA0003905246280000038
分别为Dm、sm、I和zm的离散傅里叶变换，将该式转化为集合形式为：
Figure FDA0003905246280000039
其中，定义：
Figure FDA00039052462800000310
此时，可得稀疏系数解的线性形式：
Figure FDA00039052462800000311
其中，矩阵
Figure FDA00039052462800000312
由M个串联的R×R对角矩阵组成。
8.根据权利要求2所述的一种基于卷积稀疏表示的脑影像融合方法,其特征在于，
所述步骤S202中的所述辅助变量ym通过收缩软阈值法求解，公式如下：
Figure FDA0003905246280000041
其中，收缩阈值算子满足
Figure FDA0003905246280000042
⊙表示矩阵元素乘法操作。
9.根据权利要求2所述的一种基于卷积稀疏表示的脑影像融合方法,其特征在于，
所述步骤S203中的初始字典数组dm求解的具体过程如下，将所述步骤S203中的字典迭代优化代价函数中的初始字典数组dm迭代优化过程表示为：
Figure FDA0003905246280000043
并化简为：
Figure FDA0003905246280000044
此时上式通过DFT卷积定理可转化为：
Figure FDA0003905246280000045
其中，满足
Figure FDA0003905246280000046
Figure FDA0003905246280000047
和
Figure FDA0003905246280000048
分别为
Figure FDA0003905246280000049
Ir和zm的离散傅里叶变换；将上式转化为集合形式为：
Figure FDA00039052462800000410
其中，定义：
Figure FDA00039052462800000411
此时，可得字典数组解的线性形式：
Figure FDA00039052462800000412
10.根据权利要求2所述的一种基于卷积稀疏表示的脑影像融合方法,其特征在于，
所述步骤S204中的所述辅助变量gm通过prox函数求解，公式如下：
Figure FDA00039052462800000413
且需满足：
Figure FDA00039052462800000414
CN202211305327.0A 2022-10-24 2022-10-24 一种基于卷积稀疏表示的脑影像融合方法 Pending CN115578361A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211305327.0A CN115578361A (zh) 2022-10-24 2022-10-24 一种基于卷积稀疏表示的脑影像融合方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211305327.0A CN115578361A (zh) 2022-10-24 2022-10-24 一种基于卷积稀疏表示的脑影像融合方法

Publications (1)

Publication Number Publication Date
CN115578361A true CN115578361A (zh) 2023-01-06

Family

ID=84586874

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211305327.0A Pending CN115578361A (zh) 2022-10-24 2022-10-24 一种基于卷积稀疏表示的脑影像融合方法

Country Status (1)

Country Link
CN (1) CN115578361A (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116910450A (zh) * 2023-07-18 2023-10-20 南昌大学 一种基于卷积稀疏编码的心电信号去噪方法

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116910450A (zh) * 2023-07-18 2023-10-20 南昌大学 一种基于卷积稀疏编码的心电信号去噪方法

Similar Documents

Publication Publication Date Title
Xiang et al. FISTA-Net: Learning a fast iterative shrinkage thresholding network for inverse problems in imaging
Hou et al. Brain CT and MRI medical image fusion using convolutional neural networks and a dual-channel spiking cortical model
Ben Yedder et al. Deep learning for biomedical image reconstruction: A survey
Zhang et al. Reducing uncertainty in undersampled MRI reconstruction with active acquisition
Wang et al. Single slice based detection for Alzheimer’s disease via wavelet entropy and multilayer perceptron trained by biogeography-based optimization
CN109325985B (zh) 磁共振图像重建方法、装置和计算机可读存储介质
Singh et al. Multimodal medical image sensor fusion model using sparse K-SVD dictionary learning in nonsubsampled shearlet domain
Ravishankar et al. Data-driven learning of a union of sparsifying transforms model for blind compressed sensing
CN106056647B (zh) 一种基于卷积稀疏双层迭代学习的磁共振快速成像方法
Yousif et al. An improved approach for medical image fusion using sparse representation and Siamese convolutional neural network
Rajalingam et al. A novel approach for multimodal medical image fusion using hybrid fusion algorithms for disease analysis
Zhang et al. Analysis-synthesis dictionary pair learning and patch saliency measure for image fusion
Zhou et al. Volume upscaling with convolutional neural networks
CN111487573B (zh) 一种用于磁共振欠采样成像的强化型残差级联网络模型
Rajalingam et al. Review of multimodality medical image fusion using combined transform techniques for clinical application
Li et al. DFENet: A dual-branch feature enhanced network integrating transformers and convolutional feature learning for multimodal medical image fusion
Zhang et al. A separation–aggregation network for image denoising
Wahid et al. Multi-layer basis pursuit for compressed sensing MR image reconstruction
Das et al. NSST domain CT–MR neurological image fusion using optimised biologically inspired neural network
Rajalingam et al. Multimodal medical image fusion using various hybrid fusion techniques for clinical treatment analysis
Ye et al. Unified supervised-unsupervised (super) learning for x-ray ct image reconstruction
CN115578361A (zh) 一种基于卷积稀疏表示的脑影像融合方法
CN112927235B (zh) 基于多尺度超像素与核低秩表示的脑肿瘤图像分割方法
Ramírez et al. Analysis of empirical Bayesian methods for neuroelectromagnetic source localization
Pasha et al. A computational framework for edge-preserving regularization in dynamic inverse problems

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination