CN101964118A - 一种血管内超声图像序列的三维分割方法 - Google Patents

Abstract

一种血管内超声图像序列的三维分割方法，用于提高图像序列的分割处理效率。其技术方案是：它首先对原始图像进行滤除噪声和抑制环晕伪像的预处理，然后获取IVUS图像序列的四个纵向视图，并从中提取出血管内腔边界和中-外膜边界，再通过将这些边界曲线映射到各帧IVUS图像中，得到横向视图中的初始边界，最后通过使能量函数最小，初始边界不断变形，最终得到各帧IVUS图像中的血管内腔和中-外膜边界。与传统方法相比，本发明具有以下优点：第一，可以利用整个图像序列的信息；第二，可同时完成对各帧图像的分割，实现对整个图像序列的并行处理，从而大大提高了处理效率，缩短了处理时间。

Description

一种血管内超声图像序列的三维分割方法

技术领域

本发明涉及一种血管内超声图像序列的三维分割方法，该方法可实现对IVUS图像序列的并行分割，属医学成像技术领域。

背景技术

血管内超声(Intravascular Ultrasound，IVUS)成像技术能够反映血管内腔的变化以及含斑块在内的管腔横断面结构、血管壁厚度、形态以及斑块成分等。血管壁在IVUS图像上呈现回声不同的三层环形结构，且图像序列中的前后帧之间非常相似，具有很强的相关性。

对IVUS图像进行分割，提取出血管内腔边界和中-外膜边界(包括可能存在的斑块边缘)，是定量分析和血管三维重建的重要步骤，分割质量直接决定定量测量和三维重建的精度。目前临床常采用手动分割，由富有经验的医生用鼠标在屏幕上逐帧勾画。由于IVUS图像序列常包括上千帧，因而该方法不但耗时，而且可重复性差，分割结果在很大程度上取决于操作者的经验和主观因素。现有的采用数字图像处理技术的自动分割方法都是在二维空间中进行串行处理，即对IVUS序列中的各帧图像进行逐帧分割，因而处理效率较低。

发明内容

本发明的目的是克服现有技术的不足、提供一种IVUS图像序列的三维分割方法，以提高图像序列的分割处理效率。

本发明所称问题是以下述技术方案实现的：

一种血管内超声图像序列的三维分割方法，它首先对原始IVUS图像进行滤除噪声和抑制环晕伪像的预处理，然后获取IVUS图像序列的四个纵向视图，并从中提取出血管内腔边界和中-外膜边界，再通过将这些边界曲线映射到各帧IVUS图像中，得到横向视图中血管壁的初始边界，最后通过使能量函数最小，初始边界不断变形，最终得到各帧IVUS图像中的血管内腔和中-外膜边界。

上述血管内超声图像序列的三维分割方法，具体采用以下步骤：

a、对原始IVUS图像进行预处理，包括滤波去噪和去除环晕伪像：

首先采用中值滤波和高斯平滑两种通用预处理方法，减少IVUS图像中的椒盐噪声和随机噪声，然后对各帧IVUS图像进行极坐标变换，得到其极坐标视图，再按照下式去除极坐标视图中的环晕伪像：

