CN110223331B - 一种大脑mr医学图像配准方法 - Google Patents

Abstract

本发明公开一种大脑MR医学图像配准方法，包括1采用BCFCM方法对参考图像和浮动图像进行分割；2采用MSR检测方法对原始参考图像和原始浮动图像进行对称轴检测并提取对称轴方程；3对步骤1处理后图像进行阈值分割；4根据步骤2所检测的对称轴方程对步骤3处理后图像进行近似对称性约束，所得图像记为待配准的参考图像和浮动图像；5初始化刚体变换矩阵配准参数；6通过SSD相似性度量准则配准图像，得最优变换矩阵配准参数。以BCFCM方法对参考图像、浮动图像分割；分割后的图像通过二值化阈值处理分离前背景，对二值化后处理的图像进行对称性约束并用于SSD框架下的多模态图像配准。使其配准效率、精度及鲁棒性都提升。

Description

一种大脑MR医学图像配准方法

技术领域

本发明涉及一种大脑MR医学图像配准方法。

背景技术

图像配准是将不同获取时间、不同传感器、不同获取条件的同一场景或者同一目标的两幅或者多幅图像进行匹配的过程。目前，图像配准技术已经被广泛应用于遥感图像处理、医学图像处理等领域。

对于多模态脑部核磁共振MR图像配准，需要将浮动脑部MR图像与参考脑部MR图像进行匹配(浮动图像与参考图像模态不同)，由于MR各模态图像像素信息的不一致性以及各模态图像自身存在灰度不均匀和受外界噪声影响等局限性，如何将大脑MR浮动图像精准配准于参考图像成为一个有待解决的问题。

现有技术主要包括像素级配准及特征级配准两大类。其中，像素级配准技术方案，包括主轴质心法、互相关法、最大互信息法等；特征级配准技术方案，主要有基于目标边缘、纹理、曲率等图像信息的方法。

一种基于互信息与图像分割的三维多模态医学图像自动配准方法，具体步骤：步骤1：采用门限法和数学形态学方法进行预处理；步骤2：采用k-means方法进行分割；步骤3：采用优化算法通过迭代得到基于互信息的最优配准参数；步骤4：将原始的参考图像和浮动图像相叠加；步骤5：计算经过图像分割预处理的参考图像A的灰度直方图，并把具有相同灰度值的像素划分为一组；步骤6：初始化配准参数，设定六个参数初始值为零；步骤7：利用配准参数对浮动图像B进行线性插值产生变化后的浮动图像，对在迭代过程中映射到参考图像之外的浮动图像上的像素点赋零值。

由于医学图像存在灰度不均匀和仪器噪声等问题，采用k-means方法对医学图像进行分割的效果较差。对于图像中灰度值与背景相似的组织器官，易漏检；对于图像中灰度值突变的噪声，易错检。

此外，由于互信息测度需要统计配准图像及待配准图像像素值大小并计算灰度直方图，并根据灰度直方图计算互信息测度值，进一步进行图像配准，故该算法测度函数较为复杂，配准时间长，效率低。

发明内容

发明目的：为了克服现有技术的不足，本发明提供一种可有效提高配准精度的大脑MR医学图像配准方法。

技术方案：一种大脑MR医学图像配准方法，所述配准方法包括以下步骤

步骤1：采用BCFCM方法对参考图像和浮动图像进行分割；

步骤2：采用MSR目标对称轴检测方法对原始参考图像和原始浮动图像进行对称轴检测并提取对称轴方程；

步骤3：对步骤1处理后图像进行二值化阈值分割；

步骤4：根据步骤2所检测的对称轴方程对步骤3处理后图像进行近似对称性约束，所得到的图像分别记为待配准的参考图像和浮动图像；

步骤5：初始化刚体变换矩阵配准参数，设定参数初始值为0；

步骤6：通过SSD相似性度量准则配准图像，得到最优变换矩阵配准参数。

