基于互信息技术和遗传算法的数字图像配准

互信息是信息论的一个基本概念，是两个随机变量统计相关性的测度。Woods用测试图像的条件熵作为配准的测度，用于PET 到MR 图像的配准。Collignon 、Wells[1] 等人用互信息作为多模态医学图像的配准测度。以互信息作为两幅图像的相似性测度进行配准时，如果两幅基于共同解剖结构的图像达到最佳配准时，它们对应的图像特征互信息应为最大。最大互信息法几乎可以用在任何不同模式图像的配准中，特别是当其中一个图像的数据部分缺损时，所以这种方法广泛用于多模态图像的配准中。但是，当待匹配图像是低分辨率、图像包含的信息不够充分或两幅待匹配图像的重叠部分较少时，基于互信息的配准目标函数就会极不光滑，出现较多局部最优解，为目标函数最优解的搜索带来较大的难度。但由于该测度不需要对不同成像模式下图像灰度间的关系作任何假设，也不需要对图像进行分割或任何预处理，因此，该测度可以被广泛地应用于CT-MR，PET-MR 等多种图像的配准工作。

基于最大互信息的图像配准取得了很大的成功，但是它也存在一些缺陷，例如任何图像和一幅只有一种灰度值的图像（例如灰度值为255 的全黑图像）配准，无论几何变换怎样，他们的联合灰度直方图都是一样的。因此在实际情况中我们就会碰到一个问题，一幅图像可能包括大量的同类区域（例如天空、大海等），那么这样的图像就不太适合用最大互信息的方法进行配准，实际上基于图像灰度的配准方法都不是太适合这样的情况。此外，基于最大互信息的图像配准因为要进行全局参数优化搜索，配准时间也比较长。Rangarajan 等提出了一种利用互信息匹配形状特征点进行配准的策略. 该策略针对待配准的两幅图像，首先分别提取出形状特征点的集合，并定义这两个集合它们的互信息，然后使之最大化，以达到配准。

由于角点是景物轮廓线上曲率的局部极大点，对掌握景物的轮廓特征具有决定作用。一旦找到了景物的轮廓特征点也就大致掌握了景物的形状。虽然角点相对于其他的特征来说比较少，但是使用大量特征实现配准必会导致算法复杂度的提高。出于提高速度而又不降低配准精度的考虑，角点是一个很好的配准特征。因此，在本文中我们将要对基于角点特征的图像配准做一些初步的探讨。

我们将研究的范围限制在二维脑断层图像的配准. 因为脑组织受到颅骨的严密保护，所以脑部运动可以近似为刚体运动，即内部无相对运动. 同时，假设待配准的图像经过预处理后具有相同的空间比例.我们的目标是寻求空间变换T，使MI(T) 最大.。针对前面所做假设，令T = T1*T2；其中，T1 为平移矩阵，T2为旋转矩阵。

最近邻插值法的精确度很低，而双线性插值法会产生新的灰度值。这对于联合直方图的统计是不利的。因为新加入的灰度值使得随着Tα的一些小变动，联合直方图中就会增加新的象素对，或者减少象素对，从而互信息值变化比较大，也就是互信息函数曲线会不光滑，这样不利于优化。因此为了消除新产生的灰度值的不利影响，我们在配准过程中引入了另外一种插值法：PV（Partial Volume）插值法。从产生插值图像这个方面来说，PV 插值法不能算是一种插值方法，它是专门针对联合直方图的更新而设计的。和双线性插值法一样，PV 插值法也是利用点Tα(X)的四个最近邻点和权值 ，可是，不同于双线性插值法的是，PV 插值法不是根据最近邻点的加权平均所得到的灰度值，从而更新联合直方图，而是根据权值 使周围四个象素点都贡献于联合直方图的统计，如图1所示。

可用公式表示为：

(2-1)

(2-2)

图1 pv插值法

由于角点是景物轮廓线上曲率的局部极大点，对掌握景物的轮廓特征具有决定作用。一旦找到了景物的轮廓特征点也就大致掌握了景物的形状。直观的讲，角点就是图像上所显示的物体边缘拐角所在的位置点。

Harris角点检测法是一种基于图像灰度的检测方法，这类方法主要通过计算点的曲率及梯度来检测角点。该方法是由Harris和Stephen提出来的，也叫Plessey角点检测法。其基本思想与Moravec角点算子相似，但对其作了许多改进。

Moravec角点算子计算各象素沿小同方向的平均灰度变化，选取最小值作为对应象素点的角点响应函数。定义在一定范围内具有最大角点响应的象素点为角点。假设图像的灰度定义为I那么平移(x，y)所得到的灰度变化的计算公式为：

(2-3)

这里W表示图像窗口，平移(x，y)表示了四个方向：水平、垂直、对角线和反对角线，即(0，1)，(1，0)，(1，1)，(-1，1)。

Moravec角点算子简单快速，但是它存在一些缺点，Harris角点算子正是针对这些缺点做了很大的改进。

首先，Moravec角点算子是各向异性的，因为它的角点响应只计算了四个方向。故为了包含所有的方向，Harris角点算子对式(2-3)进行了展开：

(2-4)

这里一阶微分可以由下面的式子近似：

(2-5)

因此，E可以表示

(2-6)

这里

(2-7)

为了避免噪声的影响，这里w采用高斯平滑窗口：

(2-8)

其次，Moravec角点算子对强边界敏感，这是因为它的响应值只考虑了E的最小值。Harris角点算子则利用了E在平移方向上的变化。

在平移方向(x，y)上的E可以表示如下

（2-9）

这里2×2的矩阵M为

(2-10)

可以看出，E和局部自相关函数联系非常紧密。设α，β为矩阵M的特征值，则α，β与局部自相关函数的主曲率成比例。当两个曲率都低时，局部自相关函数是平坦的，那么窗口图像区域的灰度值近似为常量；当只有一个曲率高而另一个曲率低时，局部自相关函数呈脊形，那么E只有当沿山脊移动时变化小，这就表示是边缘；当两个曲率都高时，局部自相关函数是尖峰，那么E在任意方向上移动都会增加，这就表示是角点。因此我们可以由α，β的值判断是否是角点。为了不对M进行分解求特征值，可以采用Tr(M)和Det(M)来代替α，β，其中

(2-11)

从而形成对矩阵M与旋转无关的描述：

(2-12)