$I^{\′}(r,\mathrm{\θ})=\left\{\begin{array}{cc}{I}_{\mathrm{catheter}},& r\≤\mathrm{\ϵ}\×\mathrm{ImageWidth}/2\\ I(r,\mathrm{\θ}),& \mathrm{else}\end{array}\right.$

其中，I(r，θ)和I′(r，θ)分别为原极坐标视图和去除环晕伪像后的极坐标视图中像素点(r，θ)处的灰度值；r为像素点的极径；θ为像素点的极角；I_catheter是原IVUS图像中导管区域的像素灰度值；ImageWidth为以像素为单位的图像宽度；ε为权重参数；最后经过极坐标逆变换，即可得到直角坐标系下去除环晕伪像后的IVUS图像；

b、分别取IVUS图像序列的沿血管长轴方向的四个纵向视图，即垂直切面A、水平切面B、左对角线切面C和右对角线切面D；

c、从IVUS纵向视图中提取出血管内腔边界和中-外膜边界：

从纵向视图的中轴线开始，分别向左和向右逐行遍历四个IVUS纵向视图中的各像素，用I(i，j)表示坐标为(i，j)的像素的灰度值，I(i，j+1)表示坐标为(i，j+1)的像素的灰度值，若I(i，j+1)-I(i，j)≥η，η为阈值，则(i，j)为目标边界点，否则不是；将每行左右两部分的像素中，第一个符合上述条件的像素记为内腔边界点，第二个记为中-外膜边界点；其中，η为阈值；

d、将四个纵向视图中的边界曲线映射到各帧IVUS图像中，纵向视图中的边界线对应到横向视图中为边界点，依次连接各边界点，得到各帧IVUS图像中血管内腔的初始边界和中-外膜的初始边界；

e、各帧IVUS图像中的初始边界通过snake变形获得血管内腔边界和中-外膜边界：

将各帧IVUS图像中血管内腔初始边界和中-外膜初始边界作为snake模型的初始形状，将初始边界离散成由N个点组成的有序点集，则snake能量函数的离散表达式为：

$E=\underset{i=0}{\overset{N-1}{\mathrm{\Σ}}}[E_{\mathrm{int}}\left(i\right)+E_{\mathrm{ext}}\left(i\right)]$

$E_{\mathrm{int}}\left(i\right)=\mathrm{\α}\frac{|\stackrel{\‾}{d}-|c_i-c_{i-1}||}{\max d}+\mathrm{\β}\frac{|c_{i-1}-2c_i+c_{i+1}|}{2\max d},$ $E_{\mathrm{ext}}\left(i\right)=\mathrm{\γ}\frac{I(x_i,y_i)}{255}-\mathrm{\λ}\frac{|\▿I(x_i,y_i)|}{255\sqrt{2}}$

其中，E_int是内部能量；E_ext是外部能量；c_i(x_i，y_i)(i＝1，2，…，N-2)是第i个snake点，(x_i，y_i)是它的坐标；c_i-1和c_i+1分别是第i-1个和第i+1个snake点；

和maxd分别是相邻snake点之间的平均距离和最大距离；I(x_i，y_i)和

分别是坐标为(x_i，y_i)的像素的灰度和灰度梯度值；α，β，γ，λ∈[0，1]都是权重参数。

通过使能量函数E最小，snake模型从初始形状开始不断变形，最终停留在能量函数取得全局最小值的最优位置，即为目标边界。

上述血管内超声图像序列的三维分割方法，所述权重参数ε的取值区间为[0.1，0.35]，阈值η的取值区间为[10，20]。

上述血管内超声图像序列的三维分割方法，能量函数的全局最优化采用Williams贪婪算法完成。

本发明采用三维分割方法对IVUS图像序列进行处理，与二维分割方法相比具有以下优点：第一，可以利用整个图像序列的信息；第二，可同时完成对各帧图像的分割，实现对整个图像序列的并行处理，从而大大提高了处理效率，缩短了处理时间。大量临床图像数据的实验结果证明，该方法可自动、快速、可靠地完成IVUS图像序列的分割。

附图说明

下面结合附图对本发明作进一步详述。

图1是根据本发明方法的图像处理流程图；

图2是本发明对一帧IVUS图像的环晕伪像去除过程和结果；

图3是本发明获取IVUS图像序列的纵向视图示意图；

图4是本发明逐行遍历IVUS纵向视图像素的示意图；

图5是本发明的IVUS纵向与横向视图血管壁边界对应关系示意图；

图6是完成了血管壁边界提取的IVUS图像序列纵向视图；

图7是图6所示图像序列中4帧IVUS图像中血管壁的初始边界；

图8和图9是IVUS横向视图中初始边界变形的两个实例。

文中所用符号为：I(r，θ)、I′(r，θ)、原极坐标视图和去除环晕伪像后的极坐标视图中像素点(r，θ)处的灰度值；r、像素点的极径；θ、像素点的极角；I_catheter、原IVUS图像中导管区域的像素灰度值；ImageWidth、以像素为单位的图像宽度；ε、权重参数；I(i，j)、坐标为(i，j)的像素的灰度值；I(i，j+1)、坐标为(i，j+1)的像素的灰度值；η、阈值；N、snake点的个数；E、能量函数；E_int、内部能量；E_ext、外部能量；c_i(x_i，y_i)、第i个snake点；(x_i，y_i)、c_i的坐标；c_i-1、c_i+1、第i-1个和第i+1个snake点；

max d、相邻snake点之间的平均距离和最大距离；I(x_i，y_i)、坐标为(x_i，y_i)的像素的灰度和灰度梯度值；α、β、γ、λ、权重参数。

具体实施方式

如附图1所示，本发明方法的步骤包括：

(1)对原始IVUS图像进行预处理，包括滤波去噪和去除环晕伪像：

首先，采用中值滤波和高斯平滑两种通用预处理方法，减少IVUS图像中的椒盐噪声和随机噪声。

然后采用极坐标变换法去除环晕伪像，具体方法如下：

首先，对各帧IVUS图像进行极坐标变换，得到其极坐标视图，如附图2所示。可见，环晕伪像固定的位于极坐标视图的上部。然后，按照下式去除极坐标视图中的环晕伪像：

$I^{\′}(r,\mathrm{\θ})=\left\{\begin{array}{cc}{I}_{\mathrm{catheter}},& r\≤\mathrm{\ϵ}\×\mathrm{ImageWidth}/2\\ I(r,\mathrm{\θ}),& \mathrm{else}\end{array}\right.---\left(1\right)$

其中，I(r，θ)和I′(r，θ)分别为原极坐标视图和去除环晕伪像后的极坐标视图中像素点(r，θ)处的灰度值；r为像素点的极径；θ为像素点的极角；I_catheter是原IVUS图像中导管区域的像素灰度值；ImageWidth为以像素为单位的图像宽度；ε为权重参数，经实验确定其取值区间为[0.1，0.35]。在极坐标视图中消除环晕伪像后，再经过极坐标逆变换，即可得到直角坐标系下去除环晕伪像后的IVUS图像。

(2)获取IVUS图像序列的四个纵向视图：

如附图3所示，分别取IVUS图像序列的沿血管长轴方向的四个纵向视图，即垂直切面A、水平切面B、左对角线切面C和右对角线切面D。

(3)从IVUS纵向视图中提取出血管内腔边界和中-外膜边界：

分别遍历四个IVUS纵向视图中的各像素，判断识别血管壁边界点。由于目标边界分布于纵向视图的中部，所以从纵向视图的中轴线开始，分别向左和向右进行逐行遍历，如附图4所示。用I(i，j)表示坐标为(i，j)的像素的灰度值，I(i，j+1)表示坐标为(i，j+1)的像素的灰度值，若I(i，j+1)-I(i，j)≥η，η为阈值，则(i，j)为目标边界点，否则不是。将每行左右两部分的像素中，第一个符合上述条件的像素记为内腔边界点，第二个记为中-外膜边界点。经实验确定η的取值区间为[10，20]。

(4)将四个纵向视图中的边界曲线映射到各帧IVUS图像中，得到各横向视图中血管壁的初始边界：

如附图5所示，纵向视图的边界线对应到横向视图中为边界点，四条纵向视图边界线映射到横向视图中为四个边界点。依次连接各边界点，即可获得各帧IVUS图像中血管内腔的初始边界和中-外膜的初始边界。

(5)通过使预先设定的能量函数最小，各帧IVUS图像中的初始边界在内外力的共同作用下不断变形，最终得到血管内腔边界和中-外膜边界：

将各帧IVUS图像中血管壁的初始边界作为snake模型(M.Kass，A.Witkin，D.Terzopoulos.Snakes：active contour models.International Journal of Computer Vision.1987，vol.1，no.4，pp.321-331)的初始形状，通过使预先设定的能量函数最小，snake模型不断变形，当能量函数取得全局最小值时，snake模型即停留在目标边界处，从而完成对各帧图像中血管内腔边界和中-外膜边界的并行提取。

能量函数的构成如下：

将初始边界离散成由N个点组成的有序点集，则能量函数的离散表达式为：

$E=\underset{i=0}{\overset{N-1}{\mathrm{\Σ}}}[E_{\mathrm{int}}\left(i\right)+E_{\mathrm{ext}}\left(i\right)]---\left(2\right)$

其中内部能量E_int的归一化表达式为：

$E_{\mathrm{int}}\left(i\right)=\mathrm{\α}\frac{|\stackrel{\‾}{d}-|c_i-c_{i-1}||}{\max d}+\mathrm{\β}\frac{|c_{i-1}-2c_i+c_{i+1}|}{2\max d}---\left(3\right)$

其中c_i(x_i，y_i)(i＝1，2，…，N-2)是第i个snake点，(x_i，y_i)是它的坐标；c_i-1和c_i+1分别是第i-1个和第i+1个snake点；

和max d分别是相邻snake点之间的平均距离和最大距离，在每次迭代结束时，它们的值都被更新；α和β是权重参数。(3)式中的第一项保证变形过程中snake的连续性，使snake点均匀分布，不致产生收缩的现象；第二项是曲线二阶导数的离散形式，保证变形过程中snake的光滑性。由于进行了归一化，故(3)式中两项的取值范围都在[0，1]区间内，权重α和β的取值区间也是[0，1]。

外部能量E_ext定义为：

$E_{\mathrm{ext}}\left(i\right)=\mathrm{\γ}\frac{I(x_i,y_i)}{255}-\mathrm{\λ}\frac{|\▿I(x_i,y_i)|}{255\sqrt{2}}---\left(4\right)$

其中I(x_i，y_i)和

分别是坐标为(x_i，y_i)的像素的灰度和灰度梯度值。由于8位灰度图像的灰度和灰度梯度的取值范围分别为[0，255]和

故本发明对I(x_i，y_i)和

分别进行了归一化，使其取值范围为[0，1]区间。γ，λ∈[0，1]是权重参数。

通过使能量函数E最小，snake模型从初始形状开始不断变形，最终停留在能量函数取得全局最小值的最优位置，即为目标边界。能量函数的全局最优化采用Williams贪婪算法(D.J.Williams，M.Shah.A fast algorithm for active contours and curvature estimation.Computer Vision，Graphics and Image Processing.1992，vol.55，no.1，pp.14-26)完成，不仅计算速度快，所需存贮空间小，而且非常适合于编程实现。

对临床采集的冠状动脉IVUS图像序列的实验证明，在获得与现有逐帧处理的二维串行分割方法相类似的分割精度的情况下，采用本发明的三维分割方法可大大提高处理效率，缩短处理时间。附图6～9是采用本发明方法对一个冠状动脉IVUS图像序列的分割结果。

Claims (4)

1.一种血管内超声图像序列的三维分割方法，其特征是，首先对IVUS原始图像进行滤除噪声和抑制环晕伪像的预处理，然后获取IVUS图像序列的四个纵向视图，并从中提取出血管内腔边界和中-外膜边界，再通过将这些边界曲线映射到各帧IVUS图像中，得到横向视图中的初始边界，最后通过使能量函数最小，初始边界不断变形，最终得到各帧IVUS图像中的血管内腔和中-外膜边界。

2.根据权利要求1所述血管内超声图像序列的三维分割方法，其特征是，该方法包括以下步骤：

a、对原始IVUS图像进行预处理，包括滤波去噪和去除环晕伪像：

首先采用中值滤波和高斯平滑两种通用预处理方法，减少IVUS图像中的椒盐噪声和随机噪声，然后对各帧IVUS图像进行极坐标变换，得到其极坐标视图，再按照下式去除极坐标视图中的环晕伪像：

$I^{\′}(r,\mathrm{\θ})=\left\{\begin{array}{cc}{I}_{\mathrm{catheter}},& r\≤\mathrm{\ϵ}\×\mathrm{ImageWidth}/2\\ I(r,\mathrm{\θ}),& \mathrm{else}\end{array}\right.$

其中，I(r，θ)和I′(r，θ)分别为原极坐标视图和去除环晕伪像后的极坐标视图中像素点(r，θ)处的灰度值；r为像素点的极径；θ为像素点的极角；I_catheter是原IVUS图像中导管区域的像素灰度值；ImageWidth为以像素为单位的图像宽度；ε为权重参数；最后经过极坐标逆变换，即可得到直角坐标系下去除环晕伪像后的IVUS图像；

b、分别取IVUS图像序列的沿血管长轴方向的四个纵向视图，即垂直切面A、水平切面B、左对角线切面C和右对角线切面D；

c、从IVUS纵向视图中提取出血管内腔边界和中-外膜边界：

从纵向视图的中轴线开始，分别向左和向右逐行遍历四个IVUS纵向视图中的各像素，用I(i，j)表示坐标为(i，j)的像素的灰度值，I(i，j+1)表示坐标为(i，j+1)的像素的灰度值，若I(i，j+1)-I(i，j)≥η，则(i，j)为目标边界点，否则不是；将每行左右两部分的像素中，第一个符合上述条件的像素记为内腔边界点，第二个记为中-外膜边界点；其中，η为阈值；

d、将四个纵向视图中的边界曲线映射到各帧IVUS图像中，纵向视图的边界线对应到横向视图中为边界点，依次连接各边界点，得到各帧IVUS图像中血管内腔的初始边界和中-外膜的初始边界；

e、各帧IVUS图像中的初始边界通过snake变形获得血管内腔边界和中-外膜边界：

将各帧IVUS图像中血管壁的初始边界作为snake模型的初始形状，将初始边界离散成由N个点组成的有序点集，则能量函数的离散表达式为：

$E=\underset{i=0}{\overset{N-1}{\mathrm{\Σ}}}[E_{\mathrm{int}}\left(i\right)+E_{\mathrm{ext}}\left(i\right)]$

$E_{\mathrm{int}}\left(i\right)=\mathrm{\α}\frac{|\stackrel{\‾}{d}-|c_i-c_{i-1}||}{\max d}+\mathrm{\β}\frac{|c_{i-1}-2c_i+c_{i+1}|}{2\max d},$ $E_{\mathrm{ext}}\left(i\right)=\mathrm{\γ}\frac{I(x_i,y_i)}{255}-\mathrm{\λ}\frac{|\▿I(x_i,y_i)|}{255\sqrt{2}}$

其中，E_int是内部能量；E_ext是外部能量；c_i(x_i，y_i)(i＝1，2，…，N-2)是第i个snake点，(x_i，y_i)是它的坐标；c_i-1和c_i+1分别是第i-1个和第i+1个snake点；和maxd分别是相邻snake点之间的平均距离和最大距离；I(x_i，y_i)和

分别是坐标为(x_i，y_i)的像素的灰度和灰度梯度值；α，β，γ，λ∈[0，1]都是权重参数。

通过使能量函数E最小，snake模型从初始形状开始不断变形，最终停留在能量函数取得全局最小值的最优位置，即为目标边界。

3.根据权利要求2所述血管内超声图像序列的三维分割方法，其特征是，所述权重参数ε的取值区间为[0.1，0.35]，阈值η的取值区间为[10，20]。

4.根据权利要求3所述血管内超声图像序列的三维分割方法，其特征是，能量函数的全局最优化采用Williams贪婪算法完成。

Images