所述BCFCM方法对参考图像和浮动图像进行分割，其目标函数为：

其中，y_i为观测到的图像像素，表示该像素受偏离场影响；β_i代表偏离场，x_i＝y_i-β_i为去偏离场后的图像像素；α为比例因子，控制像素邻域对图像分割的影响；

表示以第i个像素点为邻域中心的邻域方形窗像素点集合，不包括i本身，其中R为邻域方形窗的边长；y_r及β_r分别表示领域像素值及该像素值的偏离场；N表示像素总个数；c_j表示聚类簇心，C表示待划分像素的类别总数，待划分像素类别总数C≥2；其中一个聚类簇心初始值设为0，记为c_j＝0＝0；m为一个大于0的任意实数；u_ij为隶属度，表示第i个像素点属于第j类的概率，满足

在图像分割中，该值组成的矩阵U称为分割矩阵；利用(1)式进行优化得到最优隶属度

并组合成最优分割矩阵U^*。

选择U^*中类别为背景灰度值的分割子矩阵

作为原始图像的分割矩阵来分割图像，其格式为:

设分割图像的像素矩阵大小为x×y，则：

其中，N表示像素总数；

为最优隶属度，表示第i个像素点属于第j类的最优概率。

步骤2所述的MSR对称轴检测方法流程如下：

(1)将过图像中心点的垂线置为初始对称轴；

(2)将原始图像关于该对称轴进行镜像变换，得到镜面图像；

(3)将原始图像作为参考图像，镜面图像作为浮动图像进行刚性配准，得到刚性变换矩阵参数；

(4)根据刚性变换矩阵参数得到原图像对称轴。

将步骤1处理后的图像进行阈值分割，经过二值化阈值分割处理后的图像分割矩阵记为P，格式如下所示：

其中p(g)为阈值分割函数，δ为阈值，在理想情况下δ＝1。

根据步骤2所检测的对称轴方程，对步骤3处理后图像进行近似对称性约束，所得到的图像分别记为待配准的参考图像和浮动图像，对称轴方程记为l，记矩阵P经对称性约束后的矩阵为S，公式如下所示：

实际使用中令s(h)＝s_sparse(h)或者s(h)＝s_dense(h)，其中s_sparse(h)和s_dense(h)分别为稀疏性和稠密性对称性约束函数，其定义式为：

h′为关于图像对称轴l的对称点，h′_a为对称点邻域像素值(包括h′)，A²为对称点邻域方形窗大小，表示对称点邻域方形窗内像素点个数。

步骤5所述初始化刚体变换矩阵配准参数，设定参数初始值为0；所述刚体变换的过程为：

将具有6个自由度的刚体变换作为人脑多模态MR医学图像配准的几何变换参数矩阵，所述6个自由度的刚体包括t_x、t_y、t_z三个位移自由度，θ_x、θ_y、θ_z三个旋转自由度。

步骤6，通过SSD相似性度量准则配准图像，得到最优变换矩阵配准参数，其中

SSD相似性度量优化公式为：

其中，T表示由6个自由度组成的刚体变换参数矩阵；M、N代表图像像素位置，即M×N代表图像像素总个数；I和I′分别代表待配准浮动图像和待配准参考图像像素值大小矩阵。

配准过程的优化目标为：

在对配准浮动图像进行刚性变换(计算T(i,j))时，由于其值并不一定为整数坐标，故通过双线性插值法，使其映射到整数坐标；该配准通过高斯牛顿法求解(9)式最优值获得最佳刚体变换矩阵参数。

有益效果：本发明的一种大脑MR医学图像配准方法，具有以下优点：

1.本发明通过采用偏离场矫正模糊C均值聚类分割方法对参考图像和浮动图像进行分割；可将偏离场对图像影响加入考虑范围并在优化算法模型中加入可调的邻域约束项，以此减少突变噪声对分割效果的影响。

