CN102034115A - 基于马尔可夫随机场模型与非局部先验的图像配准方法 - Google Patents

Abstract

一种基于马尔可夫随机场模型与非局部先验的图像配准方法，包括以下步骤：(1)分别读入待配准的目标图像和浮动图像；(2)计算待配准的目标图像和浮动图像的差的平方和，将其均方距离作为相似性测度；(3)计算位移场的非局部先验信息作为正则项对位移场进行平滑约束；(4)将相似性测度与非局部先验信息相加建立马尔可夫随机场模型，将配准转化为求解马尔可夫随机场的能量函数最小的问题；(5)采用序列加权树信息传递算法求解上述能量函数的最小值；(6)以马尔可夫随机场的能量函数作为目标函数，搜索目标函数的最小值，马尔可夫随机场的能量函数最小时，完成配准。本发明在图像空间分辨率较低，有噪声影响等情况下，该方法配准精度高、鲁棒性强。

Description

基于马尔可夫随机场模型与非局部先验的图像配准方法

技术领域

本发明涉及一种图像配准方法，具体涉及基于马尔可夫随机场模型与非局部先验的图像配准方法。

背景技术

图像配准是现代图像处理技术应用的一个重要方面，它是指对不同时间、不同视场、不同成像模式的两幅或多幅图像进行空间几何变换，以使代表相同位置、结构的像素或体素在几何上能够匹配对应起来。图像配准的主要目的是去除或者抑制待配准图像和参考图像之间几何上的不一致，包括平移、旋转、缩放和弹性变形。它是图像分析和处理的关键步骤，是图像对比、数据融合、变化分析和目标识别的必要前提。配准技术主要应用在遥感图像处理、医学图像处理、制图学、计算机视觉、目标识别和军事目的等。

对于在不同时间、不同视场、不同成像模式等不同条件下获取的两幅图像进行配准处理，就是要定义一个配准测度函数，即相似性测度，寻找一个空间变换关系，使得经过该空间变换后，两幅图像间的相似性达到最大(或者差异性达到最小)，即两幅图像得到空间几何上的一致。常用的空间变换形式主要有仿射变换、透视变换、曲线变换和弹性变换(如流体模型、光流模型和由B样条构成的变形模型)等几种。配准时通过搜索相似性测度函数的全局最大值来得到变换参数，常见的相似性测度函数有相关系数、互相关函数、互信息量等；对于差异性测度函数，是通过搜索相似性测度函数的全局最小值来得到变换的参数，常见的差异性测度函数有距离函数、总绝对差函数、总平方差函数等。由于空间变换包含多个参数，图像配准是一个多参数最优化问题。常用的有Powell算法、Gauss-Newton法、L-M算法和最速下降法等。如果相似性测度函数的性质不好，存在大量的局部极值，经典的优化算法就会终止在局部极值，得到错误的配准参数。解决的办法是采用模拟退火、遗传算法、最大流/最小割算法和粒子群优化算法等现代优化算法，以克服局部极值问题。但是这些算法的收敛速度往往很慢，为了得较好的结果，就不得不以牺牲实时性为代价。到为了加快速度，有时也使用多分辨率的由粗到细的优化搜索策略，以减少搜索空间，提高搜索效率。

目前已有多种图像配准方法，大体可以分为基于图像特征和基于图像灰度的配准方法。基于图像特征的方法一般通过寻找图像中比较明显的解剖结构特征来计算变换参数，提取的特征包括点、线、边缘和轮廓等。该方法计算效率比较高。特别是基于轮廓的方法，尽管不同模态图像像素的灰度分布特性之间有很大的差异，但物体的一些明显轮廓在两幅图像中均能得到较好的保持，这些轮廓特征可以用来作为参照来进行多模图像间的配准。但是基于图像特征的方法的配准精度取决于特征提取的准确性与否。大部分情况下，医学图像特征点的位置比较复杂，很难进行准确地提取。基于图像灰度的配准方法不需要提取特征点，而是利用图像的灰度信息进行配准。由于充分利用了全部的灰度信息，因此一般更容易得到较精确的配准结果。其中互相关法、最小平方差法、互信息量法是几种常用的方法。互相关法一般实现起来简单，但是计算代价庞大；最小平方差法对图像数据的部分缺失和图像的背景噪声不是非常敏感，具有较强的鲁棒性，但由于需要通过迭代搜寻最优参数，故其执行效率也不是很高，而且其解严重依赖于初始值的选择，不恰当的初始值往往会使解陷入局部极值。互信息量是信息论中的一个测度，用来度量两个随机变量之间的相似性。在多模态医学图像配准问题中，虽然两幅图像来源于不同的成像设备，但是它们基于人体共同的解剖信息，所以当两幅图像的空间位置完全一致时，它们的对应象素的灰度互信息达到最大值。作为一种相似性测度，互信息量取得了巨大的成功，特别在医学图像领域。互信息配准方法目前被公认为是配准精度和鲁棒性最好的回溯性配准方法之一。但互信息量的方法没有考虑图像的空间和方向信息，在图像空间分辨率较低，有噪声影响等情况下，该算法会陷入局部极值。

