基于图像特征的CT-MR图像配准
作者：赵富强 于红玉 马斌荣 李坤成
单位：赵富强 于红玉 马斌荣(首都医科大学 北京 100054)；李坤成(首都医科大学附属宣武医院 北京 100053)
关键词：配准；图像特征；相关性；二阶导数
北京生物医学工程000404 摘 要 本文采用了基于图像特征的相关性算法对3-D CT-MR图像进行配准。首先使用基于局部坐标系的二阶导数算子对CT和MR图像进行处理，获取由图像的“脊”或“谷”组成的特征图像，然后在进行刚体变换条件下对特征图像进行相关计算，使其相关系数达到最大，从而获取配准后的图像。该方法不需在成像过程中人为增加标记物，后续处理中不需人为干预，是一种自动的图像配准方法。最后本文给出了这一方法用于实际数据的结果。
CT-MRI Registration Based on Images' Geometrical Features
, http://www.100md.com
Zhao Fuqiang，Yu Hongyu，Ma Binrong，Li Kuncheng
(Capital University of Medical Science(100054)
Abstract In this paper is described an automated approach to register CT and MRI brain images based on their geometrical features. Local second order differential operators in scale space are applied to each type image date, so as to produce feature images, depicting ridges or vallegs. Then correlation calculation is carried out for feature images under the condition of rigidity, so that the coeffiuent of correlation is maximum. The resulted CT with MRI feature images are then obtained. This method does not need human interaction and interactive operation is also not needed during post-processing. Results obtained using this method are presented. It is proved that this method is practically useful.
, 百拇医药
Key words:Registration, Images' geometrical feature, Correlation, Second differential operator
0 引 言
医学图像配准是指对一3-D医学图像寻求一种(或一系列)空间变换，使它与另一3-D医学图像上的对应点达到空间上的一致。这种一致是指位于人体上的同一解剖点在两个匹配图像上的空间位置相同。
医学图像配准具有很重要的临床应用价值。对各种使用不同或相同的成像手段所获的医学图像进行配准不仅可以用于医疗诊断，还可用于手术计划的制定、放射治疗计划的制定、病理变化的跟踪和治疗效果的评价等各个方面。例如，在制定放射治疗计划时，需要用X-CT图像进行放射剂量分布的计算并决定射线的投照区域，而病灶区域的轮廊通常在MRI(磁共振成像)中能很好地体现出来，在实际中用经过配准的MRI图像确定病灶区域后把这一区域转换为相应的CT图像区域，从而在具有CT跟踪条件下的放射治疗系统很好地完成治疗过程。在神经外科手术中，配准后的CT/MRI图像能为手术医生提供更多的病变区域及周围相关的解剖结构信息[1，2]。
, http://www.100md.com
本文首先分析了基于二维图像场的梯度场在图像的“脊”或“谷”上的表现特征，然后导出了对图像场沿与梯度向量垂直方向的二维导数的计算公式，用这一计算公式对图像进行处理就获取了由图像的“脊”或“谷”组成的特征图像。对特征图像在一定变换条件下进行相关计算，在相关系数最大的情况下就获取了配准后的结果图像。为了运行速度的考虑，本文同时给出了一种多分辨率分解技术，既首先在低图像分辨率条件下选择大致的配准参数，然后分辨率逐渐升高，而获取最终准确的配准参数。最后本文给出了对临床获取的实际数据的配准结果。
1 原理和方法
1.1 图像特征—“脊”和“谷”的检测方法
在成像过程中，成像体的一些特征，比如骨组织，根据成像原理的不同，将在结果图像中表现为山脊一样的特征，我们简称其为“脊”，或表现为山谷一样的特征，我们简称其为“谷”。在XCT成像中，由于骨组织对X线有较大的吸收系数，因此在结果图像中，骨组织表现为亮区，也就是“脊”；而在MR成像中，由于骨组织含有较低的质子密度，因此在结果图像中表现为暗的区域，也就是“谷”的特征。若能将CT图像和MR图像中的脊和谷检测出来，使这些特征能够有效的配准，然后将与之对应的配准参数应用于整幅图像，就可得到配准后的三维图像。
, 百拇医药
检测脊和谷的方法有很多，这里我们推导并选择与图像的等高线特征密切相关的局部特征拉普拉斯特征作为检测方法。我们知道，一幅医学图像可以被认为是在二维空间分布的一个标量场。在这一标量场中，由各点处与等高线垂直的并指向增大方向的矢量组成的场为梯度场，各点处的矢量值为标量函数f(x,y)在该点处的梯度，即
(1)
图1 图像特征—脊的提取
左侧是一个二维图像的强度立体图。B为脊上的一点，A是非脊点。右边的矩形框里分别是点A、点B处的w、v的强度图。
i→和j→分别为指向x方向和y方向的单位矢量。下面我们用图示的方法来研究这种梯度矢量在脊或谷上的表现。
, http://www.100md.com
从图1中可以看出，在脊上，梯度向量w的方向与脊的方向一致，而与梯度向量w垂直的矢量v则指向图像灰度变化率最大的方向。我们定义与w垂直的矢量为v，则沿v方向的对图像的二阶导数(拉普拉斯算子)即必将在脊处或谷处取最大值。设w与x轴的夹角为α。w与y轴的夹角为β，由公式(1)得
(2)
(3)
由于向量v与向量w垂直，向量v与x轴的夹角为θ，向量v与y轴的夹角为φ，
(4)
(5)
我们若定义在点p处由向量w和向量v组成的坐标系为一新的平面直角坐标系。根据坐标变换理论，图像平面中的任一点p(x,y)在新坐标系中的坐标为
, http://www.100md.com
x=wcosα+vcosθ+xp (6)
y=wcosβ+vcosφ+yp (7)
根据求导法则，f在v方向的二阶导数2f/v2为
(8)
(9)
而
(10)
, http://www.100md.com
由于
(11)
(12)
在CT图像中，脊处的2f/v2值通常为较大的负值，在MR图像中，谷处的2f/v2值为较大的正值。为了更好的提取脊或谷，避免图像中高频分量及噪声对二阶导数的严重影响，先对原始图像进行平滑，这里我们选用高斯函数与图像进行卷积而达到平滑的目的。
1.2 相关性
相关性的计算采用下面的公式：
(13)
, 百拇医药
其中，f、g为两待配准特征图像。X为图像f上的一点(x,y,z)。f(X)为点X的灰度值。T(X)表示对X进行刚体变换，T(x)=RX+t=Rx*Ry*Rz*X+t(tx,ty,tz)。R是3*3的旋转矩阵，t是3*1的平移矩阵。(x,y,z)为整数坐标，但T(X)不一定为整数坐标，因此需要采用插值技术来得到点T(X)的灰度值。插值方法为三线性插值。整个配准过程就是寻找一刚体变换T′，使C具有最大值C*，这时即认为已经配准。
2 结果与讨论
由于原始图像中含有很多高频噪声，不利于“脊”或“谷”的检测，为了更好地提取“脊”或“谷”，需要对图像进行平滑，方法为高斯滤波。图像及平滑后的三维立体图如图2所示。
图2 (a)CT原始图像 (b)CT图像的三维显示 (c)CT图像经高斯滤波(σ=4)后的三维显示
, http://www.100md.com
图像特征—“脊”或“谷”的检测结果如图3所示。
从图3可以看出，CT图像的特征提取优于MR图像。如果在MR图像中的某一区域为谷，而在CT图像相应的区域为非脊点，一个较大的值与一个接近于零的值相乘，结果仍近似为零，因此，几乎不会对C的值产生影响。另外，这种算法对数据缺失不敏感，只要特征图像中有足够的相似结构，就可以实现配准。为了简化计算，我们把CT特征图像取绝对值，再将负值置为0，在MR特征图像中把负值置为0。
图3 图像特征的检测结果
(a)CT图像 (b)MR图像 (c)CT图像的“脊”特征图像 (d)MR图像的“谷”特征图像
对于待配准的图像定义一个统一的立体坐标系统：X轴沿着图像的行扫描方向，Y轴沿着图像的列扫描方向，Z轴为从颅顶到颅底的方向。对一组CT-MR图像进行配准，配准变换参数为T(tx,ty,tz,rx,ry,rz)=(-2.5，-3.0，0.0，0.0，0.0，-6.0°)。其中tx、ty、tz分别为沿x、y、z方向的平移，rx、ry、rz分别为绕x、y、z方向的旋转。配准结果显示如图4所示。
, http://www.100md.com
配准结果表明，采用基于图像特征的相关算法是可行的。基于图像特征的相关性配准算法的特点如下：
1、算法简单，不需要人为干涉。
2、对数据缺失不敏感。
将配准后的图像在计算机屏幕上叠加显示，可以使医生更清楚地了解病变组织在体内的解剖位置，病变组织在对应的CT/MR图像上的表现为临床诊断和治疗提供更多的信息。比如对于颅内肿瘤病人，肿瘤组织是否钙化及瘤体内基质含量对治疗和预后都有指导意义，通过对配准后的图像进行观察，因为CT对钙化敏感，而MR对软组织基质敏感，因此能更准确地获取这些信息。我们课题组正在对这种方法对于临床不同疾病的诊断和治疗的指导意义进行实验。
图4 CT与MR的配准结果，从左到右分别为第2层、第3层、第5层、第7层
, 百拇医药
(a)原始CT图像 (b)MR图像 (c)经过配准变换后的CT图像
作者简介：赵富强(1965-)男，生物医学工程博士，副教授，现为首都医科大学工程系电子学教研室主任。
参考文献
[1] Petra A. den Elsen, et al. Medical image matching-a review with classification. IEEE Engineering in Medical and Biology,1993,3:26-39
[2] Petra A. Van den Elsen, et al. Automatic registration of CT and MR brain images using correlation of geometrical features. IEEE trans Medical Imaging,1995,7:384-396
, 百拇医药
[3] Beolchi L, Kuhn M H. Medical imaging-analysis of multimodality 2D/3D images. IOS Press, 1995,19.
[4] Marcel van Herk, Hanne M. Kooy.Automatic three-dimensional correlation of CT-CT,CT-MRI, and CT-SPECT using chamfer matching. Med Phys,1994,7:1163-1179
[5] Antoine J B. Maintz, Petra A. van den Elsen, et al. Evalution of ridge seeking operators for multimodality medical image matching. IEEE trans On Pattern Anal and INTELL. 1996,4:353-363
[6] Charles A. Pelizzari, George T Y.Chen, et al. Accurate three-dimensional registration of CT, PET, and/or MR images of the brain. Journal of Computer Assisted Tomography,1989,13(1):20-26
(1999-09-21收稿，1999-10-25修回), 百拇医药
单位：赵富强 于红玉 马斌荣(首都医科大学 北京 100054)；李坤成(首都医科大学附属宣武医院 北京 100053)
关键词：配准；图像特征；相关性；二阶导数
北京生物医学工程000404 摘 要 本文采用了基于图像特征的相关性算法对3-D CT-MR图像进行配准。首先使用基于局部坐标系的二阶导数算子对CT和MR图像进行处理，获取由图像的“脊”或“谷”组成的特征图像，然后在进行刚体变换条件下对特征图像进行相关计算，使其相关系数达到最大，从而获取配准后的图像。该方法不需在成像过程中人为增加标记物，后续处理中不需人为干预，是一种自动的图像配准方法。最后本文给出了这一方法用于实际数据的结果。
CT-MRI Registration Based on Images' Geometrical Features
, http://www.100md.com
Zhao Fuqiang，Yu Hongyu，Ma Binrong，Li Kuncheng
(Capital University of Medical Science(100054)
Abstract In this paper is described an automated approach to register CT and MRI brain images based on their geometrical features. Local second order differential operators in scale space are applied to each type image date, so as to produce feature images, depicting ridges or vallegs. Then correlation calculation is carried out for feature images under the condition of rigidity, so that the coeffiuent of correlation is maximum. The resulted CT with MRI feature images are then obtained. This method does not need human interaction and interactive operation is also not needed during post-processing. Results obtained using this method are presented. It is proved that this method is practically useful.
, 百拇医药
Key words:Registration, Images' geometrical feature, Correlation, Second differential operator
0 引 言
医学图像配准是指对一3-D医学图像寻求一种(或一系列)空间变换，使它与另一3-D医学图像上的对应点达到空间上的一致。这种一致是指位于人体上的同一解剖点在两个匹配图像上的空间位置相同。
医学图像配准具有很重要的临床应用价值。对各种使用不同或相同的成像手段所获的医学图像进行配准不仅可以用于医疗诊断，还可用于手术计划的制定、放射治疗计划的制定、病理变化的跟踪和治疗效果的评价等各个方面。例如，在制定放射治疗计划时，需要用X-CT图像进行放射剂量分布的计算并决定射线的投照区域，而病灶区域的轮廊通常在MRI(磁共振成像)中能很好地体现出来，在实际中用经过配准的MRI图像确定病灶区域后把这一区域转换为相应的CT图像区域，从而在具有CT跟踪条件下的放射治疗系统很好地完成治疗过程。在神经外科手术中，配准后的CT/MRI图像能为手术医生提供更多的病变区域及周围相关的解剖结构信息[1，2]。
, http://www.100md.com
本文首先分析了基于二维图像场的梯度场在图像的“脊”或“谷”上的表现特征，然后导出了对图像场沿与梯度向量垂直方向的二维导数的计算公式，用这一计算公式对图像进行处理就获取了由图像的“脊”或“谷”组成的特征图像。对特征图像在一定变换条件下进行相关计算，在相关系数最大的情况下就获取了配准后的结果图像。为了运行速度的考虑，本文同时给出了一种多分辨率分解技术，既首先在低图像分辨率条件下选择大致的配准参数，然后分辨率逐渐升高，而获取最终准确的配准参数。最后本文给出了对临床获取的实际数据的配准结果。
1 原理和方法
1.1 图像特征—“脊”和“谷”的检测方法
在成像过程中，成像体的一些特征，比如骨组织，根据成像原理的不同，将在结果图像中表现为山脊一样的特征，我们简称其为“脊”，或表现为山谷一样的特征，我们简称其为“谷”。在XCT成像中，由于骨组织对X线有较大的吸收系数，因此在结果图像中，骨组织表现为亮区，也就是“脊”；而在MR成像中，由于骨组织含有较低的质子密度，因此在结果图像中表现为暗的区域，也就是“谷”的特征。若能将CT图像和MR图像中的脊和谷检测出来，使这些特征能够有效的配准，然后将与之对应的配准参数应用于整幅图像，就可得到配准后的三维图像。
, 百拇医药
检测脊和谷的方法有很多，这里我们推导并选择与图像的等高线特征密切相关的局部特征拉普拉斯特征作为检测方法。我们知道，一幅医学图像可以被认为是在二维空间分布的一个标量场。在这一标量场中，由各点处与等高线垂直的并指向增大方向的矢量组成的场为梯度场，各点处的矢量值为标量函数f(x,y)在该点处的梯度，即
(1)图1 图像特征—脊的提取
左侧是一个二维图像的强度立体图。B为脊上的一点，A是非脊点。右边的矩形框里分别是点A、点B处的w、v的强度图。
i→和j→分别为指向x方向和y方向的单位矢量。下面我们用图示的方法来研究这种梯度矢量在脊或谷上的表现。
, http://www.100md.com
从图1中可以看出，在脊上，梯度向量w的方向与脊的方向一致，而与梯度向量w垂直的矢量v则指向图像灰度变化率最大的方向。我们定义与w垂直的矢量为v，则沿v方向的对图像的二阶导数(拉普拉斯算子)即必将在脊处或谷处取最大值。设w与x轴的夹角为α。w与y轴的夹角为β，由公式(1)得
(2) (3)由于向量v与向量w垂直，向量v与x轴的夹角为θ，向量v与y轴的夹角为φ，
(4) (5)我们若定义在点p处由向量w和向量v组成的坐标系为一新的平面直角坐标系。根据坐标变换理论，图像平面中的任一点p(x,y)在新坐标系中的坐标为
, http://www.100md.com
x=wcosα+vcosθ+xp (6)
y=wcosβ+vcosφ+yp (7)
根据求导法则，f在v方向的二阶导数2f/v2为
(8) (9)而
(10), http://www.100md.com
由于
(11) (12)在CT图像中，脊处的2f/v2值通常为较大的负值，在MR图像中，谷处的2f/v2值为较大的正值。为了更好的提取脊或谷，避免图像中高频分量及噪声对二阶导数的严重影响，先对原始图像进行平滑，这里我们选用高斯函数与图像进行卷积而达到平滑的目的。
1.2 相关性
相关性的计算采用下面的公式：
(13), 百拇医药
其中，f、g为两待配准特征图像。X为图像f上的一点(x,y,z)。f(X)为点X的灰度值。T(X)表示对X进行刚体变换，T(x)=RX+t=Rx*Ry*Rz*X+t(tx,ty,tz)。R是3*3的旋转矩阵，t是3*1的平移矩阵。(x,y,z)为整数坐标，但T(X)不一定为整数坐标，因此需要采用插值技术来得到点T(X)的灰度值。插值方法为三线性插值。整个配准过程就是寻找一刚体变换T′，使C具有最大值C*，这时即认为已经配准。
2 结果与讨论
由于原始图像中含有很多高频噪声，不利于“脊”或“谷”的检测，为了更好地提取“脊”或“谷”，需要对图像进行平滑，方法为高斯滤波。图像及平滑后的三维立体图如图2所示。
图2 (a)CT原始图像 (b)CT图像的三维显示 (c)CT图像经高斯滤波(σ=4)后的三维显示
, http://www.100md.com
图像特征—“脊”或“谷”的检测结果如图3所示。
从图3可以看出，CT图像的特征提取优于MR图像。如果在MR图像中的某一区域为谷，而在CT图像相应的区域为非脊点，一个较大的值与一个接近于零的值相乘，结果仍近似为零，因此，几乎不会对C的值产生影响。另外，这种算法对数据缺失不敏感，只要特征图像中有足够的相似结构，就可以实现配准。为了简化计算，我们把CT特征图像取绝对值，再将负值置为0，在MR特征图像中把负值置为0。
图3 图像特征的检测结果
(a)CT图像 (b)MR图像 (c)CT图像的“脊”特征图像 (d)MR图像的“谷”特征图像
对于待配准的图像定义一个统一的立体坐标系统：X轴沿着图像的行扫描方向，Y轴沿着图像的列扫描方向，Z轴为从颅顶到颅底的方向。对一组CT-MR图像进行配准，配准变换参数为T(tx,ty,tz,rx,ry,rz)=(-2.5，-3.0，0.0，0.0，0.0，-6.0°)。其中tx、ty、tz分别为沿x、y、z方向的平移，rx、ry、rz分别为绕x、y、z方向的旋转。配准结果显示如图4所示。
, http://www.100md.com
配准结果表明，采用基于图像特征的相关算法是可行的。基于图像特征的相关性配准算法的特点如下：
1、算法简单，不需要人为干涉。
2、对数据缺失不敏感。
将配准后的图像在计算机屏幕上叠加显示，可以使医生更清楚地了解病变组织在体内的解剖位置，病变组织在对应的CT/MR图像上的表现为临床诊断和治疗提供更多的信息。比如对于颅内肿瘤病人，肿瘤组织是否钙化及瘤体内基质含量对治疗和预后都有指导意义，通过对配准后的图像进行观察，因为CT对钙化敏感，而MR对软组织基质敏感，因此能更准确地获取这些信息。我们课题组正在对这种方法对于临床不同疾病的诊断和治疗的指导意义进行实验。
图4 CT与MR的配准结果，从左到右分别为第2层、第3层、第5层、第7层
, 百拇医药
(a)原始CT图像 (b)MR图像 (c)经过配准变换后的CT图像
作者简介：赵富强(1965-)男，生物医学工程博士，副教授，现为首都医科大学工程系电子学教研室主任。
参考文献
[1] Petra A. den Elsen, et al. Medical image matching-a review with classification. IEEE Engineering in Medical and Biology,1993,3:26-39
[2] Petra A. Van den Elsen, et al. Automatic registration of CT and MR brain images using correlation of geometrical features. IEEE trans Medical Imaging,1995,7:384-396
, 百拇医药
[3] Beolchi L, Kuhn M H. Medical imaging-analysis of multimodality 2D/3D images. IOS Press, 1995,19.
[4] Marcel van Herk, Hanne M. Kooy.Automatic three-dimensional correlation of CT-CT,CT-MRI, and CT-SPECT using chamfer matching. Med Phys,1994,7:1163-1179
[5] Antoine J B. Maintz, Petra A. van den Elsen, et al. Evalution of ridge seeking operators for multimodality medical image matching. IEEE trans On Pattern Anal and INTELL. 1996,4:353-363
[6] Charles A. Pelizzari, George T Y.Chen, et al. Accurate three-dimensional registration of CT, PET, and/or MR images of the brain. Journal of Computer Assisted Tomography,1989,13(1):20-26
(1999-09-21收稿，1999-10-25修回), 百拇医药