基于帧间相关性的乳腺MRI病灶的自动分割方法与流程

本发明涉及医学影像分割技术领域，特别涉及一种乳腺磁共振图像分割方法

背景技术：

乳腺癌是女性最常见的恶性肿瘤之一，严重威胁着女性身心健康及生命。乳腺癌的早期发现诊断和治疗有助于提高患者的年生存率和生存质量。磁共振成像技术(MRI)有公认的高敏感性(超过90％)，对各类软组织结构的分辨能量强，提供多个序列，多个方向上的图像，为明确病因提供了组织结构，病灶的形态，肿块分布等更加丰富的图像信息。乳腺MRI对乳腺癌早期诊断，治疗期观察，拟定合理的手术方案以及后期跟踪都有极大的辅助作用。

然而，由于在核磁共振成像的过程中存在电子噪声、偏移场失真与容积效应，导致图像中存在噪声与偏场效应，乳腺MRI图像中各种组织之间的边界往往模糊不清。此外，通常乳腺MRI的诊断大多依靠临床医生的个人经验与专业知识。计算机辅助诊断技术，是一种可以帮助放射科医生对医学影像进行筛查的新兴技术。它在处理大量病人影像数据时具有很高的实用性，通过大量的医学影像处理方法对某种疾病进行分析和建模，将病灶信息通过参数模型与正常组织区分开来，甄别出真正的病灶，给出标注作为医师进一步诊断的参考。计算机辅助诊断技术的出现和发展，大幅降低了医师的工作量，而更重要的是提高了乳腺癌诊断的准确率和灵敏度，从而挽救千万乳腺癌患者的生命，同时也对利用计算机准确识别和分割乳腺病灶提出了更高的要求。

病灶分割作为乳腺计算机辅助诊断系统中的一个研究热点，其主要分割目标有：1)为肿瘤性质做定量分析的病变体的分割以及后续的分类任务；2)相关组织结构(皮肤边缘、胸壁、心脏等)的自动轮廓获取。针对这两类分割，国内外很多学者做了相关研究并取得了一定成果，包括基于形态学的分割方法^[1-3]、基于聚类的分割方法^[4-6]、基于边缘的分割方法^[7，8]和基于特定理论的分割方法^[9，10]。但是各种分割方法均有其局限性，例如阈值法和区域增长法效率较低且对噪声较敏感，FCM和MRF等算法过于依赖数据，不能保证分割结果在解剖学上正确。最新的研究提出了两种算法相结合的分割方法，如结合SLIC超像素和水平集的分割算法^[11](SLIC+与DRLSE)，该方法用于乳腺X线图像的肿瘤分割，取得了较好的分割效果.然而由于MRI图像的区域灰度不均匀，该算法在MRI图像上的实验结果并不理想。首先，该算法采用的简单线性迭代聚类算法(simple linear iterative clustering，SLIC)产生的超像素形状差异较大，会对后续分割产生影响。细分割步骤采用了距离正则化水平集方法(DRLSE)，但是该方法要根据初始曲线的位置人为确定常量演化速度的符号，需要人工干预，不能实现自适应分割。另外，该方法只考虑了二维分割，没有充分利用乳腺MRI扫描提供的三维影像空间信息。在MRI图像上，尤其当起始和终止帧上肿瘤较小且灰度与周围组织较接近时，该方法的实验结果并不理想，甚至无法分割出病灶。一些结合帧间相关性的文献^[12-13]，虽然结合了帧间信息，但是需要医师手动标注，不适合处理大量数据。

技术实现要素：

针对现有技术的上述不足，本发明提供一种可以不依赖人工标注，自动准确分割起始和终止帧上乳腺病灶的乳腺MRI病灶的分割方法，本发明根据乳腺MRI图像的特点，改进现有的C-V水平集算法，将基于帧间相关性的SLIC0超像素和C-V水平集相结合，粗分割和细分割相结合，在粗分割部分采用SLIC聚类算法；在细分割部分结合改进的C-V水平集模型，更加快速准确地确定肿瘤轮廓，提高分割精度；在二维分割过程中结合三维帧间相关信息，提高分割的稳定性和灵敏度，为显示病灶的三维立体结构打下基础。技术方案如下：

一种基于帧间相关性的乳腺MRI病灶的自动分割方法包括：

A.读取MRI图像；

B.预处理图像，包括以下两个步骤：

C.对预处理后图像粗分割，确定病灶的初始轮廓，包括以下步骤：

D.采用改进的C-V水平集模型方法对步骤C得到的粗分割图像I(x,y)进行细分割，在粗分割轮廓基础上进一步细化肿瘤轮廓，改进之处如下：

