CN103632345B - 一种基于正则化的mri图像非均匀性校正方法 - Google Patents

Abstract

本发明提供一种基于正则化的MRI图像非均匀性校正方法，通过对低通滤波方法进行改进，选取灰度值相差不大的一部分区域进行低通滤波、以获得非均匀场的初始估计，然后通过一种正则化的多项式拟合方法进行外推，获得全图像的非均匀场估计。本发明确定感兴趣区域进行滤波，以避免全局图像滤波方法中的边界效应；正则化方法拟合，能够保证在提高拟合阶次的情况下不会过拟合。

Description

一种基于正则化的MRI图像非均匀性校正方法

技术领域

本发明涉及核磁共振、图像处理的技术领域，具体涉及一种基于正则化的MRI图像非均匀性校正方法。

背景技术

近十几年来核磁共振成像由于其高分辨率、无电离辐射损伤以及任意角度成像的优点被广泛地应用于医疗诊断。随着对更高分辨率扫描图像的需求，扫描仪的磁场强度越来越高、磁场梯度也越来越精细。但随之而来的问题是核磁共振图像所遭受的非均匀场的干扰也越来越严重。非均匀场是由发射的空间磁场不均匀或接收线圈的不均匀灵敏度所引起的偏差场，这种偏差场一般被假设为一种平滑的、缓慢变化的乘性偏差场，会导致图像的灰度值与真实值之间存在一定的偏差。非均匀场会对临床的定性诊断造成一定的影响，并且对图像配准、图像分割等定量的计算机辅助诊断带来很大的干扰。图像配准与图像分割等技术假设同一人体组织内所有体素的灰度值应该是相同的，而非均匀场的存在破坏了这一假设。因此，非均匀场的校正对每一幅核磁共振图像都是必不可少的(参见文献UrosVovk,FranjoPernus,andBostjanLikar.AReviewofMethodsforCorrectionofIntensityInhomogeneityinMRI,IEEETrans.onMedicalimaging,Vol.26,No.3,405-420,2007)。

相关现有技术如下：

方案一

方案名称：基于特殊非均匀场采集装置的非均匀性估计(参见文献RoemerPB,EdelsteinWA,HayesCE,SouzaSP,MuellerOM.TheNMRphasedarray.MagnResonMed1990；16:1992-225.)

算法思想：先进行一次预扫描，利用体线圈获得一幅较低分辨率均匀的参考图像。然后用正式扫描图像与参考图像进行图像配准，再用配准后的图像除以参考图像，得到非均匀性信息，最后进行平滑处理。

算法缺点：二次扫描需要额外的硬件设备，增加了扫描时间；并且两次扫描时病人的位移等造成图像不匹配，需要图像配准，使得估计非均匀场的存在误差，造成伪影。并且，低信噪比区域的估计值不可靠，需要外推。

方案二

方案名称：基于低通滤波的线圈灵敏度估计(参见文献BrinkmannBH,ManducaA,RobbRA.OptimizedhomomorphicunsharpmaskingforMRgrayscaleinhomogeneitycorrection.IEEETransMedImaging1998；17:161-171.)

算法思想：将非均匀场看作是一个低频的、缓慢变化的偏差场，对非均匀图像进行低通滤波，提取出图像的低频信息，然后通过样条插值或多项式拟合等做平滑处理。

算法缺点：非均匀场与图像的低频信息存在混叠，低通滤波的频率阈值难以确定。在图像的高对比度区域存在边界效应，边界处灰度值变化被过度抬高，造成图像信息丢失。

方案三

方案名称：基于图像分类的非均匀场校正(R.GuillemaudandM.Brady,“EstimatingthebiasfieldofMRimages,”IEEETrans.MedImag.,vol.16,no.3,pp.238–251,Jun.1997.)

算法思想：利用图像的灰度变化的统计特性对图像进行分类，在分类的同时获取像素点的非均匀性信息。

算法缺点：普适性较差，只适用于特殊部位的图像。分类时需要事先指定多项参数。

