CN106886986A - 基于自适应组结构稀疏字典学习的图像融合方法 - Google Patents

Abstract

本发明公开了一种基于稀疏表示自适应学习组结构化字典的图像融合方法，具体为：从输入待融合图像并进行滑窗分块，提取特征向量得到训练样本训练出自适应的组结构字典，对每个图像块计算稀疏表示向量，取相同位置对应的稀疏表示向量采用基于组结构的L1范数最大化的方法得到融合后的系数表示向量矩阵，最后将其加上相应的均值转换为图像表示，得到图像块，使用滑窗的逆操作输出最终融合图像。本发明的融合后的图像更加精确、细节更加丰富，能够减少融合图像空间上的不连续性。

Description

基于自适应组结构稀疏字典学习的图像融合方法

技术领域

本发明属于图像处理技术领域，特别涉及图像的融合，可用于图像处理、计算机视觉、医学图像处理及军事等领域。

背景技术

近年来，稀疏表示作为一种有效的方法受到了越来越多的关注，它被应用到了图像处理的许多问题中，例如：图像去噪、图像融合、图像压缩等。一幅场景中含有丰富的信息，用单成像传感器系统不能捕获场景中所有的信息，为了解决此问题，多传感器图像融合系统应运而生，多传感器图像融合系统有效利用不同成像传感器之间的互补性，消除多传感器图像之间的冗余信息，对不同成像传感器获取的图像进行综合，形成对目标更加完整、清晰、准确的描述，从而有效地增加对象的信息量、提高图像信息的利用率、提高系统的精度和鲁棒性、降低系统的不确定性、扩大系统的应用范围、多传感器图像融合系统在军事、医学、遥感和公共安全等领域具有广泛的应用前景。

基于稀疏表示的图像融合主要分为三个阶段：构造字典、稀疏编码和对源图像的稀疏表示系数设计融合规则。基于稀疏表示的图像融合的框架中，字典和融合规则的设计对图像融合的质量有着非常重要的影响。

当前主要构造字典的方法是通过样本学习得到字典。字典的学习方法有最优方向法即MOD(Method of Optimal Directions),K-SVD(Singular Value Decomposition)和任务驱动(Task-Driven)的字典学习方法的基础。其中K-SVD算法是最常用的字典学习方法，也是其他许多字典学习方法的基础。但是K-SVD算法构造字典的缺点：单个字典的表示能力有限，对图像的细节的表示能力不足，这样会导致融合后的图像会有所平滑，对细节的重构能力比较差。自然图像中有丰富的结构信息，而K-SVD的字典学习方法没有利用图像中的结构信息。在文献(Eldar,Y.C.,P.Kuppinger and H.Bolcskei,Block-Sparse Signals:Uncertainty Relations and Efficient Recovery.Signal Processing,IEEETransactions on,2010.58(6):p.3042-3054.)中提出了组结构字典，该文献中证明了组结构字典重构信号的理论可行性并且说明了组结构字典考虑了信号中的结构信息后在重构信号的能力上表现更好。文献(Shutao,L.,Y.Haitao and F.Leyuan,Group-SparseRepresentation With Dictionary Learning for Medical Image Denoising andFusion.Biomedical Engineering,IEEE Transactions on,2012.59(12):p.3450-3459.)中利用了组结构化字典进行图像融合。但是这里字典的组结构是预先固定了的，获得的组结构字典不能够自适应的根据信号调整自己的分组，在该文章中约束组内字典原子是通过增加字典中的相似性，而这样会导致同一组的字典过于相似，使一个组的字典原子的表示能力受限。以上的原因会导致获得的字典在表示信号时的能力不够好，从而导致融合的结果不够好。

基于稀疏表示的融合规则的设计，一般是L1范数最大或者取平均。取平均的融合规则会使融合后的图像的对比度下降，还会有比较严重的平滑现象；L1范数最大的融合规则是根据L1范数大的那个图像块被选中作为融合图像中对应的图像块，而L1范数较小的对应的整个图像块就被抛弃了，这样会导致以下两个问题：

(1)L1范数较小的那个图像块中也会有一些很显著的信息需要保留，而这样的融合规则会使在融合图像中并没有得到这一部分的信息，使融合图像的信息损失较大、细节不够丰富。