1)C-V水平集模型是通过最小化能量函数，确定最终分割轮廓的，但其能量函数模型中只利用全局信息，而无法正确分割亮度分布不均匀图像，在C-V水平集模型的能量函数中加入局部能量项，用局部的灰度的统计信息来提高图像分割的能量，用一个窗口大小为k×k的均值卷积算子，对图像I(x,y)做卷积差值计算，以使图像的灰度分布差异减小，在细分割步骤中准确找到肿瘤轮廓；

2)为提高C-V水平集模型的收敛速度，以提高细分割效率，在C-V水平集模型的能量函数中引入能量惩罚项，对依据闭合曲线C构造的水平集函数φ做梯度计算，并对梯度值的平方再积分，以加快计算速度；

3)为控制零水平集函数的光滑度，将进化曲线的长度补偿项添加到C-V水平集模型的能量函数中，对Dirac函数与水平集函数的梯度乘积做积分运算，避免在细分割结果中出现孤立的小区域；

E.将步骤D得到的细分割结果，结合帧间相关性，进行优化，包括以下步骤：

1)计算细分割后得到的每一帧的疑似肿瘤区域面积，自动选取肿瘤面积最大帧作为迭代关键帧，得到关键帧的肿瘤轮廓R；

2)将R作为初始轮廓，利用上述改进的C-V水平集模型分别对关键帧的上、下一帧图像进行分割调整，得到更精确的分割结果，再以此结果为初始轮廓向前、后迭代，直到检测不到更细小的轮廓，即可认为该帧中不包含肿块。

作为优选实施方式，

B中的预处理图像，包括以下两个步骤：

1)截取ROI，取出包含乳腺的最小图像矩形区域；

2)采用顶帽(top-hat)运算进行形态学增强，得到预处理后图像。

对预处理后图像粗分割，确定病灶的初始轮廓，包括以下步骤：

1)超像素分割：采用SLIC聚类算法进行超像素分割；

2)可疑区域筛选，采用将相近灰度值的超像素区域合并，只保留灰度值偏高的区域作为可疑区域的方法，进行筛选，得到粗分割后的图像。

采用本发明，对预处理后的乳腺MRI图像进行超像素粗分割和C-V细分割之后，再结合帧间相关性进行分割和优化，相比于单一的分割方法，准确率更高，即使MRI图像起始和终止帧上的病灶很小、病灶灰度值和背景灰度值比较接近的情况，也能较好的分割出病灶，为后续病灶的三维立体显示，打下了良好基础。同时，整个分割过程为全自动分割，减少了人工交互。

附图说明

图1：基于帧间相关性的分割算法系统框图。

图2：由磁共振扫描仪获得的乳腺MRI图像。

图3：各步骤产生的结果：3(a)为自动截取的ROI图像；3(b)为形态学处理后的结果；3(c)为SLIC0超像素分割结果；3(d)为水平集细分割结果。

图4：本发明对乳腺MRI病灶图像一个完整序列的分割结果，每张小图下面的数值是其在高分辨MRI序列中实际的图片序号。

具体实施方式

参照图1所示，其中包含以下执行步骤：首先读取MRI图像10；然后对所得图像进行预处理20；接下来对预处理后图像粗分割30，确定病灶的初始轮廓；对粗分割后的图像细分割40，细化病灶边缘；最后采用基于帧间相关性的分割50，提高分割准确率。

在上述步骤中读取乳腺MRI图像10，获得的图像如图2所示。上述图像的采集装置为利浦Intera Achieva 1.5T磁共振扫描仪。采用轴位快速容积采集动态增强脂肪抑制序列(dyn_eTHRIVE)扫描，成像的相关参数为：重复时间TR＝4.4ms，回波时间TE＝2.2ms，翻转角FA＝12°。层厚2mm，FOV为100×100cm，图像灰度均为12位，矩阵大小为352×352(单位：像素)，重复扫描8次，每个动态扫描时间为60s。

步骤20对得到的MRI图像进行预处理，具体步骤如下：

自动截取ROI21：由于乳腺病灶的尺寸相对整幅影像的尺寸来说比较小，所以要自动选择出包含病灶的ROI，尽量使病灶位于ROI的中心。选取ROI的操作使分割对象更具有针对性，可以排除心脏等具有明显增强的组织对分割过程的影响，降低分割过程的复杂性。