发明内容

本发明目的为：1)能够较为精确地估计出非均匀场；2)估计方法不依赖于其他设备和先验知识；3)估计方法运算速度满足实时要求。

本发明采用的技术方案为：一种基于正则化的MRI图像非均匀性校正方法，利用一种正则化的拟合方法获取平滑的非均匀场估计，具体流程如下：

步骤一、感兴趣区域的确定：

首先获得去掉背景区域的图像直方图，找到峰值p，然后计算出图像噪声的方差σ²；最后，感兴趣区域RoI，即标识感兴趣区域的二维矩阵，每个元素对应的坐标如果是感兴趣区域，标记为1，否则标记为0，通过以下公式得到：

$RoI\left(x\right)=\left\{\begin{array}{cc}1& p-5\mathrm{\&sigma;}<v\left(x\right)<p+5\mathrm{\&sigma;}\\ 0& other\end{array}\right.$

其中RoI(x)是点x是否属于感兴趣区域的指示器，如果x属于感兴趣区域，则RoI(x)＝1，如果不属于则RoI(x)＝0，区间[p-5σ,p+5σ]即代表感兴趣区域的图像灰度变化范围；

步骤二、低通滤波：

确定出感兴趣区域之后，对区域内的点进行同态域低通滤波，假设待校正的原始图像用v表示，则点x的灰度值为v(x)，滤波步骤如下：

1.将RoI内的图像变换到对数域：

I_log＝log(v(RoI))

2.对I_log进行低通滤波：

I_filt＝LPF(I_log)

3.由于低通滤波在边缘处包含大量空白区域，为了防止边界效应，对感兴趣区域的二元掩膜RoI进行滤波：

RoI_filt＝LPF(RoI)

4.最后通过以下公式得到滤波后的粗略的非均匀度：

$\hat{u}\left(RoI\right)=\exp (I_{filt}/{\mathrm{RoI}}_{filt})$

由于滤波校正的缺点，肯定多余地包含了许多图像低频信息，并存在边界效应，并且只对感兴趣区域内的点进行了校正；

步骤三、正则化多项式拟合：

一旦估计出感兴趣区域内的局部非均匀场，全局的非均匀场就能够通过多项式拟合，为了防止过拟合现象，利用了正则化的拟合方法来达到较为精确的估计，选择n阶多项式作为拟合曲面，对估计的感兴趣区域内的非均匀度进行拟合，假设选择k阶多项式，则包含项其中l+m≤k,l≥0，m≥0并且r_x与r_y分别为图像像素的水平方向和竖直方向的坐标，对于给定的k阶多项式，一共有K＝(k+1)(k+2)/2个项，设为F_i(r)(i＝1,…,K)，假设各项系数为w_i,i＝1,...,K,则多项式拟合通过最小化下列代价函数来实现：

$J\left(W\right)=\underset{W}{min}\left\{\frac{1}{2}\underset{n=1}{\overset{N}{\&Sigma;}}{(\hat{U}-FW)}^2+\frac{\mathrm{\&lambda;}}{2}\left|\right|W|{|}^2\right\}$

其中 $W=\left[\begin{array}{c}{\mathrm{\&omega;}}_1\\ {\mathrm{\&omega;}}_2\\ .\\ .\\ .\\ {\mathrm{\&omega;}}_K\end{array}\right],$ F＝[F₁(r)F₂(r)…F_K(r)]，λ为惩罚因子，是感兴趣区域中的每个点的非均匀度组成的向量，通过加入使得拟合的多项式参数不至于过大，从而抑制了多项式的过拟合，使得最终曲面为平滑缓慢变化的，更符合精确的非均匀场，

上述公式的求解如下：

在此步骤中，λ是一项决定非均匀场平滑度的重要参数，通过比较不同λ值拟合出的非均匀场来调节获得一个最佳的λ；

获得多项式系数之后，就能够通过外推得到全图像的非均匀估计：

g(r)＝[F₁(r),F₂(r),…,F_K(r)]W

步骤四、迭代：

当得到最终估计的非均匀场后，对原始图像进行校正，假设原始图像的灰度值为v(x)，那么有：

I_correct(x)＝v(x)/g(x)

在校正时不包含空白的背景区域，即x只代表图像区域的坐标，获得校正的图像，将I_correct(x)赋值给v(x)，并带入下一次迭代重复以上校正过程，对图像进行下一轮校正，这样能有效提高非均匀场估计的准确度，当估计出的非均匀场的变化范围在感兴趣区域小于一定的阈值时，代表图像已经足够均匀，这时停止迭代，最终获得了较为均匀的校正图像。

本发明技术方案的优点和积极效果为：

(1)、本发明选取一个灰度变化范围较小的区域进行低通滤波以减小边界效应。

(2)、本发明通过正则化的拟合方法有效避免了过拟合，可以适当提高拟合阶次，达到更精确的估计，从而能得到平滑的、又足够精确的非均匀场估计。

(3)、本发明方法易于实现，复杂度低，实际运用时耗时极低。

附图说明

图1为本发明一种基于正则化的MRI图像非均匀性校正方法的流程图；

图2为本发明用于校正头部图像的效果图。

具体实施方式

下面结合附图以及具体实施方式进一步说明本发明。

本发明通过对低通滤波方法进行改进，选取灰度值相差不大的一部分区域进行低通滤波、以获得非均匀场的初始估计，然后通过一种正则化的多项式拟合方法进行外推，获得全图像的非均匀场估计。具体流程如下：

一、感兴趣区域的确定

直接对图像进行滤波的方法容易造成边界效应，即在高对比度区域的非均匀度失真，在随后的插值时就容易造成较大的误差，而实际上在插值拟合时大部分的数据又都是冗余的、不必要的。因此在本发明中我们通过选取部分灰度值接近的区域、而不是整个图像区域进行低通滤波来避免边界效应。感兴趣区域的确定是通过在图像的直方图中确定一个阈值来实现的。为了确定感兴趣区域，我们首先获得去掉背景区域的图像直方图，找到峰值p,然后计算出图像噪声的方差σ²(通过计算背景区域的像素点灰度值平方和的一半获得)；最后，感兴趣区域RoI(即标识感兴趣区域的二维矩阵，每个元素对应的坐标如果是感兴趣区域，标记为1，否则标记为0)通过以下公式得到：

$RoI\left(x\right)=\left\{\begin{array}{cc}1& p-5\mathrm{\&sigma;}<v\left(x\right)<p+5\mathrm{\&sigma;}\\ 0& other\end{array}\right.$

其中RoI(x)是点x是否属于感兴趣区域的指示器，如果x属于感兴趣区域，则RoI(x)＝1，如果不属于则RoI(x)＝0.区间[p-5σ,p+5σ]即代表感兴趣区域的图像灰度变化范围。

二、低通滤波

确定出感兴趣区域之后，对区域内的点进行同态域低通滤波。

假设待校正的原始图像用v表示，则点x的灰度值为v(x)。滤波步骤如下：

1.将RoI内的图像变换到对数域：

I_log＝log(v(RoI))

2.对I_log进行低通滤波：

I_filt＝LPF(I_log)

3.由于低通滤波在边缘处包含大量空白区域，为了防止边界效应，对感兴趣区域的二元掩膜RoI进行滤波：

RoI_filt＝LPF(RoI)

4.最后通过以下公式得到滤波后的粗略的非均匀度：

$\hat{u}\left(RoI\right)=\exp (I_{filt}/{\mathrm{RoI}}_{filt})$

由于滤波校正的缺点，肯定多余地包含了许多图像低频信息，并存在边界效应，并且只对感兴趣区域内的点进行了校正。下面就通过拟合对非均匀场进行进一步的平滑处理和外推。

三、正则化多项式拟合

一旦估计出感兴趣区域内的局部非均匀场，全局的非均匀场就可以通过多项式拟合。为了防止过拟合现象，我们利用了正则化的拟合方法来达到较为精确的估计。选择n阶多项式作为拟合曲面，对估计的感兴趣区域内的非均匀度进行拟合。假设我们选择k阶多项式，则包含项其中l+m≤k,l≥0，m≥0并且r_x与r_y分别为图像像素的水平方向和竖直方向的坐标。对于给定的k阶多项式，一共有K＝(k+1)(k+2)/2个项，设为F_i(r)(i＝1,…,K)，例如二阶二维多项式包含[1xyxyx²y²]项，假设各项系数为w_i,i＝1,...,K,则多项式拟合通过最小化下列代价函数来实现：

$J\left(W\right)=\underset{W}{min}\left\{\frac{1}{2}\underset{n=1}{\overset{N}{\&Sigma;}}{(\hat{U}-FW)}^2+\frac{\mathrm{\&lambda;}}{2}\left|\right|W|{|}^2\right\}$

其中 $W=\left[\begin{array}{c}{\mathrm{\&omega;}}_1\\ {\mathrm{\&omega;}}_2\\ .\\ .\\ .\\ {\mathrm{\&omega;}}_K\end{array}\right],$ F＝[F₁(r)F₂(r)…F_K(r)]，λ为惩罚因子，是感兴趣区域中的每个点的非均匀度组成的向量。通过加入使得拟合的多项式参数不至于过大，从而抑制了多项式的过拟合，使得最终曲面为平滑缓慢变化的，更符合精确的非均匀场。

上述公式的求解如下：

在此步骤中，λ是一项决定非均匀场平滑度的重要参数，在本发明中，我们通过比较不同λ值拟合出的非均匀场来调节获得一个最佳的λ。

获得多项式系数之后，就可以通过外推得到全图像的非均匀估计：

g(r)＝[F₁(r),F₂(r),…,F_K(r)]W

四、迭代

当得到最终估计的非均匀场后，我们对原始图像进行校正，假设原始图像的灰度值为v(x)，那么有：

I_correct(x)＝v(x)/g(x)

(在校正时不包含空白的背景区域，即x只代表图像区域的坐标)。获得校正的图像，我们将I_correct(x)赋值给v(x)，并带入下一次迭代重复以上校正过程，对图像进行下一轮校正，这样能有效提高非均匀场估计的准确度。当估计出的非均匀场的变化范围在感兴趣区域小于一定的阈值时，代表图像已经足够均匀，这时停止迭代，最终我们获得了较为均匀的校正图像。

图2是本发明用于处理头部图像的效果及直方图对比。由图的均匀度可以看出，头部图像中心的低信号区域有所抬高，图像细节更加明显，由图像直方图可以看出，校正后的头部图像灰质与白质波峰更加明显，更容易区分。

本发明未详细公开的部分属于本领域的公知技术。

尽管上面对本发明说明性的具体实施方式进行了描述，以便于本技术领域的技术人员理解本发明，但应该清楚，本发明不限于具体实施方式的范围，对本技术领域的普通技术人员来讲，只要各种变化在所附的权利要求限定和确定的本发明的精神和范围内，这些变化是显而易见的，一切利用本发明构思的发明创造均在保护之列。

Claims (1)

1.一种基于正则化的MRI图像非均匀性校正方法，其特征在于利用一种正则化的拟合方法获取平滑的非均匀场估计，具体流程如下：

步骤一、感兴趣区域的确定：

首先获得去掉背景区域的图像直方图，找到峰值p，然后计算出图像噪声的方差σ²；最后，感兴趣区域RoI，即标识感兴趣区域的二维矩阵，每个元素对应的坐标如果是感兴趣区域，标记为1，否则标记为0，通过以下公式得到：

$RoI\left(x\right)=\left\{\begin{array}{cc}1& p-5\mathrm{\&sigma;}<v\left(x\right)<p+5\mathrm{\&sigma;}\\ 0& other\end{array}\right.$

其中RoI(x)是点x是否属于感兴趣区域的指示器，如果x属于感兴趣区域，则RoI(x)＝1，如果不属于则RoI(x)＝0，区间[p-5σ,p+5σ]即代表感兴趣区域的图像灰度变化范围；

步骤二、低通滤波：

确定出感兴趣区域之后，对区域内的点进行同态域低通滤波，假设待校正的原始图像用v表示，则点x的灰度值为v(x)，滤波步骤如下：

1.将RoI内的图像变换到对数域：

I_log＝log(v(RoI))

其中，v代表扫描获得的图像，RoI代表感兴趣区域；

2.对I_log进行低通滤波：

I_filt＝LPF(I_log)

其中，LPF(·)代表低通滤波器；

3.由于低通滤波在边缘处包含大量空白区域，为了防止边界效应，对感兴趣区域的二元掩膜RoI进行滤波：

RoI_filt＝LPF(RoI)

4.最后通过以下公式得到滤波后的粗略的非均匀度：

$\hat{u}\left(RoI\right)=\exp (I_{filt}/{\mathrm{RoI}}_{filt})$

其中，即代表粗略估计的非均匀度；

由于滤波校正的缺点，肯定多余地包含了许多图像低频信息，并存在边界效应，并且只对感兴趣区域内的点进行了校正；

步骤三、正则化多项式拟合：

一旦估计出感兴趣区域内的局部非均匀场，全局的非均匀场就能够通过多项式拟合，为了防止过拟合现象，利用了正则化的拟合方法来达到较为精确的估计，选择n阶多项式作为拟合曲面，对估计的感兴趣区域内的非均匀度进行拟合，假设选择k阶多项式，则包含项其中l+m≤k,l≥0，m≥0并且r_x与r_y分别为图像像素的水平方向和竖直方向的坐标，对于给定的k阶多项式，一共有K＝(k+1)(k+2)/2个项，设为F_i(r)(i＝1,…,K)，假设各项系数为ω_i,i＝1,…,K,则多项式拟合通过最小化下列代价函数来实现：

$J\left(W\right)=\underset{W}{min}\{\frac{1}{2}\underset{n=1}{\overset{N}{\&Sigma;}}{(\hat{U}-FW)}^2+\frac{\mathrm{\&lambda;}}{2}|\left|W\right|{|}^2\}$

其中 $W=\left[\begin{array}{c}{\mathrm{\&omega;}}_1\\ {\mathrm{\&omega;}}_2\\ \&CenterDot;\\ \&CenterDot;\\ \&CenterDot;\\ {\mathrm{\&omega;}}_K\end{array}\right],$ ||W||代表向量W的欧几里得范数，F＝[F₁(r)F₂(r)…F_K(r)]，λ为惩罚因子，是感兴趣区域中的每个点的非均匀度组成的向量，通过加入使得拟合的多项式参数不至于过大，从而抑制了多项式的过拟合，使得最终曲面为平滑缓慢变化的，更符合精确的非均匀场，

上述公式的求解如下：

其中F^t代表矩阵F的转置，在此步骤中，λ是一项决定非均匀场平滑度的重要参数，通过比较不同λ值拟合出的非均匀场来调节获得一个最佳的λ；

获得多项式系数之后，就能够通过外推得到全图像的非均匀估计：

g(r)＝[F₁(r),F₂(r),…,F_K(r)]W

步骤四、迭代：

当得到最终估计的非均匀场后，对原始图像进行校正，假设原始图像的灰度值为v(x)，那么有：

I_correct(x)＝v(x)/g(x)

在校正时不包含空白的背景区域，即x只代表图像区域的坐标，获得校正的图像，将I_correct(x)赋值给v(x)，并带入下一次迭代重复以上校正过程，对图像进行下一轮校正，这样能有效提高非均匀场估计的准确度，当估计出的非均匀场的变化范围在感兴趣区域小于一定的阈值时，代表图像已经足够均匀，这时停止迭代，最终获得了较为均匀的校正图像。