(2)L1范数最大的规则会造成在融合图像中空间上的不连续性。这个现象在异质图像的融合中会经常发生。

发明内容

本发明的发明目的在于：针对上述存在的问题，提出一种基于稀疏表示自适应学习组结构化字典的图像融合方法。本发明学习到的字典是自适应的组结构化字典，这种字典的表示能力比较强，对细节的重构能力比较好。在融合规则上充分利用结构化字典的信息，使融合的图像不会失去源图像中的显著信息，并且减少了融合图像空间上的不连续性。本发明的具体实现步骤如下：

1、基于自适应组结构稀疏字典学习的图像融合方法，其特征在于，包括下列步骤：

步骤1：对待融合图像A和B分别进行滑窗分块处理，得到t个大小为的局部图像块其中图像块标识符i＝{1,2,,…,t},N为预设值；

所述滑窗分块处理为：采用滑窗步长为1，窗口大小为的滑动窗口自左到右从上到下进行重叠滑窗；

步骤2：把图像块的个像素值排成一列，得到列向量i＝{1,2,3,…,t}，列向量减去各自对应的列向量的均值得到零均值的列向量即其中Q是N×1的全为1的列向量；

步骤3：判断图像块数目t是否大于阈值T，若是，则分别从t个中随机选择T个列向量作为训练样本集Y；否则，直接将t个作为训练样本集Y；用参数L表示训练样本集Y的个数；

步骤4：随机初始化一个大小为N×M的字典D₀，其中N等于步骤1所述的图像块的边长的平方，M为预设值且M>N；

将字典D₀作为组结构字典D初始值，基于字典D₀、训练样本集Y，采用自适应组结构化字典学习方法学习得到组结构字典D，在字典学习时，约束组结构字典D的组稀疏度为预设值k，且每个组的最大字典原子数为预设值s；

步骤5：基于步骤4得到的组结构字典D，分别对各采用组稀疏度为k的块正交匹配追踪算法(B-OMP)计算稀疏系数向量

步骤6：对稀疏系数向量采用基于组结构的L1范数最大的方法得到对应的融合稀疏系数向量X_i，i＝{1,2,3,…,t}；

步骤7：根据公式计算第一融合图像列向量即第一融合图像列向量为没有源图像亮度信息的融合图像列向量；

根据列向量的均值的均值得到融合列向量均值即再根据公式得到第二融合图像列向量即二融合图像列向量为含有源图像对应亮度信息的融合图像列向量；

基于步骤2中得到列向量的逆方式，将第二融合图像列向量转化成大小为的融合图像块

步骤8：将每个融合图像块按照步骤1滑窗分块时的顺序进行排列，对存在融合图像块重叠的像素点的像素值进行累加、再除以各像素点的累加次数，得到的融合图像F。

进一步的，步骤4所述的自适应组结构化字典学习方法得到组结构字典D的具体步骤为：

401：采用正交匹配追踪算法(OMP)计算训练样本集Y的每个样本基于当前组结构字典D的稀疏表示系数矩阵X，所述稀疏表示系数矩阵X的稀疏度为s*k，矩阵大小为M×L，即X＝argmin_X‖Y-DX‖_F且||X_j||₀≤s*k，其中j＝1,…,L，||X_j||₀表系数矩阵X的第j列不为零的元素的数目，符号||·||_F表示Frobenius范数，组结构字典D的初始值为D₀；

402：对稀疏表示系数矩阵X进行更新，先对矩阵X按行进行标准化，再对矩阵按列进行标准化后得到更新后的稀疏表示系数矩阵X；

403：基于当前稀疏表示系数矩阵X计算组结构字典D的分组信息d：

403-1：初始化矩阵X′＝X，分组向量d为稀疏表示系数矩阵X的行索引值，即d＝[1,2,3,…,M]；

403-2：计算矩阵X′中任意两行的内积，并内积大小进行排序，初始化i＝1，并执行步骤403-2-1：

403-2-1：判断内积为第i大的两行的原子数目之和是否大于s，若是，则执行步骤403-2-2，否则执行步骤403-3；