由于同一设备得到的乳房的位置基本不变，所以可以设定一个矩形窗，包括整个乳房图像，去掉心脏等部分。在自动选取ROI的过程中，ROI的矩形窗大小可以根据病灶实际大小调节。为方便统计分析分割结果，除个别巨大病灶外，其余的病灶本发明统一采用大小为80×80像素的ROI窗，该尺寸远大于一般病灶的大小。

形态学增强22：

与普通图像相比，MRI初始图像容易受到影像设备和生成条件或者其他成份的影响，图像的质量可能会出现退化的状况，甚至有可能出现伪迹。为了解决这些问题，必须对图像进行预处理，去除图像的噪音，增强图像的对比度。本发明应用了一种双重形态学方法。可以有效的抑制背景区域，并突出可能区域，为之后的分割检测打好基础。

形态学处理方法包括开闭运算和顶帽运算。形态学开运算用于删除不包含结构元素的区域，使得目标对象的轮廓平滑，断开狭窄的链接，去除细小的突起。形态学的闭运算主要是用来填充物体内细小的空洞、连接临近的物体的。选择合适的结构元素，对图像先开运算后闭运算能够有效去除图像噪声。结构元素的选择直接决定形态学去除噪声的效果。

形态学处理方法中另外一种非常常用的方法是顶帽(top-hat)运算，它能够增强图像细节。顶帽运算将原始图像减去经过开运算的结果图像。顶帽运算能够增强图像的阴影细节，检测出图像的峰值信号。顶帽运算公式如下所示：

这种增强方法的具体步骤如下：

f(i,j)代表原图，B₁代表所选用的第一个结构元素，则首先对原图进行一次Top-hat运算：

此时r₁(i,j)就是原图与原图和结构元素B₁开运算后图像的差分图像。这里要选择结构元素B₁稍大于肿块的尺寸，那么通过这一步就能够去除图像中尺寸大于肿块的腺体组织，保留尺寸与肿块相当以及小于肿块的区域。

若另B₂代表所选用的第二个结构元素，为了去除尺寸很小的噪声区域及其它小尺寸的无用区域，还需要对得到的结果再进行一次开运算：

这种双重的形态学运算，分别去除了小尺寸的噪声杂质和大尺寸的腺体组织，仅保留图像中肿块可能大小范围内的区域。这种方法不仅可以突出疑似肿块区域，还可以免去在一些基于阈值方法中为避免干扰而在预处理中去除胸大肌的步骤。

粗分割步骤30，主要包括以下两个操作：

超像素分割31，利用像素之间特征的相似程度将像素分组，可以获取图像的冗余信息，在很大程度上降低了后续图像处理任务的复杂度。本发明采用简单的线性迭代聚类即SLIC聚类算法的零参数版本进行粗分割。

1)初始化种子点：假设图像有N个像素点，预分割为K个相同尺寸的超像素，那么每个超像素的大小为N/K，且相邻种子点的距离近似为为了避免种子点处在图像的边缘位置而对后续的聚类过程造成干扰，需要将种子点在以自身为中心的3×3窗口内移动到梯度值最小的位置，同时为每个种子点分配一个单独的标签。

2)相似度衡量：对于每个像素点，分别计算与之相邻的种子点之间的相似程度，将最相似种子点的标签赋给该像素；不断迭代该过程，直到收敛。定义(l_k,a_k,b_k,x_k,y_k)为在Lab空间里任意点的五维坐标，(l_i,a_i,b_i,x_i,y_i)为种子点的坐标值，则相似度的衡量关系为：

D_S是Lab距离和空间像素距离在网格间距上的归一化距离的和。d_lab代表Lab颜色距离，d_xy代表空间距离，S是类内最大空间距离.最大Lab颜色距离既随图片不同而不同，也随聚类不同而不同，所以取一个取值范围为[1,40]的固定常数m代替。m越大，空间距离相似性在聚类过程中的影响就越大。为了提高算法的运算速度，SLIC算法在对每个种子点聚类时，只在以种子点为中心的2S×2S区域内搜索相似像素点，而不是在整张图像中寻找。

SLIC0算法是SLIC的零参数版本。在SLIC算法中，m和S均为单一常量，而SLIC0算法则是动态地归一化每次迭代得到的最大空间和颜色距离，作为下一次迭代的空间和颜色距离。改进后的距离和公式为：

m_lab和m_xy为前一次迭代中得到的最大颜色和空间距离。由上述公式可知，SLIC0不再需要设置参数m和S，算法将根据图像不同区域的纹理复杂度，动态地调整合适参数值，使得SLIC0算法在保持速度不变的情况下，所产生的超像素更加规整统一，有利于减少后续处理时超像素形状差异对算法的影响。

