2. 清华大学 生物医学工程系, 北京 100084, 中国;
3. 美国纽约州立大学石溪分校, 石溪镇NY 11794, 美国
2. Department of Biomedical Engineering, Tsinghua University, Beijing 100084, China;
3. Stony Brook University, Stony Brook NY 11794, USA
膀胱癌是最常见的泌尿系统恶性肿瘤，在中国和美国其发病率和死亡率都居于恶性肿瘤前列。部分膀胱癌中晚期患者需要进行膀胱切除手术，而早期诊断发现后可以通过保留膀胱手术方案进行治疗[1]。目前，世界上通用的光学膀胱镜检查对于膀胱癌早期诊断与术后复查有重要的意义，但其缺点在于检查过程有创，甚至存在膀胱穿孔的危险[1]。结合计算机技术和医学影像的虚拟膀胱镜技术成为一种有效的替代方法。通过定量化测量膀胱壁厚度等参数，对可能存在膀胱肿瘤的区域进行检测和定位[2]。
段侪杰等研究[1]表明，T1加权的核磁共振(magnetic resonance,MR)图像更适合用于虚拟膀胱镜。MR图像中的膀胱壁分割是虚拟膀胱内窥镜的构建过程中关键的一环，分割结果的精确与否直接影响虚拟膀胱内窥镜的效果。不过MR图像信噪比一般比较低，并且存在运动伪影，且边缘不清晰，影响分割的精度，影响虚拟内窥镜的效果。所以在分割前需要对MR图像进行预处理，提高信噪比，减弱伪影，锐化边缘。改善图像信噪比一个简单的方法就是改变MR扫描中的平均参数[3]。为平缓膀胱壁边缘伪影和提高膀胱壁信噪比，使用短时采集，将采集中运动(长时)变成采集间运动(短时)是有效方法。
在众多医学图像分割方法中，水平集方法受限于几何约束，结果的精确度达不到虚拟内窥镜的要求。而图割方法[4]在精度上要远高于水平集方法[5]，在细节方面有着非常卓越的性能，应用在虚拟内窥镜上更具优势。
本文把几套短时采集核磁图像中的1套作为基准，将其它短时图像向它配准，然后将配准结果与基准图像求平均，再对平均结果利用图割方法进行分割。
1 原理 1.1 配准图像配准就是将一幅图像(浮动图像)映射到另一幅图像(基准图像)上同源点的空间变换过程。
从配准的基本框架图(图 1)中可以看到，浮动 图像经过变换得到的数据通过插值器(interpolator) 的插值，得到与基准图像像素位置对应的图像，插值后的结果经过测度(metric)与基准图像进行匹配度的评估，优化器利用评估结果调整变换的参数，完成一次迭代。根据优化器的参数设置，迭代结束后将得到使得测度值最优化的变换参数，此时变换得到的数据经过插值得到的图像即最终配准结果。在实验中，采用仿射变换粗配准，分层B样条变换细配准的策略，测度采用均方差，采用边界约束BFGS优化算法(L-BFGS-B)，插值器采用线性插值。
|
| 图 1 配准框架图 |
同一人利用一台MRI机器在相同环境条件下采集到的几套数据，对应人体组织的相同位置应采集到相同的灰度信息。经过配准求平均后，相同的信息增强，随机的噪音和伪影减弱。
1.2 分割由于虚拟内窥镜对分割精度要求较高，分割算法采用图割方法中的闭集模型[6]。图割方法建立在图论中解决最大流-最小割问题的一套方法。最大流-最小割问题是网络流中的一种最优化问题，这一问题的解决最终可以将网络中所有节点分成2个集合。图割方法则是将图像分割问题转换成最大流-最小割问题，求解该问题，最终得到分割结果。转换的方法是将图像映射成一个图(graph)，图中每个节点与图像中的像素或者亚像素对应，节点之间根据一定的规则连接，每个连接的边被赋予图像对应位置信息的代价值，最后求解这个图的最大流-最小割问题。闭集模型是一种节点具有代价值的有向图，该模型由若干列(column)组成，列上有相同的等间距分布的节点，列与列间距相等，互相平行。列内节点从上至下连接，相邻列之间节点按照固定间隔从上向下连接。这种模型存在一个最小闭集，即节点代价值总和最小的集合，最小闭集中节点指向的节点属于最小闭集。求解最小闭集问题后，有向图被分为上下2个部分，交界面与每一列只有一个交点，如果图像被映射为该模型，那么交界面就是最终图像分割结果。这种模型可以转换成边具有权重的新图，原图的最小闭集问题求解则等价[7]于新图的最大流-最小割问题求解[8]。将图像映射成一个闭集模型，解决其最小闭集问题，最后得到分割结果[6]。
2 实验方法 2.1 配准 2.1.1 模拟实验仿真数据(图 2)为199像素×199像素的图像，圆心坐标(100,100)，外圆半径66像素，膀胱径向厚度 10像素，以(140,140)为圆心，半径为8像素画实心圆，模拟小肿瘤(图 2箭头所指处)。该图像作为基准原始图像，通过平移和形变，将原始图像位移变形得到6个浮动图像。图像上有若干条纹状突起，代表伪影(图 2中肿瘤以外的突起均为伪影)，膀胱壁和伪影以及息肉的灰度均为100，最后加上Gauss噪声，方差为5，膀胱壁内部噪声方差为10。模拟数据的灰度设定近似真实数据的灰度分布。
|
| 图 2 用于实验的模拟图像 |
利用ITK(Insight Segmentation and Registration Toolkit)函数包在VS2010下编程，测度选择均方差，优化算法选择L-BFGS-B，插值方法选择线性插值，变换方法选择先仿射变换进行粗配准，后采用2层B样条变换，第一层在x方向上为5个网格点(y和z方向上网格点个数按照像素数比例进行自适应调整)，第二层在x方向上为12个网格点，该层控制点间隔大概为伪影尺寸的2倍。
2.1.2 临床MR图像配准实验数据由美国纽约州立大学石溪分校医学健康中心的3.0 T核磁共振全身扫描仪采集，采用三维T1高分辨率各向同性容积激发脉冲序列，视野域21 cm×21 cm，TE3.873 ms，TR6.246 ms，像素0.937 mm×0.937 mm，层厚2 mm，层间距1 mm，反转角10°，信号采集为1，NSA为1。
实验利用一个被试的连续7套短时采集数据，按照采集时间先后顺序选取中间1套数据作为基准，剩余其它数据作为待配准数据，按照上述模拟实验步骤进行配准。由于多次短时采集存在时间跨度，浮动图像中的膀胱可能与基准图像中的膀胱运动变形差异较大，利用仿射变化进行粗配准可以减少整体配准的时间，而且可以提高后续B样条变换配准在变形较大时的性能。多层B样条变换也可以减少B样条变换配准的时间且提高配准效果。 最后一层网格点数在x方向上仍然选择为12。如果网格过于密集，会使得浮动图像在配准后被引入基准图像的伪影信息； 如果网格点过于稀疏，过少的网格点将难以处理局部的变形。
2.2 分割利用配准结果平均后的数据进行分割，基于一个圆心沿径向进行采样[9]，这样的采样方式对于充分膨胀并且没有肿瘤的膀胱MR图像是适用的，因为符合闭集模型[6]。但是对于膨胀不完全或者存在较大肿瘤的膀胱，一个采样列会穿过期望分割曲面2次，不符合最小闭集模型的基本假设，是不可能得到准确的分割结果，所以在第一层利用一个初始曲线，这个初始曲线尽量接近待分割区域，并且该曲线上每个点的沿采样方向往内侧的采样列只穿过膀胱壁内边缘一次。这样的初始曲线利用其它算法得到，或者直接手动勾勒得到。初始曲线上所有像素点经过一定的插值之后，沿法向往膀胱内侧采样，然后按照相关方法[9]进行构图(graph)，利用Boykov的函数包[8]对该图求解最小割。采样方向选用法向，但由于沿法向采样非常容易造成采样列发生交叉现象，需要将最小割得到的点重新连接，得到的结果边与边之间不存在交叉。当然可以通过一些方法选定在一定距离内不存在交叉的采样方向，并且采样区域不存在重复，这样更加利于分割，但是这样理想的直线采样方式实际上很难实现，最好的方式是使用电场线的方法，采样路径由直线变成均匀分布并且没有交叉的曲线。
3 实验结果 3.1 配准分别对模拟数据和真实临床MR数据的配准增强结果和基准图像的比较，通过绘制灰度分布曲线和测量信噪比来进行评估。
3.1.1 灰度分布曲线图 3给出了模拟图像配准平均前后效果图，图 4给出了图 3的2幅图像上某一行像素的灰度值折线图。从图 3可以清楚地看到，增强图像的信噪比和膀胱壁与膀胱腔间的对比度明显强于基准图像，膀胱壁附近的噪声得到很好的抑制，膀胱壁的分辨情况也远好于基准图像。图 5为临床MR图像配准前 后效果图，图 6为图 5的2幅图像上某行像素的灰度值折线图，可见配准平均的方法对抑制临床MR图像噪声有效。
|
| 图 3 模拟图像配准前后比较 |
|
| 图 4 模拟图像配准前后灰度值折线图比较 |
|
| 图 5 MR图像配准前后比较 |
|
| 图 6 MR图像配准前后灰度值折线图比较 |
为了测量平均之后的信噪比增益，选择核磁图像中膀胱腔体作为信号区，人体外的空气区域作为无信号噪声区。基于Henkelma[10]和Kaufma[11]的方法，利用式(1)进行计算。计算结果见表 1。
模拟图像配准后信噪比增大了2.26倍，MR图像配准后信噪比增大了1.17倍。有
| $SNR = \frac{{{u_s}}}{{1.5{\delta _N}}}$ | (1) |
其中 ：us表示信号的平均值，σN表示噪声的标准差。
3.2 分割分别分割模拟配准增强数据和MR配准增强数据。模拟数据采用径向采样方式[9]进行全局采样，这由于模拟数据的原始图是中心对称的，这与闭集模型[6]是对应的，从中心展开后，每一列与期望分割表面只有一个交点。对于真实数据，利用 2.2 中提到的方法获得采样方向，并在每层切片分割之后更新采集的方向。
3.2.1 模拟数据图 7为模拟图像配准前后结肠内壁分割结果，从图中可以看到，原图由于存在噪声和伪影，分割结果并不理想。配准增强后的分割完全不受伪影影响，虽然伪影肉眼依然可以识别，但是平均后淡化。同时原来图像中由于受到噪音干扰影响的分割结果(图 7箭头所指处)，在配准增强之后得到改善。从图 7可以非常直观地看出，配准增强可以对后续分割起到很好的帮助。
|
| 图 7 模拟图像配准增强前后分割结果 |
图 8 第一行为MR图像基准图第43层分割情况，第二行为配准增强图分割情况。图 8c是基于基准图像第43层切片手描初始曲线，初始曲线并不一定完全在膀胱壁上，但尽量贴近膀胱壁区域，后面其它层的分割依赖于这一层的分割结果，若第一层分割出现严重问题，后面的分割也会随之出现问题。由于设定的采样方向在一定距离内比较均匀且极少交叉，曲线上所有点的采样距离设定相同长度，初始曲线需要尽可能靠近膀胱壁内边缘，否则采样路径可能会无法经过膀胱壁内边缘。
|
| 图 8 临床MR数据分割情况 |
从图 8可以清楚地看到，由于噪音和伪影的存在，单次采集数据的分割结果容易最后收敛在内部黑色噪点处，在肿瘤区域也存在这样的噪点，本文中使用的分割算法存在将这些噪点或者伪影当作边缘的倾向，所以最终会影响分割的效果。在图 8b中左侧箭头所指处有一个黑色噪点，最终分割到该点附近区域，而图 8e中，配准增强图像中的噪声和伪影得到很好的抑制，最后分割结果正确。对于临近层，可以采用同样的思路，利用之前层的分割结果来作为初始曲线进行分割，图 8f是利用第43层分割结果对第44层分割的结果。
4 结 论本文对T1加权的膀胱MR图像提出了一套完整的从增强到分割的策略，旨在提高图像的信噪比并进一步提高对膀胱壁的分割精度，为构建虚拟膀胱镜及提高膀胱病灶的检出率打下基础。本文对多套MR单次激发短时采集数据进行配准增强，并提出图割算法的改进策略。使用仿真数据和真实MR数据进行实验，结果表明单次激发图像受运动伪影影响较轻，配准增强方法可以大大提升MR图像信噪比，使得膀胱壁有更好的边缘特性。在此基础上，使用改进的图割分割策略，可以得到精确的膀胱壁分割边缘。
| [1] | 段侪杰, 田珍, 梁正荣, 等. 基于医学影像的虚拟膀胱镜技术 [J]. 中国医学物理学杂志, 2010, 27(2): 1712-1715.DUAN Chaijie, TIAN Zhen, LIANG Zhengrong, et al. A review on virtual cystoscopy techniques [J]. Chinese Journal of Medical Physics, 2010, 27(2):1712-1715. (in Chinese) |
| [2] | Jaume S, Ferrant M, Macq B, et al. Tumor detection in the bladder wall with a measurement of abnormal thickness in CT scans [J]. IEEE Transactions on Bio-medical Engineering, 2003, 50(3): 383-390. |
| [3] | LIN Qin, LIANG Zhengrong, DUAN Chaijie, et al. Motion correction for MR cystography by an image processing approach [J]. IEEE Transactions on Biomedical Engineering, 2013, 60(9): 2401-2410. |
| [4] | Boykov Y, Veksler O, Zabih R. Fast approximate energy minimization via graph cuts [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2001, 23(11): 1222-1239. |
| [5] | Osher S, Sethian J A. Fronts propagating with curvature dependent speed: algorithms based on Hamilton-Jacobi formulation [J]. Journal of Computational Physics, 1988, 79(1): 12-49. |
| [6] | LI Kang, WU Xiaodong, CHEN Ziyi, et al. Optimal surface segmentation in volumetric images—A graph-theoretic approach [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2006, 28(1): 119-134. |
| [7] | Hochbaum D. A new-old algorithm for minimum-cut and maximum-flow in closure graphs [J]. Networks, 2001, 37(4): 171-193. |
| [8] | Boykov Y, Kolmogorov V. An experimental comparison of min-cut/max-flow algorithms for energy minimization in vision [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2004, 26(9): 1124-1137. |
| [9] | ZHANG Yijia, DUAN Chanjie, YE Datian, et al. The segmentation of MR bladder wall in 3D based on minimum closed set model [C]//2013 IEEE International Conference on Medical Imaging Physics and Engineering. Shenyang, China: IEEE, 2013: 319-323. |
| [10] | Henkelman R. Measurement of signal intensities in the presence of noise in MR images [J]. Medical Physics, 1985, 12(2): 232-233. |
| [11] | Kaufman L, Kramer D, Crooks L, et al. Measuring signal-to-noise ratios in MR imaging [J]. Radiology, 1989, 173(1): 265-267 |