发明内容

本发明的目的是提供一种基于马尔可夫随机场模型与非局部先验的图像配准方法，在图像空间分辨率较低，有噪声影响等情况下，该方法配准精度高、鲁棒性强。

本发明实现上述目的的技术解决方案是：

一种基于马尔可夫随机场模型与非局部先验的图像配准方法，包括以下步骤：

(1)分别读入待配准的目标图像和浮动图像；

(2)计算待配准的目标图像和浮动图像的差的平方和，将其均方距离作为相似性测度；

(3)计算浮动图像的位移场的非局部先验信息作为正则项对位移场进行平滑约束；

(4)将(2)中的相似性测度与(3)中浮动图像的位移场的非局部先验信息相加建立马尔可夫随机场模型，将配准转化为求解马尔可夫随机场的能量函数最小的过程；

(5)采用序列加权树信息传递算法求解该能量函数的最小值；

(6)以马尔可夫随机场的能量函数作为目标函数，搜索目标函数的最小值，马尔可夫随机场的能量函数最小时，完成配准。

本发明中所述的步骤1中，通过线性变换将目标图像和浮动图像所有像素点的灰度值变换到0～255的范围内，并分别采用中值滤波对输入图像进行降低噪声的预处理。

所述步骤3中，非局部先验信息是通过非局部滤波器得到。

所述步骤5中，采用序列加权树信息传递算法求解能量最小问题，序列树中每棵树的权重相同，假设每一条边(edge)只属于一条链(chain)，引入辅助变量

分别表示向前和向后信息，初值设为0，下边界LB(θ)＝∑_imin_XE(X|θⁱ)初值设为0，对下边界进行向前和向后消息迭代更新。

本发明的有益效果是：

(1).本发明将图像的位移场视为马尔可夫随机场，建立基于马尔可夫随机场肝脏图像配准模型。以固定图像与浮动图像相减后的平方和为数据项，以非局部先验作为正则项对位移场进行平滑约束。这就将配准转化为寻找Bayesian最大后验的最优化问题，为了加强求解的速度与稳定性，采用序列加权树信息传递算法对该问题进行求解。在提高配准精度的同时，减少了配准所需的时间，减少了陷入局部极值的可能。实验结果表明在图像空间分辨率较低，有噪声影响等情况下，该算法具有精度高、鲁棒性强的特点，在无需人工干预和预处理条件下能自动实现图像的精确配准，特别适合于医学、卫星、航空等图像的配准。

(2).二次先验的过平滑效应和非二次先验所导致的阶梯状伪影将就都是由于无法有效的区分边缘信息和噪声。而这些缺陷正是由于局部邻域内像素点的简单加权仅能提供十分有限的先验信息。非局部先验利用了较大区域内，甚至是整体的图像信息，因此最终能达到比较高的配准精度。

附图说明

图1为本发明的步骤流程图；

图2为本发明的非局部先验滤波器图；

图3a为人体腹部注入造影剂前扫描的CT图像；

图3b为人体腹部注入造影剂后扫描的CT图像。

具体实施方式

本发明的具体实施方式如图1～图3所示，在本实施例中，结合一套人体腹部注入造影剂前、后扫描的CT图像(参见图3a～3b)，详细阐述本发明的工作步骤。

