一种基于瑞利噪声分布的超声与磁共振影像融合方法与流程

本发明涉及医学影像图像处理技术领域，具体是一种基于瑞利噪声分布及非负矩阵分解的超声与磁共振影像融合方法。

背景技术：

超声成像由于其实时性好、成本低等优点，在临床上有着广泛的应用。和磁共振和CT成像磁共振比较，超声成像分辨率不高。磁共振是临床上另一种被广泛使用的成像模式，对软组织成像效果好，并有着较高的空间分辨率。除此之外，磁共振还能够提供临床上需要的功能信息。但是磁共振成像实时性不好。目前，没有任何一种成像模式能够满足临床上得所有成像需求。各种成像模式在临床上有着各自的有点和缺点，在很多情况下，超声和磁共振在临床上是互相补充的，为了能够更好的辅助临床医生进行诊断，有必要对超声和磁共振进行融合。

由于超声影像斑点噪声突出、对比度低，直接对超声影像与其他模态影像(包括磁共振)进行融合将导致融合后影像噪声突出、清晰度降低。目前融合算法的研究集中于以下4个方面:基于数学/统计学的融合算法，包括加权平均法、NMF方法、PCA方法、贝叶斯估计融合和EM融合方法；基于多分辨率分析的图像融合算法，包括多分辨金字塔法、小波变换法；基于彩色空间的图像融合方法，包括HSI变换法、视觉神经生理学法；智能图像融合算法，主要包括模糊理论方法、神经网络法等。上述融合算法本质都是构建和优化数学模型，对两幅或多幅参考图像的信息进行综合，使融合后的图像具有更高的清晰度和信息量。上述融合算法通过降低图像清晰区域的图像质量为代价提升图像模糊区域的图像质量，以得到目标均清晰的融合结果，导致融合图像中一些重要边缘信息的丢失。

技术实现要素：

为了克服上述现有技术的缺点，本发明提供一种基于瑞利噪声分布及非负矩阵分解的超声与磁共振影像融合方法，提高融合图像的对比度。

本发明是以如下技术方案实现的：一种基于瑞利噪声分布的超声与磁共振影像融合方法，其具体步骤如下：

步骤1：基于自适应加权非负矩阵分解A＝WUV，首先根据行优先，将超声图像和磁共振图像分别转化为列向量US和MRI，构建非负矩阵A，并令A＝[US,MRI]，构建非负的权重矩阵W＝[WUS,WMRI]，并且非负的权重矩阵W的每一个元素值都被初始化为0.5；W和A均为k×2矩阵，这里k＝n×m；U、V分别为非负的k×r矩阵和非负的r×2矩阵，其中r＝1，即U为按照行优先转换后的列向量；

步骤2：计算超声图像和磁共振图像对应的梯度能量矩阵并记为EOGUS和EOGMRI，并对上述两个梯度能量矩阵对应的元素进行比较，对应元素的较大值构成梯度能量矩阵EOGUS～MRI；

步骤3：对U和V进行随机初始化；

步骤4：对U和V进行迭代求解；

步骤5：更新非负的权重矩阵W；

步骤6：如果迭代满足结束条件||A-WUV||2＜0.01，那么转到步骤7；如果迭代不满足结束条件，转到步骤4；

步骤7：基于行优先，将向量U转换为矩阵，即为融合结果，融合过程结束；

步骤5其更新非负的权重矩阵W具体方法如下：

步骤5.1：迭代一次后，将向量U按照行优先转化为矩阵U，并计算U对应的梯度能量矩阵，记为EOGU；

步骤5.2：将WUS和WMRI按照行优先转化为矩阵；

步骤5.3对于超声和磁共振影像的每一个元素，如果EOGU(i,j)≧EOGUS～MRI(i,j)，进入步骤5.4，如果EOGU(i,j)<EOGUS～MRI(i,j)，则对应的W自适应更新算法如下：

SUM(i,j)＝WUS(i,j)+WMRI(i,j)；

步骤5.4：将WUS和WMRI按照行优先转换分别为列向量，并令W＝[WUS,WMRI]。

本发明的有益效果：将自适应非负矩阵算法用于本发明中，获得的融合图像具有更强的对比度。

附图说明

图1是本发明融合流程图；

图2是基于不同融合方法的融合结果比较图。

具体实施方式

一种基于瑞利噪声分布的超声与磁共振影像融合方法，具体步骤如下：

步骤1：基于自适应加权非负矩阵分解A＝WUV，首先根据行优先，将超声图像和磁共振图像分别转化为列向量US和MRI，构建非负矩阵A，并令A＝[US,MRI]，构建非负的权重矩阵W＝[WUS,WMRI]，并且非负的权重矩阵W的每一个元素值都被初始化为0.5；W和A均为k×2矩阵，这里k＝n×m；U、V分别为非负的k×r矩阵和非负的r×2矩阵，其中r＝1，即U为按照行优先转换后的列向量；

步骤2：计算超声图像和磁共振图像对应的梯度能量矩阵并记为EOGUS和EOGMRI，并对上述两个梯度能量矩阵对应的元素进行比较，对应元素的较大值构成梯度能量矩阵EOGUS～MRI；

步骤3：对U和V进行随机初始化；

步骤4：对U和V进行迭代求解；

步骤5：更新非负的权重矩阵W；

步骤5.1：迭代一次后，将向量U按照行优先转化为矩阵U，并计算U对应的梯度能量矩阵，记为EOGU；

步骤5.2：将WUS和WMRI按照行优先转化为矩阵；

步骤5.3对于超声和磁共振影像的每一个元素，如果EOGU(i,j)≧EOGUS～MRI(i,j)，进入步骤5.4，如果EOGU(i,j)<EOGUS～MRI(i,j)，则对应的W自适应更新算法如下：

