2. 东北大学 医学影像计算教育部重点实验室,辽宁 沈阳 110169
2. Key Laboratory of Medical Image Computing, Ministry of Education, Northeastern University, Shenyang 110169, China
MRI(magnetic resonance imaging)在脑功能疾病诊断方面被广泛应用于临床医学和医学研究, 其中, 脑组织的提取是脑功能疾病诊断中的重要环节.由于MRI影像通常存在偏移场效应、局部体效应[1-3], 导致MRI影像组织边缘模糊, 灰度分布不均匀;此外,MRI脑影像还存在个体差异性.这些问题的存在, 使得脑组织的提取成为医学影像分析领域极具挑战性的课题.目前, 在脑组织提取方面, 文献[4-5]结合各向异性扩散方法、形态学方法和阈值分割方法, 最终实现了脑组织的提取.分水岭技术是一种基于区域的分割方法.文献[6]提出一种基于水平集方法的分水岭变换和模糊C值聚类的脑图像分割算法, 该方法通过灰度梯度判断子区域的连通性, 并采用水平集方法得到最终的分割结果.文献[7]提出三维分水岭方法, 采用简单的合并规则避免过分割问题.文献[8]提出可变形点阵模型, 但是不能直接处理三维脑影像.文献[9-10]提出一种基于活动轮廓模型来提取脑白质、脑灰质的方法, 但计算量较大.文献[11]提出基于模板的分割方法, 该方法迭代地将模板与脑组织进行匹配, 最终获得脑组织, 与区域分割的方法相比, 该方法鲁棒性较好, 且人工交互少, 对影像质量不敏感.
BET(brain extraction tool)[12]算法是当前常用的基于可变形点阵模型的脑组织提取方法.该方法首先计算MRI脑影像的灰度直方图, 通过灰度直方图估计三个值, 分别是区分脑组织和非脑组织的影像灰度的灰度阈值、影像灰度最大值和最小值;然后粗略估计脑组织的重心, 并依据脑与非脑的灰度值获得初始脑组织;最后在脑组织内通过三维三角面片构建脑初始表面, 每个三角面片建立切向力和平滑力, 在三角面片上两个力的驱使下使得初始面保持一定的距离和平滑性, 直到脑表面足够平滑且稳定下来, 分割结束.该算法具有较高的鲁棒性和准确性, 同时速度也非常快, 是一种非常好的脑组织提取算法.但该算法对初始脑重心位置选取依赖较大, 处理真实全脑图像时往往对重心的估计不够准确, 效果不好, 并且在真实全脑图像中, 脑组织与相邻组织的灰度值相近, 脑组织提取时很容易发生泄漏或过度分割[13].
针对上述问题, 本文引入图像梯度信息, 进一步改进脑重心的选取, 在脑表面形变过程中引入边缘力以减少脑组织在边界的泄漏和过度分割.
1 实验方法改进的BET脑组织自动提取算法主要由两部分组成:①参数估计和脑表面初始化;②脑表面形变力构建.首先在参数初估计中, 引入图像梯度信息, 提高脑重心(center of gravity, COG)的计算准确性,其次构建新的脑表面形变力, 在垂直于脑表面切线的扩张力中引入边缘力, 很好地抑制了脑组织的边界泄漏和过度分割问题.
1.1 脑重心估计和脑表面初始化为了估计脑重心, 首先计算区分脑与非脑组织的阈值t.在原始BET算法中t=t2+0.1(t98-t2), 其中t98是脑影像直方图中灰度值等于灰度分布98%处的影像灰度值, t2是脑影像直方图中灰度值等于灰度分布2%处的影像灰度值.根据先验知识, 整个脑影像中灰度值大于t的位置属于脑组织, 其他位置为非脑组织部分.获得脑组织区域后, 将计算脑组织区域像素位置的平均值作为脑重心COG; 但对全脑数据而言, 在z轴方向, 由于阈值t无法准确地对头颈部灰度值与脑部组织灰度值进行区分, 导致在z轴方向的COG坐标均值不准确, 在x和y轴方向的坐标均值比较准确.图 1为MRI脑图像.
针对z轴坐标值不准确的问题, 引入图像梯度信息以提高脑重心的计算准确性.沿着z轴分别向上和向下求取每个像素点的梯度值zu(i)和zd(i), 并取得z轴上方向最大梯度值zu和z轴下方向最大梯度值zd.这里认为zu和zd所在位置是z轴方向上的脑组织边界, 则脑重心在z轴方向上的坐标为zu和zd所在位置的1/2处, 如式(1)所示:
(1) |
假设脑组织是一个三维球体, 半径为R.将R作为初始脑的半径:
(2) |
图 2是脑表面模型示意图.初始脑表面是一个二十面体, 为了使表面上每个顶点到二十面体中心的距离尽可能相等, 将二十面体中的每个小三角分成4个更小的三角.在新的表面中, 原来二十面体中的每个顶点都有5个或是6个相邻的顶点.按照这样的方法, 将表面上新的三角面片继续分解成更小的三角, 直到这个二十面体近似于一个球面为止.
BET算法的主要思想是通过3种力将相互作用的轮廓点推到脑组织边缘, 从而得到脑组织.这3种力分别是平行于脑表面切线的拉力u1, 该力使脑表面点保持间距; 垂直于脑表面切线的平滑力u2, 该力和脑表面的曲率有关, 使得轮廓平滑; 垂直于脑表面切线的扩张力u3, 该力使得脑表面向外或向内演化(变形).
脑组织轮廓演化的过程就是轮廓线上的轮廓点根据受力情况更新轨迹的过程.最后受力接近平衡时, 轮廓点在一定的小范围内运动, 这时的轮廓点连线为脑组织轮廓.图 3描述的是在二维平面中求向量s及其两个分解向量sn和st.如图 3所示, 假设有顶点O, A和B是与O相邻的两个点, A和B中间的位置为D, 用向量s代表顶点O到点D的向量, s的分解向量为sn和st,“^”表示单位向量:
(3) |
(4) |
sn和st两个正交向量的变化是推动顶点O演化的三个力u1, u2, u3更新的基础, u1由脑表面的切向量st及其权值组成, u1限制顶点只在脑表面中移动.u2和u3的方向和sn的方向平行.
(5) |
(6) |
式中:
原始BET算法中, 脑组织在边界处会出现泄漏和过度分割, 为此在u3引入边缘力f:
(7) |
式中:d是当前点O到脑重心的距离; R为半径估计值; gi是当前点的梯度值, gμ是当前点周围26个邻域内所有点的平均梯度值, gm是当前点所在的边界梯度值, 即gm =zu或者zd.yμ为当前点周围26个邻域内所有点的平均灰度值, yi为当前点的灰度值, η是一个极小常数, 则
(8) |
式中:Imin表示从顶点O沿着平行于它的局部表面法向量方向,向里前进距离d1获得的最小灰度值; Imax表示从顶点O沿着平行于它的局部表面法向量方向,向里前进距离d2获得的最大灰度值.通常, d1=20 mm, d2= d1/2;t1=(Imax-t2)bt+t2, bt的默认值是0.5.
最终的形变力u为
(9) |
式中a1, a2, a3为比例系数.
2 结果与讨论北京宣武医院为本文的实验提供健康人活体脑部MRI影像数据, 这些数据符合DICOM 3.0标准, 脑影像的分辨率是512×512, 层数为192层, 层厚为1 mm, 像素间距为0.5 mm.
采用本文提出的算法对MRI脑数据进行测试.图 4是本文算法的实验结果.从图中可以看出, 使用本文算法提取的脑组织结构清晰完整, 表面光滑, 没有明显的脑组织泄漏和过度分割现象, 较好地完成了脑组织的提取工作.
选取三组数据对本文算法与原始BET算法进行对比, 实验结果如图 5所示:每组图中红色部分标定为脑组织提取结果, 蓝色标记为脑组织过度分割位置, 黄色标记为脑组织泄漏位置.
从图 5可以看出, 采用BET方法提取脑组织, 会出现一些脑组织泄漏和过度分割现象.而采用本文方法能很好地解决这一问题, 从而能较好地提取脑组织.
图 6是BET算法和本文算法的局部对比, 从图中可以看出, 本文算法在头颈部区域、脑组织边缘区域的提取结果都优于原始BET算法.从实验结果可以看出, 与BET方法相比, 本文方法能够有效提取脑组织, 并能够有效减少脑组织的泄漏和过度分割.
1) 本文提出的MRI脑组织自动提取算法能够提高初始脑重心的定位精准度, 也能有效解决脑组织泄漏和过度分割的问题.
2) 本文算法能够自动获得更加准确的脑组织提取结果, 特别是与BET算法相比, 能够获得更准确的脑组织边缘提取结果.
3) 实验结果表明本文算法具有有效性和准确性, 优于BET方法.
[1] |
Cheng H D, Jiang X H, Sun Y.
Color image segmentation:advances and prospects[J]. Pattern Recognition, 2001, 34(12): 2259–2281.
DOI:10.1016/S0031-3203(00)00149-7 |
[2] |
Min R, Wu G, Cheng J, et al.
Multi-atlas based representations for Alzheimer's disease diagnosis[J]. Human Brain Mapping, 2014, 35(10): 5052–5071.
|
[3] |
Tijms B M, Yeung H M, Sikkes S A.
Single-subject grey matter graph properties and their relationship with cognitive impairment in early-and late-onset Alzheimer's disease[J]. Brain Connectivity, 2014, 4(5): 337–346.
DOI:10.1089/brain.2013.0209 |
[4] |
Atkins M S, Mackiewich B T.
Fully automatic segmentation of the brain in MRI[J]. IEEE Transactions on Medical Imaging, 1998, 17(1): 98–107.
|
[5] |
Meegama R G N.
Fully automated peeling technique for T1-weighted, high-quality MR head scans[J]. International Journal of Image and Graphics, 2004, 4(2): 141–156.
DOI:10.1142/S0219467804001348 |
[6] |
Saikumar T, Shashidhar B, Harshavardhan V, et al.
MRI brain image segmentation algorithm using watershed transform and kernel fuzzy c-means clustering on level set method[J]. International Journal on Computer Science and Engineering, 2011, 3(4): 1591–1598.
|
[7] |
Hahn H K, Peitgen H O. The skull stripping problem in MRI solved by a single 3D watershed transform[C]//International Conference on Medical Image Computing and Computer-Assisted Intervention. Pittsburgh, 2000: 134-143.
|
[8] |
贾迪, 杨金柱, 张一飞, 等.
序列磁共振颅脑影像的脑组织自动提取方法[J]. 仪器仪表学报, 2011, 32(8): 1781–1787.
( Jia Di, Yang Jin-zhu, Zhang Yi-fei, et al. Automatic brain tissue extraction approach of serial magnetic resonance head images[J]. Chinese Journal of Scientific Instrument, 2011, 32(8): 1781–1787. ) |
[9] |
Atkins M S, Mackiewich B T.
Fully automatic segmentation of the brain in MRI[J]. IEEE Transactions on Medical Imaging, 1998, 17(1): 98–107.
DOI:10.1109/42.668699 |
[10] |
贾迪, 杨金柱, 张一飞, 等.
自适应脑组织影像分割[J]. 吉林大学学报(工学版), 2012, 42(1): 161–165.
( Jia Di, Yang Jin-zhu, Zhang Yi-fei, et al. Self-adapting segmentation for brain tissue[J]. Journal of Jilin University(Engineering and Technology Edition), 2012, 42(1): 161–165. ) |
[11] |
Dale A M, Fischl B, Sereno M I.
Cortical surface-based analysis I:segmentation and surface reconstruction[J]. Neuroimage, 1999, 9(2): 179–194.
DOI:10.1006/nimg.1998.0395 |
[12] |
Smith S M.
Fast robust automated brain extraction[J]. Human Brain Mapping, 2002, 17(3): 143–155.
DOI:10.1002/(ISSN)1097-0193 |
[13] |
Yang J, Tan W, Ma S, et al.
Automatic MRI brain tissue extraction algorithm based on three-dimensional gray-scale transformation model[J]. Journal of Medical Imaging & Health Informatics, 2014, 4(6): 907–911.
|