1　引言

       前列腺癌已成为严重危害老年男性健康的疾病，而且据相关报道，中国前列腺癌发病率正呈上升趋势^[1,2]。前列腺癌的早期诊断和检测对疾病的治疗预防非常重要。考虑到传统前列腺癌早期检测方法的局限性，一些借助于MRI技术的人工检测方法，以及结合前列腺癌的病理学知识与MR影像特征的计算机辅助诊断方案也在研究之中^[3,4,5,6,7]。

       这些借助于MR影像的诊断方法对于影像的质量有较高要求，然而医学图像中普遍存在着一种灰度不均匀现象，也即偏场。偏移场的存在会影响病灶区的视觉效果以及计算机对病灶组织的判断。根据文献[8,9]，偏移场的产生是由于成像设备和患者特异性的综合作用。目前对偏移场的校正方法大致可以分为前瞻性方法和回顾性方法^[9]。前瞻性方法是从MR图像获得过程的角度来考虑校正方法，这类方法虽然也能取得一些效果，但是往往不能解决由于患者特异性所带来的影响。回顾性方法只依赖于获得图像的信息，因而能更全面地解决偏场，当然也得到了更多的关注和研究。Johnston等^[10]最先提出一种基于同态滤波的方法，这个方法假设偏移场具有低频特性，但是图像本身也具有很多低频成分，因此高通滤波器也不能取得理想的结果。Miles等^[11]提出一个结合直方图信息和B样条曲面拟合的方法。Li等^[12]在乘法固有成分优化的思想下，提出一个能量最小化模型，同时实现偏场矫正和多组织分割。

       目前对于偏场矫正的方法很多，本文在已有算法的基础上作了一些改进，提出一种针对前列腺组织MR图像的偏场校正方法。首先根据前列腺磁共振图像的灰度特征进行自适应的预处理，然后利用峰值检测技术从图像的直方图信息中自动获得前列腺组织的类别参数，接着根据图像模型构造能量最小化模型，最后提出一种由最简多项式函数与三角函数组合而成的曲面拟合方法，实现对偏场的平滑拟合。

2　相关阐述

2.1　MR图像模型

       目前一种普遍认可的磁共振图像模型如式(1)所示：

       其中x代表图像像素点，v(x)表示获得的图像，u(x)代表没有偏场与噪声作用的真实图像，b(x)是偏场，n(x)则指图像噪声。通常认为图像中的背景与前景中的各个目标在灰度上具有一致性，而图像中的灰度不均匀变化是一个平滑过程。所以在这个图像模型中，笔者认为代表真实图像的u(x)具有分段常量特性，即图像中不同的组织可以近似为不同的灰度常量区域，而表征图像不均匀现象的b(x)具有平滑特性。

点击查看大图

2.2　能量函数

       令Ω代表图像区域，根据2.1中的图像模型得到能量函数：

       针对u(x)的分段常量特性，假设图像中存在n个不同的组织，每个组织的灰度值为常量ci(i=1,...,n)。隶属度函数μi(x)表示像素点x属于第i类组织的概率。为了更真实地反映实像，这里参照模糊C均值(fuzzy c-means algorithm, FCM)算法采用一种模糊表示方法()。那么实像u(x)表示如下：

       其中q是加权指数。对于偏场，笔者用一系列连续函数的加权和表示：

       m的取值决定着曲面的拟合精度。通常对于偏场的拟合会用多项式函数作为连续函数，但是笔者发现三角函数有着等同高次函数的平滑特性，而且计算更简。所以用低次最简多项式函数与三角函数共同作为拟合曲面的基函数。

       令a=(a1,...,am)^T，G(x)=(g1(x),...,gm(x))^T，得到偏场的向量表达形式：

       将(3)、(5)代入(2)中得到能量的表达：

       根据文献[11]中的能量优化过程，得到能量的变化形式：

点击查看大图

点击查看大图

点击查看大图

点击查看大图

点击查看大图

点击查看大图

2.3　能量优化

       2.2部分提出的能量函数实际是噪声的表达。当噪声取得最小值时，可以得到最接近真实前列腺磁共振图像的校正图像。由于能量函数是凸函数，所以通过一种类似梯度下降法的迭代方法便可使能量达到全局最小。

       对于变量c、u、a分别固定其中两个，然后让能量函数对另外一个求导，重复多次，直到收敛。计算得到变量的迭代式如下：

点击查看大图

3　实施

       前列腺磁共振图像偏移场校正的流程图见图1，详细步骤将在后续部分介绍。

3.1　预处理

       采用Pual等^[13]开发的ITK-SNAP软件(www.itksnap.org)对前列腺磁共振图像进行分割获得感兴趣组织。然后通过一个直方图线性拉伸，使感兴趣组织获得更好的视觉效果，预处理结果见图2。

点击查看大图

图2  预处理过程。A：原始图像；B：分割结果；C：直方图拉伸结果
Fig. 2  Preprocessing. A: Original image; B: Segmentation results; C: Results of histogram stretch.