2.本发明对经BCFCM分割后的图像通过二值化阈值处理分离前背景，将原始图像经基于配准的镜像对称性检测方法处理提取对称轴，对二值化后处理的图像进行对称性约束并最终用于SSD框架下的多模态图像配准。改进的方法在配准效率、配准精度及鲁棒性上都大大提高。

附图说明

图1为本发明的流程图。

图2为使用s_sparse(h)约束函数的效果示意图；

图3使用s_dense(h)约束函数的效果示意图。

具体实施方式

下面结合附图对本发明的技术方案作进一步说明。

本发明的图像配准方法关注的是基于图像分割及对称性检测下的图像配准。为了克服MR多模态医学图像存在偏离场及噪声等干扰因素对图像配准性能的影响，本发明先后进行了如图1所示的操作来处理图像，并采用图像配准指标目标配准误差(TargetRegistration Error,TRE)在数据集RIRE和BrainWeb上对本发明的方法效果作出评价。

目标配准误差(TRE)为配准后图像和物理目标物位置之间的距离，其格式如下式表示：

其中，D表示已配准图像的像素点，

表示配准后图像中像素点坐标位置，

表示与

对应的目标物理坐标位置，N表示图像维度大小。为了便于研究分析，对于仅涉及平移旋转变换的配准来说，上式可转化为如下形式：

TRE(D)＝||T₀-T_D|| (12)

其中，T₀为待配准图像转化为目标参考图像的标准变换矩阵，T_D为待配准图像转化为经算法配准后获得的已配准图像的变换矩阵。实际仿真实验中，先将已配准的两幅图像分别作为浮动图像和参考图像，然后对参考图像作刚体变换(平移旋转变换)，该变换矩阵即为T₀。

为了测试配准图像该如何选择对称性约束条件，以及不同对称性约束条件对配准效果的影响，本发明对比了稀疏性约束条件和稠密性约束条件对配准性能的影响。选用图像大小为256×256像素的RIRE数据集MR大脑切片T1、T2、PD模态图像分别作为浮动图像和参考图像进行实验，其中PD以及T2模态的图像经所提方法处理流程如图2和图3所示；图2为使用s_sparse(h)约束函数的效果示意图，图2(a)为RIRE数据集真实病人核磁共振PD模态大脑图像的原始图像；图2(b)为经过MSR检测后的图像，绿线为对称轴；图2(c)为经过BCFCM图像分割后，根据背景类分割后的图像

图2(d)为根据阈值处理后的二值化图像P，其中δ置为0.8；图2(e)为经过稀疏对称性约束后的图像S，A置为1。

图3为使用s_dense(h)约束函数的效果示意图，图3(a)为RIRE数据集真实病人核磁共振T2模态大脑图像的原始图像；图3(b)为经过MSR检测后的图像，绿线为对称轴；图3(c)为经过BCFCM图像分割后，根据背景类分割后图像

图3(d)为根据阈值处理后的二值化图像P，其中δ置为0.8；图3(e)为经过稠密对称性约束后的图像S，A置为1。

实验分为两部分：第一部分仅考虑图像绕图像中心旋转所带来的多模态大脑MR图像配准的影响，其中旋转角度范围为[-20°,20°]；

第二部分仅考虑图像位移所带来的多模态大脑MR图像配准的影响，其中水平及垂直位移范围为[-20,20]像素。

实验对比结果如表1所示。其中”-”前为浮动图像，”-”后为参考图像。譬如，表中“稀疏-稠密”表示浮动图像选用稀疏性约束，参考图像选用稠密性约束。T2模态的分割图像受噪声及偏移场的影响过大，故使用稀疏性约束条件后的图像效果较差，但T1及PD模态图像均不存在此问题，故稀疏-稠密条件下的对比实验仅考虑T1及PD稀疏约束，T2稠密约束。

表1对称性约束对比实验