可疑区域筛选32，SLIC0得出初始分割结果后，将对所得的超像素进行筛选，留下可疑区域。由于背景区域超像素的灰度值可见明显低于肿块区域，且形态学增强后背景区域超像素的灰度值明显低于原图像相同区域的灰度值，而且肿块区域超像素形状更规则，类圆度较高。因此可以将相近灰度值的超像素区域合并，只保留灰度值偏高的区域作为可疑区域。

对粗分割后的图像细分割40，采用C-V水平集方法细化肿瘤轮廓。Chan和Vese提出的C-V模型是对Mumford-Shan能量泛函的优化，利用内外灰度均值来促进水平集曲线的演化，可分割无边缘模糊、无梯度意义的图像，且对噪声不敏感。该模型假设定义域为Ω的图像I(x,y)被闭合曲线C划分为目标inside(C)和背景outside(C)两个同质区域。其能量泛函表示如下:

其中μ,ν是各项系数，均取非负值；c₁,c₂分别为目标和背景区域的均值。φ设为依据闭合曲线C构造的水平集函数，Ω为闭合曲线C的长度，▽φ为φ的梯度，得到水平集函数表达式为：

其中，定义正则化Heaviside函数H(φ(x,y))和Dirac函数δ(φ(x,y))为：

对轮廓C求导，得到以水平集函数φ表示的偏微分方程为

φ(0,x,y)＝φ₀(x,y)

其中：t表示时间，为的梯度产生的单位向量场散度，其意义是曲面运动时单位体积的改变率。c₁和c₂计算式为

C-V模型主要解决没有明显边缘的图像分割问题，但其存在两种主要缺陷：1)只利用了全局信息，无法正确分割亮度分布不均匀的图像。2)在计算中效率不高，所需时间较长。为克服这两种主要缺陷，对C-V模型进行改进。

基于C-V模型只利用到了全局信息，而无法正确分割亮度分布不均匀图像这个问题，给C-V模型中的加入局部项。局部能量项用到局部的灰度的统计信息来提高图像分割的能量，局部能量函数E(d₁,d₂,C)描述如下:

其中，g_k为窗口大小为k×k的均值卷积算子。*为卷计算子，d₁,d₂分别是差值图像(g_k*I(x,y)-I(x,y))在目标和背景区域的灰度均值。

C-V模型的第二个主要缺陷是在计算过程中效率不高，所需要的计算时间相当长，相当耗时，解决方法是引入能量惩罚项P(φ)：

▽为求梯度符号，▽φ为对φ求梯度。

另外，为了控制零水平集函数的光滑度及避免在分割结果中出现孤立的小区域，将进化曲线的长度补偿项添加到调整能量项中。长度补偿项定义如下：

正则化Heaviside函数H(φ(x,y))和Dirac函数δ(φ(x,y))的定义同上。

改进后的C-V模型，其基于水平集的能量泛函为：

其中

本发明中参数设置为:λ₁＝1，λ₂＝1，μ＝2，υ＝1.5，α＝0.04，ε＝1.0。改进后的C-V算法具有无重新初始化的曲线演化模型优点，收敛速度更快，且收敛性更好，分割结果更加贴近肿瘤轮廓。

利用帧间相关性的分割，首先要选取中间帧。一般方法采用的中间帧大多需要医师手动标注，分割过程并不是全自动。为降低人机交互的需要，提高检测效率，本文提出一种自动选取中间帧的算法。根据细分割后的结果计算每一帧中疑似病灶的面积，将包含病灶的面积最大的一帧认定为中间帧，来约束前后向的分割。这是因为肿瘤面积较大的帧，肿瘤轮廓相对清晰，更容易得到精确的分割结果作为迭代的初始轮廓。

首先，对细分割的结果进行优化，除去错误的分割区域。一般乳腺肿瘤的最大直径大约为30mm，在MRI图像上，大约占250像素。而错误分割部分的区域面积，一般小于最大肿瘤面积的1/5，因此，除去每帧图像上面积小于最大肿瘤面积1/5(即50像素)的连通区域，即可在保留最大肿瘤的同时，除去错误分割的部分。

计算细分割后得到的每一帧的疑似肿瘤区域面积，自动选取肿瘤面积最大帧作为迭代关键帧，得到关键帧的肿瘤轮廓R。

将R作为初始轮廓，利用改进的C-V模型分别对关键帧的上、下一帧图像进行分割调整，得到更精确的分割结果，再以此结果为初始轮廓向前、后迭代。

直到检测不到更细小的轮廓，即可认为该帧中不包含肿块。