3.2　自动获取组织的类别参数

       本研究笔者直接对感兴趣组织的直方图进行处理。先对直方图进行均值滤波，消除波形中的不规则点(灰度值出现跳变的点)，平滑波形。然后对波形进行峰值检测^[14]，根据检测策略判断图像直方图峰值的个数，由于图像中不同的组织会呈现不同的灰度直方图分布，从而自动获得组织的类别参数，波形平滑过程如图3所示。

点击查看大图

图3  波形平滑过程。A为原图；B为滤波后的图形
Fig. 3  Procedure of waveform smoothing. A: Original outline of ROI; B: Outline after fitting.

3.3　能量优化

       能量优化过程是整个校正过程的核心。首先给定能量函数的初始化参数c、μ；然后分别对c、μ、a进行迭代运算，使能量函数沿着梯度下降的方向变化；接着检查能量函数的值，通过比较判断是否达到收敛；根据判断结果再次执行迭代运算(不收敛)或者评估偏场(收敛)、校正输出。

3.4　加权指数的选择

       据文献[14]中描述，加权指数是模糊聚类算法中一个重要的参数，它的选取决定着偏场矫正过程中能量优化的速度和程度。笔者通过图4、图5所示的两组对比实验来确定指数q的值。

       通过以上两组对比实验可以得出以下结论：(1)当q值取整数时能量函数的迭代速度更快；(2)随着q值的增大，能量函数能够达到更好的收敛程度。另外结合文献[15]中对于最佳q值的理论判断，笔者使用q=2作为前列腺磁共振图像偏场校正算法中的加权指数。

点击查看大图

图4  不同图像在不同q下的优化时间
图5  能量函数在不同q值下的收敛程度
Fig. 4  Show the time of four group different MR image in different exponent q.
Fig. 5  Show the different degree of convergence in different exponent q.

点击查看大图

图1  前列腺磁共振图像偏场校正过程
Fig. 1  Bias field correction of prostate MR image.

4　实验

       笔者将分别展示本研究的方法在定性和定量实验下的结果，所有的算法程序都运行在2.53GHzCPU，4G内存的Matlab2012a平台上。

4.1　定性评估

       笔者选取6幅来自不同个体的前列腺磁共振扫描数据中的横断面图像进行实验，定性的实验结果见图6。

点击查看大图

图6  定性实验结果。第一行是截取的横断面图像，第二行是预处理后的结果，第三行是偏场校正后的结果
Fig. 6  Qualitative experiment result. The three row are original prostate axial image, results of preprocessing, results of bias field correction.

4.2　定量评估

       为了定量评估算法的有效性，笔者将本研究的方法与Li等^[12]的方法在一些公认指标下进行对比。变异系数(coefficient of variation，CV)(式11)被定义为标准差与期望的比值，它反映了组织(T)的相对离散程度。通常认为真实图像的同种组织具有一致性，所以一个好的处理结果会得到一个小的CV值。本研究的方法与Li等^[12]的方法对6幅前列腺横断面图像处理的结果见表1。

       均方根(root mean square，rms)(式12)反映了处理结果与真实图像的差异。较小的rms值反映较好的处理结果。由于rms指标需要无偏场图像作为参考，所以笔者选取来自模拟脑数据库(http://brainweb.bic.mni.mcgill.ca/brainweb/)的10幅无噪声、无不均匀现象的脑部磁共振图像作为参考图像，将对应的10幅添加了3%噪声、20%偏移场的图像作为待处理图像进行实验。

       JS (jaccard similarity)(式13)反映了处理结果与标准图像的相似程度。S1、S2分别代表算法处理结果与参考图像，较大的JS值反映了更好的处理结果。

       巴氏距离(bhattacharyya distance)(式14)被认为是测量直方图相似性效果最好的度量工具。笔者通过测量校正图像p(x)与参考图像q(x)的直方图之间的巴氏距离来衡量算法的处理效果。据算法定义，较高的巴氏距离表示了更好的处理效果。

       本研究的算法与Li等^[12]的算法在各项指标下的实验情况见图7。实验结果表明，笔者的方法相比Li等^[12]的方法能取得更好的校正效果。

点击查看大图

点击查看大图

点击查看大图

点击查看大图

点击查看大图

图7  定量实验结果。左、中、右分别为在rms、JS、BD下的对比实验结果
Fig. 7  Quantitative experiment result. The left, middle, and right columns show the result on index rms, JS, BD respectively

点击查看表格

表1  本文方法与Li等^[12]的方法处理结果的CV值
Tab.1  Contrast of CV between ours and Li.et al’s method

5　总结

       本文提出了一个针对前列腺磁共振图像的偏移场校正方法。一系列实验和峰值检测算法被用来自动获取能量优化过程的初始参数，使得算法具有更好的准确性和稳定性。一组由三角函数和多项式函数组成的基函数保证了拟合偏场的平滑性。最后通过定性实验和与相关算法在定量指标下的对比，验证了本文提出的算法的有效性。