由表1可看出，浮动图像和参考图像均使用稀疏约束时，仅T1-PD可配准，但无论是平移还是旋转，配准性能最好，原因在于此数据集T2模态图像使用稀疏性约束后的分割图像缺失严重，达不到配准要求，但T1、PD模态下图像不存在此问题；此外，由于使用稠密性约束可能会造成图像噪声的扩大，导致图像误分割严重，该问题对配准造成了一定的负面影响。当T2模态使用稠密性约束，TI或者PD模态使用稀疏性约束配准时，其效果并没有配准图像和参考图像都使用稠密性约束配准效果好。故经综合考虑，如所配准图像的参考图像及浮动图像的灰度不均匀性及偏离场对分割效果影响较小，宜采用稀疏性对称约束，反之则宜都采用稠密性对称约束处理图像。

为了评估所提方法对大脑MR多模态图像配准的性能，我们选用BrainWeb数据集中正常大脑T1、T2、PD模态图像做了刚体配准实验。其中数据集的参数为：切片厚度为1mm，噪声强度为3％(基于最大图像灰度值)，灰度不均匀性为40％。在旋转实验中，±20°为旋转角度范围；平移实验中，平移范围为[-20,20]mm。本次实验所选用的对比方法有互信息(mutual information,MI)配准、基于上下文自相似(self-similarity context,SSC)配准、基于模态独立局部二值纹理描述子(modality independent local binary pattern,miLBP)配准、基于鲁棒自相似描述子(robust self-similarity descriptor,RSSD)配准。本实验中，所提方法对各模态图像均使用稀疏性对称约束，其中根据阈值处理的分割图像，其阈值δ＝0.8；稀疏性对称约束函数中邻域大小为3×3，即A＝3。实验结果如表2所示。

表2多模态刚性配准实验

由表2验可看出，提出的方法的配准性能均明显好于其他方法。其中提出的方法T1-T2及T1-PD的配准性能相似，PD-T2几乎可达到无误差配准。由此可见，使用提出的方法对多模态MR大脑图像刚体配准具有优良的性能。此外，综合表1和表2可发现，在RIRE数据集上，T1、PD模态图像用提出的方法得到的待配准图像配准效果较好，T2模态图像则相对较差；在BrainWeb数据集上，T2、PD模态图像用提出的方法得到的待配准图像配准效果较好，T1模态图像则相对较差。由此可见，此方法受各模态拍摄环境造成的噪声影响较大。

在SSD相似性度量准则配准框架下，提出了一种基于图像分割及对称自相似的大脑MR多模态医学图像配准方案，该方法将不同模态图像转换为统一的人工模态，很好地保持了图像的一致性。提出的方法较之基于SSD相似性度量准则下用自相似描述子(如SSC、RSSD、miLBP)来处理各模态图像的配准方案，在配准精度上有明显优势。

应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。本实施例中未明确的各组成部分均可用现有技术加以实现。

Claims (9)

1.一种大脑MR医学图像配准方法，其特征在于：所述配准方法包括以下步骤

步骤1：采用BCFCM方法对参考图像和浮动图像进行分割；

步骤2：采用MSR目标对称轴检测方法对原始参考图像和原始浮动图像进行对称轴检测并提取对称轴方程；

步骤3：对步骤1处理后图像进行二值化阈值分割；

步骤4：根据步骤2所检测的对称轴方程对步骤3处理后图像进行近似对称性约束，所得到的图像分别记为待配准的参考图像和待配准的浮动图像；

步骤5：初始化刚体变换矩阵配准参数，设定参数初始值为0；

步骤6：通过SSD相似性度量准则配准步骤4得到的待配准的参考图像和待配准的浮动图像，得到最优变换矩阵配准参数。

2.根据权利要求1所述的一种大脑MR医学图像配准方法，其特征在于：所述BCFCM方法对参考图像和浮动图像进行分割，其目标函数为：

其中，y_i为观测到的图像像素，表示该像素受偏离场影响；β_i代表偏离场，x_i＝y_i-β_i为去偏离场后的图像像素；α为比例因子，控制像素邻域对图像分割的影响；