403-2-2：若i<L则令i＝i+1后，继续执行步骤403-2-1；若i＝L则分组聚类算法结束，执行步骤404；

403-3：将内积为第i大的两行合并，得到更新后的矩阵X′，其中合并的方式为：两个行向量相加除以2；

并基于内积为第i大的两行的行索引值对分组向量d进行更新后，继续执行步骤403-2；

其中对分组向量d的更新具体为：用p、q表示内积为第i大的行索引值，将分组向量d中的行索引值p、q全部修改为p与q中的较小索引值；

404：将分组信息d中同一取值的列作为一组，得到组结构字典D的分组方式，并基于组结构字典D的分组结果更新稀疏表示系数矩阵X，即X＝argmin_X‖Y-DX‖_F且s.t.||X_j||_0,d≤k,j＝1,…L，其中||X_j||_0,d表示稀疏表示系数矩阵的第j列不为0的元素对应的分组数目；405：逐组更新组结构字典D；

406：判断是否满足迭代收敛条件，若是，则基于当前组结构字典D执行步骤5；否则继续执行步骤401；

其中迭代收敛条件为：达到最大迭代次数或最近两次迭代得到的组结构字典D的重构信号的误差小于或等于预设阈值。

进一步的，步骤6采用基于组结构的L1范数最大的方法进行融合的具体为：

根据权利要求1的步骤4学习得到组结构字典D，基于组结构的L1范数最大的融合规则为：

若则否则其中j＝1…M,k＝1…B，B表示分组向量d中元素的最大值，分别表示稀疏系数向量的第j列中属于第k组元素的绝对值之和，X_i,d(j)＝k表示融合系数向量X_i的第j列中属于第k组的元素。

综上所述，由于采用了上述技术方案，本发明的有益效果是：

(1)本发明采用自适应的组结构字典，使得图像稀疏表示系数更为稀疏，因此融合后的图像更加精确，且融合处理速度快；

(2)本发明采用基于组结构的L1范数最大的融合规则，能够减少融合图像中的信息损失，使融合图像中的细节更加丰富；并且能够减少融合图像空间上的不连续性。

附图说明

图1是本发明的流程图

图2是待融合的红外和可见光图像

图3是用本发明的融合方法和基于K-SVD的融合方法对图2的融合结果

图4是待融合的脑部医学图像

图5是用本发明的融合方法和基于K-SVD的融合对图4的融合结果。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面结合实施方式和附图，对本发明作进一步地详细描述。

参照图1，本发明的具体实施过程如下：

(1)从输入的待融合图像A和B(图像A、B的尺寸相同)，分别对图像A和B执行滑窗分块操作：采用滑窗步长为1的滑动窗口为自左到右、从上到下进行滑窗分块，从而分别得到t个大小为的图像块其中：图像块标识符i＝{1,2,,…,t}，N根据待融合图像的大小选定，本具体实施方式中，选取N＝64。

(2)把图像块分别拉成列向量即将每个图像块的个像素分别排列成一列，生成列向量列向量减去各自对应的列向量的均值得到零均值的列向量即其中Q是N×1的全为1的列向量。

(3)若t>T,则分别从集合θ^A、θ^B中随机选择T个列向量作为训练样本集Y，否则直接用θ^A、θ^B作为训练样本集Y，令Y中样本的个数为L，其中T为预设阈值，本具体实施方式中，取T＝10000；

(4)随机初始化一个大小为N×M(M>N)的字典D₀，将字典D₀作为组结构字典D初始值，在训练样本集Y上采用自适应组结构化字典学习方法学习得到组结构字典D，在字典学习时，约束字典组结构D的组稀疏度为k＝4，最大的组能含有的字典原子数为s＝4：

401：采用正交匹配追踪算法计算训练样本集Y的每个样本基于当前组结构字典D的稀疏表示系数矩阵X，所述稀疏表示系数矩阵X的稀疏度为s*k，矩阵大小为M×L，即X＝argmin_X‖Y-DX‖_F且||X_j||₀≤s*k，其中j＝1,…,L，||X_j||₀表系数矩阵X的第j列不为零的元素的数目，符号||·||_F表示Frobenius范数，组结构字典D的初始值为D₀；