SUM(i,j)＝WUS(i,j)+WMRI(i,j)；

步骤5.4：将WUS和WMRI按照行优先转换分别为列向量，并令W＝[WUS,WMRI]。

步骤6：如果迭代满足结束条件||A-WUV||2＜0.01，那么转到步骤7；如果迭代不满足结束条件，转到步骤4；

步骤7：基于行优先，将向量U转换为矩阵，即为融合结果，融合过程结束；

本发明原理如下：

设A为非负的待分解的n×m矩阵，U、V分别为非负的n×r矩阵和非负的r×m矩阵，非负矩阵分解的计算模型为：A＝UV。

非负矩阵分解即为求解一个非负的n×r矩阵U和非负的r×m矩阵V，在此基础上构建代价函数L(A|U,V)。L(A|U,V)作为能够度量A和UV逼近程度的代价函数，当非负矩阵分解算法迭代到一定次数，代价函数L(A|U,V)最小，迭代结束。本发明拟基于超声影像的噪声分布，构建一种新型的代价函数，基于上述代价函数，提出了自适应加权非负矩阵分解算法并进行超声与磁共振的融合。下面给出新型代价函数及自适应加权非负矩阵分解算法(A＝WUV)的构建过程，其中W为权重矩阵。

超声影像的噪声分布符合瑞利分布，瑞利分布定义如下：

在(1)式中x≥0，a为常数且a＞0.为了简化后续计算，a被设置为常数1.对(1)式进行求偏导：

求得

当x等于时，px(x)获得最大值。由于超声的斑点噪声符合噪声分布，本发明提出了一种新的架构如下：

在式(2)中，A为一个n×m非负矩阵。A被近似的分解为两个非负的矩阵U(n×r)和V(r×m)。

设

那么(3)式对应的最大似然解，即为对如下损失函数进行求最小极值：

为了对(4)式进行化简，将预先加到矩阵A的所有元素上，(4)式可以化简为：

通过最小化(6)式构建的加权代价函数，我们能够推导出自适应加权非负矩阵分解算法。

在(6)式中W＝{Wij}＞0为非负的权重矩阵，XοY代表矩阵X和矩阵Y的Hadamard积(即矩阵X和矩阵Y点乘)，这里我们定义(6)式为加权瑞利测度。

在本发明中我们仅给出V的迭代公式的推导过程，因为我们可以用同样的方法获得U的迭代式。A,V和W中任意一列向量分别表示为a,v和w，我们引入部分代价函数F(v)：

在(7)式中：

DW＝diag(w), (8)

其中DW为基于w构建的对角方阵。对(9)式和(10)式分别进行一阶和二阶求偏导，分别得到式(11)、(12)、(13)和(14)：

式(13)和(14)中的Ui代表U的第i列。对(7)式进行迭代，可以得到v的迭代式(15)

同理能够得到u的迭代式：

(16)式中，A^T,U^T和W^T分别为A,U和W对应的转置矩阵,a^T,u^T和w^T分别为A^T,U^T和W^T对应的列向量，其中的一阶偏导和二阶偏导定义见式(17)、(18)、(19)和(20)：

在式(17)、(18)、(19)和(20)中，(V^T)i代表V^T的第i列,DwT代表基于w^T构建的对角方阵。

在图像融合过程中，随着迭代的进行，融合图像像素灰度也是动态变化的，图像各区域的灰度分布也会发生变化，因此图像各区域的权重也应该发生相应的变化。采用固定权重矩阵不能体现这种动态变化过程。针对以上问题，本文提出一种动态迭代策略，随着迭代的进行，根据最新迭代结果对加权系数进行动态调整。对于迭代结果的特定指标改善已经达到目标阈值的情况，则停止对对应加权系数的调整；反之，加权系数朝有利于特定指标改善的方向进行调整。WNMF的加权系数体现了对应元素的重要程度，通过对加权系数进行持续动态的调整，即可实现融合图像特定指标的改善。

本发明采用平均梯度(average gradient，AG)、边缘强度(edge intensity，EI)、信息熵(entropy，EN)和空间频率(space frequency，SF)作为融合图像的评价标准。AG也称为清晰度，反映了图像的细节反差与纹理变化特征，融合图像的AG值越大，表明融合图像的清晰度越高；EN表示图像所包含信息量的多少，EN值越大，表明融合图像包含的信息量越大；EI反应了边缘与邻近像素的对比强度，EI值小意味着强度变化不大而难以分辩，值越大表示边缘越清晰；SF反应了图像空间域的活跃程度，目标与背景对比明显，图像中存在剧烈的突变等情况下，SF值较大，反之较小。

为了对比本发明提出的融合方法的有效性，另外分别基于DWT(discrete wavelet transform)的融合算法、基于PCA的融合算法、基于SIDWT(shift invariance discrete wavelet transform)的融合算法对两只裸鼠的超声和磁共振影像进行了融合实验，实验结果如图2所示，直观观察，本发明提出方法获得的融合图像具有更强的对比度。为了定量比较，基于四种方法分别进行了融合指标计算，计算结果见表1和表2。分析表1和表2，使用本发明提出的方法对应的AG、EI、EN和SF值更大，说明本发明或的的融合结果在定量指标方面表现更好。

表1.基于第一只裸鼠的融合结果指标比较.

表2.基于第二只裸鼠的融合结果指标比较.

图2

图2(a)～(d)分别为基于第一只裸鼠使用本发明方法、DWT、PCA以及SIDWT进行融合的结果；

图2(e)～(h)分别为基于第二只裸鼠使用本发明方法、DWT、PCA以及SIDWT进行融合的结果。