表示以第i个像素点为邻域中心的邻域方形窗像素点集合，不包括i本身，其中R为邻域方形窗的边长；y_r及β_r分别表示领域像素值及该像素值的偏离场；N表示像素总个数；c_j表示聚类簇心，C表示待划分像素的类别总数，待划分像素类别总数C≥2；其中一个聚类簇心初始值设为0，记为c_j＝0＝0；m为一个大于0的任意实数；u_ij为隶属度，表示第i个像素点属于第j类的概率，满足

在图像分割中，该值组成的矩阵U称为分割矩阵；利用(1)式进行优化得到最优隶属度

并组合成最优分割矩阵U^*。

3.根据权利要求2所述的一种大脑MR医学图像配准方法，其特征在于：选择U^*中类别为背景灰度值的分割子矩阵

作为原始图像的分割矩阵来分割图像，其格式为：

设分割图像的像素矩阵大小为x×y，则：

其中，N表示像素总数；

为最优隶属度，表示第i个像素点属于第j类的最优概率。

4.根据权利要求1所述一种大脑MR医学图像配准方法，其特征在于：步骤2所述的MSR目标对称轴检测方法流程如下：

(1)将过图像中心点的垂线置为初始对称轴；

(2)将原始图像关于该对称轴进行镜像变换，得到镜面图像；

(3)将原始图像作为参考图像，镜面图像作为浮动图像进行刚性配准，得到刚性变换矩阵参数；

(4)根据刚性变换矩阵参数得到原始图像对称轴。

5.根据权利要求1所述的一种大脑MR医学图像配准方法，其特征在于：将步骤1处理后的图像进行阈值分割，经过二值化阈值分割处理后的图像分割矩阵记为P，格式如下所示：

其中p(g)为阈值分割函数，δ为阈值，在理想情况下δ＝1。

6.根据权利要求1所述的一种大脑MR医学图像配准方法，其特征在于：根据步骤2所检测的对称轴方程，对步骤3处理后图像进行近似对称性约束，所得到的图像分别记为待配准的参考图像和浮动图像，对称轴方程记为l，记矩阵P经对称性约束后的矩阵为S，公式如下所示：

实际使用中令s(h)＝s_sparse(h)或者s(h)＝s_dense(h)，其中s_sparse(h)和s_dense(h)分别为稀疏性和稠密性对称性约束函数，其定义式为：

h′为关于图像对称轴l的对称点，h′_a为对称点邻域像素值，A²为对称点邻域方形窗大小，表示对称点邻域方形窗内像素点个数。

7.根据权利要求1所述的一种大脑MR医学图像配准方法，其特征在于：步骤5所述初始化刚体变换矩阵配准参数，设定参数初始值为0；所述刚体变换的过程为：

将具有6个自由度的刚体变换作为人脑多模态MR医学图像配准的几何变换参数矩阵，所述6个自由度的刚体包括t_x、t_y、t_z三个位移自由度，θ_x、θ_y、θ_z三个旋转自由度。

8.根据权利要求1所述的一种大脑MR医学图像配准方法，其特征在于：步骤6，通过SSD相似性度量准则配准图像，得到最优变换矩阵配准参数，其中

SSD相似性度量优化公式为：

其中，T表示由6个自由度组成的刚体变换参数矩阵；M、N代表图像像素位置，即M×N代表图像像素总个数；I和I′分别代表待配准浮动图像和待配准参考图像像素值大小矩阵。

9.根据权利要求8所述的一种大脑MR医学图像配准方法，其特征在于：

配准过程的优化目标为：

在对配准浮动图像进行刚性变换即计算T(i，j))时，由于其值并不一定为整数坐标，故通过双线性插值法，使其映射到整数坐标；该配准通过高斯牛顿法求解(9)式最优值获得最佳刚体变换矩阵参数。

Images