步骤1，分别读入注入造影剂前、后扫描的CT图像，分别作为待配准的目标图像和浮动图像。通过线性变换(I-I_min)×255/(I_max-I_min)将目标图像和浮动图像所有像素点的灰度值变换到0～255的范围内，其中I为图像灰度值，即对注入造影剂前扫描的CT图像进行线性变换时，I为注入造影剂前扫描的CT图像的灰度值；对注入造影剂后扫描的CT图像进行线性变换时，I为注入造影剂后扫描的CT图像的灰度值。再分别采用中值滤波对输入图像进行降低噪声的预处理。

步骤2，计算待配准的目标图像和浮动图像的差的平方和，即均方距离作为相似性测度。

步骤3，计算浮动图像的位移场的非局部先验信息作为正则项对位移场进行平滑约束，利用非局部先验滤波器滤波得到非局部先验信息，具体如下：

在去噪模型中，像素点x_s滤波后的灰度值为大邻域V_s内图像灰度的加权平均

，其中v(·)表示灰度函数，则v(x_t)表示像素点x_t处的灰度值，ω(x_s，x_t)表示分配给v(x_t)的权重；对于V_s内的任意一个像素x_t，权重ω(x_s，x_t)与像素点x_s，x_t的邻域N_s，N_t之间的距离d(v(N_s)，v(N_t))有关：

其中为归一化常数；σ为估计的噪声标准差，本专利中可规定σ＝1；h为滤波器参数；d为邻域内所有像素点的高斯加权Euclidian距离，如图2所示。将图像的位移场(矢量)看作含噪声图像，计算θ_st(x_s，x_t)时采用这种非局部滤波方法必然可以引入位移场中更多的信息，从而达到更好的配准效果。只需将灰度函数v(·)改成在该像素点的位移矢量，即v(x_t)表示像素点x_t处的位移。

步骤4，将相似性测度与非局部先验信息相加建立马尔可夫随机场模型，将配准转化为求解马尔可夫随机场的能量函数最小的过程。

步骤5，采用序列加权树信息传递算法求解该能量最小问题，序列树中每棵树的权重相同，假设每一条边(edge)只属于一条链(chain)，引入辅助变量分别表示向前和向后信息，初值设为0，下边界LB(θ)＝∑_imin_XE(X|θⁱ)初值设为0，对下边界进行向前和向后消息迭代更新。

步骤6，以马尔可夫随机场的能量函数作为目标函数，在搜索目标函数的最小值时，要对浮动图像先进行刚性或者弹性变换，再计算相似性，而对浮动图像每变换一次都要进行插值，采用双线性插值方法对浮动图像的变换图像进行插值，当目标函数的最小时，完成配准。

本发明的实施方式不限于此，上述一些参数的获得也可以用类似的计算公式获得，如采用互信息量作为相似性测度，或者采用梯度下降法、粒子群优化等优化算法搜索最优解。在此不一一列举；根据本发明的上述内容，按照本领域的普通技术知识和惯用手段，在不脱离本发明上述基本技术思想前提下，本发明还可以做出其它多种形式的等效修改、替换或变更，均可实现本发明目的。

Claims (3)

1.一种基于马尔可夫随机场模型与非局部先验的图像配准方法，其特征在于包括以下步骤：

(1)分别读入待配准的目标图像和浮动图像；

(2)计算待配准的目标图像和浮动图像的差的平方和，将其均方距离作为相似性测度；

(3)计算浮动图像的位移场的非局部先验信息作为正则项对位移场进行平滑约束；

(4)将步骤(2)中的相似性测度与步骤(3)中浮动图像的位移场的非局部先验信息相加建立马尔可夫随机场模型，将配准转化为求解马尔可夫随机场的能量函数最小的过程；

(5)采用序列加权树信息传递算法求解上述能量函数的最小值；

(6)以马尔可夫随机场的能量函数作为目标函数，搜索目标函数的最小值，马尔可夫随机场的能量函数最小时，完成配准。

2.权利要求1所述的基于马尔可夫随机场模型与非局部先验的图像配准方法，其特征在于所述的步骤3中，非局部先验信息由非局部滤波器得到。

3.根据权利要求2所述的基于马尔可夫随机场模型与非局部先验的图像配准方法，其特征在于所述的步骤5中，所采用的序列加权树信息传递算法中，其序列树中每棵树的权重相同，设每一条边只属于一条链，引入辅助变量

分别表示向前和向后信息，初值设为0，下边界LB(θ)＝∑_imin_XE(X|θⁱ)初值设为0，对下边界进行向前和向后消息迭代更新。

Images