402：对稀疏表示系数矩阵X进行更新，先对矩阵X按行进行标准化，再对矩阵按列进行标准化后得到更新后的稀疏表示系数矩阵X，本具体实施方式中采用的标准化方式是Z分数标准化(Z-score标准化)；

403：基于当前稀疏表示系数矩阵X计算组结构字典D的分组信息d：

403-1：初始化矩阵X′＝X，分组向量d为稀疏表示系数矩阵X的行索引值，即d＝[1,2,3,…,M]；

403-2：计算矩阵X′中任意两行的内积，并内积大小进行排序，初始化i＝1，并执行步骤403-2-1：

403-2-1：判断内积为第i大的两行的原子数目之和是否大于s，若是，则执行步骤403-2-2，否则执行步骤403-3；

403-2-2：若i<L则令i＝i+1后，继续执行步骤403-2-1；若i＝L则分组聚类算法结束，执行步骤404；

403-3：将内积为第i大的两行合并，得到更新后的矩阵X′，其中合并的方式为：两个行向量相加除以2；

并基于内积为第i大的两行的行索引值对分组向量d进行更新后，继续执行步骤403-2；

其中对分组向量d的更新具体为：用p、q表示内积为第i大的行索引值，将分组向量d中的行索引值p、q全部修改为p与q中的较小索引值；

404：将分组信息d中同一取值的列作为一组，得到组结构字典D的分组方式，以M＝8为例，若分组信息d为[1,2,1,1,2,3,4,5]，则表示在组结构字典D中，第1、3、4列属于一组、第2、5列属于一组，第6、7、8列分别属于一组。

然后，基于组结构字典D的分组结果更新稀疏表示系数矩阵X，即X＝argmin_X‖Y-DX‖_F且s.t.||X_j||_0,d≤k,j＝1,…L，其中||X_j||_0,d表示稀疏表示系数矩阵的第j列不为0的元素对应的分组数目；

405：逐组更新组结构字典D：

组结构字典D的组数用B表示，即分组向量d中元素的最大值，则保持B-1组的字典列不变，更新第j₀组的字典列，其中j₀表示当前待更新组，用公式可以表示为：

其中，是表示稀疏表示系数矩阵X的第j行，符号(·)^T表示矩阵转置，D_j表示组结构字典D的第j列，d(j)表示分组向量d中第j个元素的值。

为了简化优化处理过程，定义参量则优化属于第j₀组的字典和对应的稀疏表示系数为：用m表示第j₀组含有的原子数，该优化问题等价于求秩为m的对近似的矩阵，可以使用SVD算法求解完成对组结构字典D的逐组更新。

406：判断是否满足迭代收敛条件，若是，则基于当前组结构字典D执行步骤5；否则继续执行步骤401；

其中迭代收敛条件为：达到最大迭代次数或最近两次迭代得到的组结构字典D的重构信号的误差小于或等于预设阈值。

(5)分别对用训练好的组结构字典D(步骤(4)学习得到的组结构字典)采用组稀疏为k的B-OMP算法求出稀疏表示的系数向量

(6)对系数向量采用基于组结构的L1范数最大的方法得到对应的融合系数向量X_i，i＝{1,2,,…,t}，其中融合规则具体为：若则否则其中j＝1…M,k＝1…B，B表示分组向量d中元素的最大值，分别表示稀疏系数向量的第j列中属于第k组元素的绝对值之和，X_i,d(j)＝k表示融合系数向量X_i的第j列中属于第k组的元素。

(7)计算融合后的图像列向量根据列向量的均值的均值得到融合列向量均值再根据公式得到含有源图像对应亮度信息的融合图像列向量

基于步骤(1)中生成列向量的逆方式，将第二融合图像列向量转化成大小为的融合图像块

(8)将每个融合图像块按照步骤1滑窗分块时的顺序进行排列，对存在融合图像块重叠的像素点的像素值进行累加、再除以各像素点的累加次数，得到的融合图像F。

本发明的技术效果由以下仿真实验进一步说明。

1、仿真条件：

1)选取4幅待融合图像，如图2、3所示：

2)仿真实验的对比方法使用传统的基于K-SVD的图像融合方法；

3)仿真实验中本发明的组稀疏度k＝4，组内最大原子数s＝4；传统的基于K-SVD的图像融合方法的组稀疏度取值为16。

2、仿真内容和仿真结果

仿真一，用本发明方法与现有的小波方法对图2进行融合的仿真实验，结果如图4所示，其中：图4(a)是本发明得到的结果，图4(b)是传统的基于K-SVD的融合方法得到的结果图。

从图4中可以看出，本发明融合出的图像包含的信息更多，融合的结果也比对比的方法更加清晰，特别是对源图像的一些显著的小特征融合的非常好，且也可以看出本发明的融合图像减少了图像在空间上的不连续性。

将用本发明方法和现有的基于K-SVD的融合方法对图2所示的两幅待融合图像进行融合的性能参数进行对比如表1下所示：

表1

其中SSIM(structural similarity index measure)是结构相似度的评价准则，是图2(a)和各自方法融合图像的结构相似度，是图2(b)和各自方法融合图像的结构相似度，是图2(a)和图2(b)和各自融合方法的平均结构相似度。CORR(correlation)是两幅图像的相关系数，是图2(a)和各自方法融合图像的相关系数，是图2(b)和各自方法融合图像的相关系数，是图2(a)和图2(b)和各自融合方法的平均结构相似度。

从表1中可以看出，由结构相似度的评价指标可以发现：对比方法的融合结果会倾向于只选择某一副源图像中的结构信息，导致另一幅源图像的一些显著结构在源图像中得不到保留，而本发明方法在选择结构信息的时候是综合考虑各自源图像中结构的显著性来进行融合的，这说明了本发明方法的基于组结构的L1范数最大的融合规则能够从源图像中获取的结构信息更为丰富，从而使融合图像的细节表现力更强。表1中相关系数的评价指标也反映了本发明方法比对比方法能更均衡的融合各自图像中的内容。

仿真二，用本发明方法与现有的基于K-SVD的融合方法对图3所示的两幅待融合图像进行融合的仿真实验，结果如图5所示，其中:图5(a)是本发明的融合结果图，图5(b)是基于现有K-SVD融合方法得到的融合结果图。

从图5可以看出，本发明的融合结果信息更为丰富，对各自源图像中的特征保留的更好。

表2用本发明方法与基于现有K-SVD融合方法进行融合的性能参数对比

表2中各项的定义与表1的相同。

从表2中可以看出，本发明的图像融合结果比对比方法的与源图像的结构相似度更高，说明本发明方法从源图像中获取的结构信息更加丰富。本发明方法的与源图像的相关系数更高，说明本发明方法的融合结果中的内容和源图像中的内容相似度更高，说明获取了源图像中的更多信息。

以上所述，仅为本发明的具体实施方式，本说明书中所公开的任一特征，除非特别叙述，均可被其他等效或具有类似目的的替代特征加以替换；所公开的所有特征、或所有方法或过程中的步骤，除了互相排斥的特征和/或步骤以外，均可以任何方式组合。

Claims (3)

1.基于自适应组结构稀疏字典学习的图像融合方法，其特征在于，包括下列步骤：

步骤1：对待融合图像A和B分别进行滑窗分块处理，得到t个大小为的局部图像块其中图像块标识符i＝{1,2,,…,t},N为预设值；

所述滑窗分块处理为：采用滑窗步长为1，窗口大小为的滑动窗口自左到右从上到下进行重叠滑窗；

步骤2：把图像块的个像素值排成一列，得到列向量i＝{1,2,3,…,t}，列向量减去各自对应的列向量的均值得到零均值的列向量i＝{1,2,,…,t}，即其中Q是N×1的全为1的列向量；

步骤3：判断图像块数目t是否大于阈值T，若是，则分别从t个中随机选择T个列向量作为训练样本集Y；否则，直接将t个作为训练样本集Y；用参数L表示训练样本集Y的个数；

步骤4：随机初始化一个大小为N×M的字典D₀，其中N等于步骤1所述的图像块的边长的平方，M为预设值且M>N；

将字典D₀作为组结构字典D初始值，基于字典D₀、训练样本集Y，采用自适应组结构化字典学习方法学习得到组结构字典D，在字典学习时，约束组结构字典D的组稀疏度为预设值k，且每个组的最大字典原子数为预设值s；

步骤5：基于步骤4得到的组结构字典D，分别对各采用组稀疏度为k的块正交匹配追踪算法计算稀疏系数向量i＝{1,2,3,…,t}；

步骤6：对稀疏系数向量采用基于组结构的L1范数最大的方法得到对应的融合稀疏系数向量X_i，i＝{1,2,3,…,t}；

步骤7：根据公式计算第一融合图像列向量i＝{1,2,3,…,t}；

根据列向量的均值的均值得到融合列向量均值即再根据公式得到第二融合图像列向量i＝{1,2,3,…,t}；

基于步骤2中得到列向量的逆方式，将第二融合图像列向量转化成大小为的融合图像块i＝{1,2,3,…,t}；

步骤8：将每个融合图像块i＝{1,2,3,…,t}按照步骤1滑窗分块时的顺序进行排列，对存在融合图像块重叠的像素点的像素值进行累加、再除以各像素点的累加次数，得到的融合图像F。

2.根据1所述的方法，其中步骤4所述的自适应组结构化字典学习方法得到组结构字典D的具体步骤为：

401：采用正交匹配追踪算法计算训练样本集Y的每个样本基于当前组结构字典D的稀疏表示系数矩阵X，所述稀疏表示系数矩阵X的稀疏度为s*k，矩阵大小为M×L，即X＝argmin_X‖Y-DX‖_F且||X_j||₀≤s*k，其中j＝1,…,L，||X_j||₀表系数矩阵X的第j列不为零的元素的数目，符号||·||_F表示Frobenius范数，组结构字典D的初始值为D₀；

402：对稀疏表示系数矩阵X进行更新，先对矩阵X按行进行标准化，再对矩阵按列进行标准化后得到更新后的稀疏表示系数矩阵X；

403：基于当前稀疏表示系数矩阵X计算组结构字典D的分组信息d：

403-1：初始化矩阵X′＝X，分组向量d为稀疏表示系数矩阵X的行索引值，即d＝[1,2,3,…,M]；

403-2：计算矩阵X′中任意两行的内积，并根据内积大小进行排序，初始化i＝1，并执行步骤403-2-1：

403-2-1：判断内积为第i大的两行的原子数目之和是否大于s，若是，则执行步骤403-2-2，否则执行步骤403-3；

403-2-2：若i<L则令i＝i+1后，继续执行步骤403-2-1；若i＝L则分组聚类算法结束，执行步骤404；

403-3：将内积为第i大的两行合并，得到更新后的矩阵X′，其中合并的方式为：两个行向量相加除以2；

并基于内积为第i大的两行的行索引值对分组向量d进行更新后，继续执行步骤403-2；

其中对分组向量d的更新具体为：用p、q表示内积为第i大的行索引值，将分组向量d中的行索引值p、q全部修改为p与q中的较小索引值；

404：将分组信息d中同一取值的列作为一组，得到组结构字典D的分组方式，并基于组结构字典D的分组结果更新稀疏表示系数矩阵X，即X＝argmin_X‖Y-DX‖_F且s.t.||X_j||_0,d≤k,j＝1,…L，其中||X_j||_0,d表示稀疏表示系数矩阵的第j列不为0的元素对应的分组数目；

405：逐组更新组结构字典D；

406：判断是否满足迭代收敛条件，若是，则基于当前组结构字典D执行步骤5；否则继续执行步骤401；

其中迭代收敛条件为：达到最大迭代次数或最近两次迭代得到的组结构字典D的重构信号的误差小于或等于预设阈值。

3.根据2所述的方法，其中步骤6采用基于组结构的L1范数最大的方法进行融合的具体为：

根据权利要求1的步骤4学习得到组结构字典D，基于组结构的L1范数最大的融合规则为：

若则否则其中j＝1…M,k＝1…B，B表示分组向量d中元素的最大值，分别表示稀疏系数向量的第j列中属于第k组元素的绝对值之和，X_i,d(j)＝k表示融合系数向量X_i的第j列中属于第k组的元素。