-
摘要: 经典的Active demons算法利用参考图像和浮动图像的梯度信息作为驱动力，并使用均化系数调节两种驱动力之间的强度.该算法克服了Demons算法单一使用参考图像的梯度信息作为驱动力的缺点，但是Active demons算法中的均化系数无法同时兼顾大形变和小形变区域的准确配准，还会导致配准的收敛速度和精确度相互制约的问题.为此，本文提出一种新的Active demons非刚性配准算法.提出的算法在Active demons扩散方程中引入一个称为平衡系数的新参数，与均化系数联合调整驱动力，不仅可以兼顾图像中同时具有的大形变和小形变区域的准确配准，而且在一定程度上缓和了收敛速度和精确度相互制约的问题.为了进一步提高配准的收敛速度和精确度，避免陷入局部极值，在新的配准算法的实现中引入由粗到细的多分辨率策略.在Checkboard测试图像、自然图像和医学图像上的实验结果表明，提出的算法较经典的Active demons算法收敛速度更快，配准精度平均提高了54.28%，接近最新的TV-L1光流场图像配准算法的配准精度，解决了Active demons算法存在的问题.
-
关键词:
- 非刚性配准 /
- Active demons算法 /
- 光流场图像配准 /
- 驱动力 /
- 多分辨率策略
Abstract: Classic active demons algorithm uses gradient information of the static image and the moving image as driving forces, and uses a homogeneous coefficient to adjust their intensities. Although the algorithm overcomes the disadvantage of the demons algorithm using the gradient information of a single static image, the homogeneous coefficient of the active demons algorithm can not accurately handle registration with both large deformation and small deformation, and will cause the mutual restraint problem of convergence speed and registration accuracy. In order to solve this problem, this paper presents a non-rigid registration algorithm based on active demons algorithm, which introduces a new parameter called balance coefficient to the active demons algorithm to adjust the driving force in combination with the homogeneous coefficient. Not only can the large deformation and small deformation be taken into account at the same time, but also the mutual restraint problem of speed and accuracy can be eased to a certain extent. In order to further improve registration accuracy and convergence speed and avoid falling into local extremes, a coarse-to-fine multi-resolution strategy is introduced into the registration process. Experiments on checkboard test images, natural images and medical images demonstrate that the proposed algorithm is faster and more accurate. The registration accuracy is improved by 54.28% on average, and is close to that of the latest TV-L1 optical flow image registration algorithm. -
图像配准是图像处理和计算机视觉等领域的热点研究问题, 广泛应用于遥感图像拼接、目标定位和疾病诊断治疗等方面.一般我们将配准过程中固定不变的图像称为参考图像, 将施加变换的图像称为浮动图像.图像配准通过进行参考图像和浮动图像的空间变换, 最终使得参考图像和浮动图像对齐.一般来说, 图像配准在不同领域存在一定差异, 比如:将多幅遥感图像拼接成一幅更大的图像, 多幅图像之间有重叠部分, 重叠部分的特征提取是配准的关键; 而在医学图像配准中, 两幅图像之间几乎是完全重叠的, 研究关注的核心问题是图像的非均匀形变.按照配准过程中图像是否存在形变来划分, 图像配准分为刚性配准和非刚性配准.其中, 非刚性图像配准是近年来的研究热点之一[1].非刚性配准算法可以分为基于特征和基于灰度的两类算法.虽然基于特征的配准算法运算速度较快, 但是特征提取和特征匹配往往又是一个困难的问题, 对那些特征不清晰的图像, 该类算法存在一定的局限性.而基于灰度的配准算法直接利用图像灰度的梯度信息, 避免了特征提取过程中可能涉及的人工干涉, 完全自动处理, 受到越来越多的青睐.
基于光流场理论[2]的Demons算法最初由Thirion[3]提出, 是一种基于灰度的全自动配准算法.该算法将配准视作扩散问题, 在图像之间的形变较小的前提下, 可利用参考图像灰度的梯度和两幅图像的差值来估计浮动图像的形变. Demons算法具有完备的数学理论基础和高效的配准性能, 广泛应用于各种图像配准中, 特别是医学图像的配准[3-10]. Hellier等[11]对常用的6种配准算法对比研究表明, Demons算法配准精度高于其他5种配准算法.然而, Demons算法存在以下不足: 1)仅采用参考图像的梯度信息驱动形变, 收敛速度慢; 2)在参考图像梯度信息接近于零时, 浮动图像的形变方向不能确定, 将导致错误的配准变换; 3) Demons算法建立在图像之间存在小形变的假设之上, 对于形变较大的图像, Demons算法基于图像灰度很难准确估计形变.
针对以上的问题, Rogelj等[12]和Wang等[13]提出Active demons算法, 将浮动图像的梯度信息也引入扩散方程, 将来自两幅图像的两个单向力叠加为一个合力, 引入均化系数α调节驱动力的强度.该算法很好地克服了Demons算法无法配准形变相对较大的图像的缺陷, 且大幅提高了配准速度和精确度.但是, Active demons算法仅通过均化系数α调节驱动力的强度, 也存在很大的缺陷:当α取值较大时, 驱动力较小, 适合小形变区域的配准, 不适合大形变的区域配准; 当α取值较小时, 驱动力较大, 收敛速度很快, 适合大形变区域的快速配准, 但对于小形变区域, 则容易产生过矫正的问题, 导致图像的整体配准精度较低.
在实际图像处理中, 尤其是在医学图像配准中, 两帧图像之间通常同时存在大形变和小形变, 这导致Active demons算法经常会陷入两难的境地.综合以上分析得出, 均化系数α调节驱动力无法兼顾大形变和小形变区域的配准.我们应该寻找一种可以解决该问题的方法, 使其既能够兼顾大形变和小形变区域的配准, 同时也能够提高配准精度和收敛速度.
为了提高Demons算法的配准精度, 许多学者对Demons算法进行了不同的改进, 这些改进算法可以分为3种策略.策略1保证形变场的拓扑保持性.最具代表性的学者是来自法国的Vercauteren等, 他们提出了一种标准的配准模型[14-15], 将基于图像灰度的Demons配准算法看成一个能量优化过程, 通过建立一个适当的目标函数, 对其进行优化搜索得到相应的配准参数.为了保证形变场的拓扑保持性, 避免图像配准时产生物理上不合理的形变, Vercauteren根据李群理论提出了Diffeomorphic demons配准算法.基于该理论, 后人提出了很多的改进算法.例如, Lorenzi等[16]提出在形变场中增加SVFs (Stationary velocity fields)约束和使用新的相似性测度函数LCC (Local correlation coefficient).另外, 林相波等[17]从矢量场的特点出发, 详细分析了该算法得到的形变场拓扑性质, 并在此基础上给出了强化形变场拓扑保持性的方法.策略2是在驱动力计算公式中增加新的信息项, 弥补仅依靠灰度梯度信息配准图像的不足; 常见的信息项有灰度梯度场的相似性[18]、对称正交梯度信息[19]、几何形状约束[20]等.策略3从分力的影响程度上开展研究, 如Rogelj等[12]提出对称梯度平均化的Symmetric demons算法; Lin等[21]对两个分力进行了取舍, 如果某个分力指向图像匹配方向, 则对其进行加强, 否则将其设为零, 取消该分力的影响.文献[21]中的算法有利于提高配准的精度, 但是分力的取舍需要较长的时间.
为了加快配准的收敛速度, Sharp等[22]研究了Demons算法基于硬件的加速方法, 利用GPU中的并行计算能力, 大幅减少算法运行时间.也有一些学者采用多分辨率策略[21, 23], 将低分辨率层的空间变换作为高分辨率层的初始形变矩阵, 按照由粗到细的方式进行逐层配准.文献[24]提出一种在配准迭代过程中采用自适应地调节均化系数的算法.
Pock等[25]及Sánchez等[26]在传统光流模型的基础上, 利用全变差正则项(Total variation regularization, TV正则项)和鲁邦性数据项(L1范数项)构建能量函数得到的TV-L1光流场图像配准算法的配准精度较高, 与Demons类算法相比运算时间缩短不小, 但是同样不能解决大形变的问题.
本文在Active demons算法的基础上, 提出了一种能够同时提高配准精度和收敛速度的非刚性配准算法.主要思路如下: 1)在Active demons算法的驱动力计算公式中引入平衡系数k, 该系数联合均化系数α精细地调整驱动力的强度, 使得Active demons算法适合形变范围更大的图像配准; 2)通过实验分析, 给出平衡系数k的一般取值范围, 增强算法的普适性; 3)引入了多分辨率策略.在实验测试分析中, 采用配准精度和收敛速度指标, 将提出的算法分别与经典的Active demons配准算法和TV-L1光流场图像配准算法进行对比.
本文剩余部分结构安排如下:在第1节, 本文回顾经典的Demons算法以及Active demons算法; 在第2节, 介绍本文提出的改进Active demons算法; 在第3节, 通过对Checkboard测试图像、自然图像和医学图像的实验验证本文算法的有效性; 第4节对本文工作进行总结.
1. 方法与原理
1.1 Demons算法
Demons算法[3]将参考图像和浮动图像看成是连续运动图像序列中的两帧, 并估计从浮动图像到参考图像的形变向量u.通常该向量理解为施加在浮动图像上的形变驱动力.假设图像在运动的过程中保持灰度不变, 对于空间中任意一点(x, y), 如果参考图像
S和浮动图像M上的灰度分别为S(x, y)和M(x, y), $\nabla S(x, y)$ 为参考图像在(x, y)处的梯度值, 根据光流场方程可以得到下式: $\begin{align} u\nabla S(x, y)=M(x, y)-S(x, y) \end{align}$
(1) 则u可以用下式估计:
$\begin{align} u = \frac {M(x, y)-S(x, y)}{|\nabla S(x, y)|^2} \nabla S(x, y) \end{align}$
(2) 显然, 当 $|\nabla S(x, y)|^2$ 很小时, u可能无限大, 与实际情况不符.为了解决该问题, 在式(2)分母上增加一新项, 变为
$\begin{align} u= \frac {M(x, y)-S(x, y)}{|\nabla S(x, y)|^2+(M(x, y)-S(x, y)^2)} \nabla S(x, y) \end{align}$
(3) 式(3)可以保证在参考图像灰度变化较小的地方对浮动图像的位移向量u接近于零, 即驱动力为零.
在Demons算法中, 采用迭代方式得到最终变换.在第n次迭代时, 对每一点(x, y), 根据式(3)计算该点的形变驱动力 $u(x, y)$ , 得到该点的瞬时位移 $v_n (x, y)=-u_n (x, y)$ , 计算形变场 $T_n (p)=T_{n-1} (p)+v_n (p)$ .此外, 为了使该变换在全局范围内连续, 在迭代过程中, 使用高斯滤波对瞬时位移进行平滑.
1.2 Active demons算法
在Demons算法中, 使浮动图像发生形变的力单纯取自参考图像的梯度信息.这样的模型只适合于处理小形变问题, 不能满足较大形变图像配准应用的需要.为了扩大该算法的应用范围, Wang等[13]根据牛顿第三定律的作用力与反作用力的原理, 提出了Active demons算法, 将浮动图像的梯度信息也作为一种正内力引入到驱动力计算公式中, 得到下式:
$\begin{split} u= \frac{M(x, y)-S(x, y)}{|\nabla S(x, y)|^2 +(M(x, y)-S(x, y))^2} \nabla S(x, y)+ \\ \frac{M(x, y)-S(x, y)}{|\nabla M(x, y)|^2 +(M(x, y)-S(x, y))^2} \nabla M(x, y) \end{split}$
(4) 为了能够调整驱动力的强度, 在上式中引入了一个均化系数α, 式(4)变为
$u= \frac{M(x, y)-S(x, y)}{|\nabla S(x, y)|^2 +\alpha^2(M(x, y)-S(x, y))^2} \nabla S(x, y)+ \\ \frac{M(x, y)-S(x, y)}{|\nabla M(x, y)|^2 +\alpha^2(M(x, y)-S(x, y))^2} \nabla M(x, y)$
(5) 文献[27]详细分析了α的取值对配准结果的影响, 该分析结果表明:较小的α可以加快算法的收敛速度, 减少配准时间, 但配准精度较差; 较大的α会增加配准时间, 但得到的配准结果较好.
2. 基于Active demons算法的非刚性配准方法
2.1 一种改进的Active demons算法
本文针对Active demons算法在图像配准中存在的问题, 提出了一种能够同时提高配准精度和收敛速度的改进算法.在Active demons扩散方程中, 为了能够调节驱动力, Wang在分母中灰度差的平方项上引入了一个均化系数α.受此启发, 将该思想应用在另一分母项上, 即在梯度的模值的平方项上引入一个新的系数k, 我们称之为平衡系数.然后通过调整均化系数α和平衡系数k的取值, 联合调节驱动力的强度, 较依靠单一系数能够更加精细地控制配准过程中的的形变程度.改进后的Active demons扩散方程如下所示:
$u= \frac{M(x, y)-S(x, y)}{k^2|\nabla S(x, y)|^2 +\alpha^2(M(x, y)-S(x, y))^2} \nabla S(x, y)+\\ \frac{M(x, y)-S(x, y)}{k^2|\nabla M(x, y)|^2 +\alpha^2(M(x, y)-S(x, y))^2} \nabla M(x, y)$
(6) 式(6)中, 如果均化系数α不变, 平衡系数k对驱动力u的大致影响可以分为三种情况:当 $k < 1$ 时, 分母项减小, u变大, 形变程度增大; 当 $k>1$ 时, 分母项增大, u变小, 形变程度减小; 当 $k=1$ 时, 式(6)等价于Active demons扩散方程的形式.下一节将具体分析平衡系数k的取值对配准结果的影响.
2.2 平衡系数的引入及取值分析
2.2.1 引入平衡系数的有效性分析
为了清晰地了解参数k对配准结果的影响, 本节采用简单的二值图像进行配准.如图 1所示, 正方形为浮动图像, 圆为参考图像.实验参数设置如下: 1)图像分辨率为256 × 256; 2) Demons估计浮动图像中所有点的形变大小; 3)初始形变矩阵为零, 空间变换为自由项变换; 4)最大迭代次数为200, 单尺度实现; 5)图像插值使用三线性插值法.
为了进行详细的分析, 我们先使用经典的Active demons算法进行配准实验, 均化系数α取10个不同的值, 依次为0.05, 0.1, 0.4, 0.5, 0.6, 1, 1.5, 2, 2.5, 3, 得到迭代过程的均方差曲线图, 如图 2所示.图 2表明:随着迭代次数的增加, 均方差基本上保持单调递减, 收敛后保持微幅振动. α越小, 配准的收敛速度越快, α较大时收敛变慢, 但是配准精度更高, 这也与文献[27]的分析相符.
实验条件不变, 采用本文提出的改进Active demons算法进行配准实验.平衡系数k取0.5, 实验结果的均方差与迭代次数的关系曲线图如图 3所示.为了比较两种算法的配准精度, 表 1为两种算法迭代200次时得到的均方差(Mean square error, MSE)值. AD为Active demons算法, IAD (Improved active demons)为引入平衡系数k之后的Active demons算法.可以很明显地看到, IAD算法的均方差值均小于AD算法, 这说明改进算法配准精度更高.为了考察收敛速度, 分别从两个实验中取出4组具有代表性的实验数据进行对比分析. 图 4是两种算法的配准曲线对比图, 其中, α取为0.05, 0.4, 1, 2, 图 4中带标志的曲线表示Active demons算法, 未带标志的曲线表示改进算法, 相同的线型代表相同的均化系数.从曲线图中可以明显看到, 改进后的Active demons算法的配准曲线与经典的Active demons算法相比更贴近y轴, 这说明本文算法的收敛速度更快.
表 1 两种算法的配准结果均方差对比Table 1 The comparison of two registration algorithms on the MSEα 0.05 0.1 0.4 0.5 0.6 AD(×10−4) 3.55 4.02 6.49 8.18 10 IAD(×10−4) 2.25 2.35 3.71 4.29 5.27 α 1.0 1.5 2 2.5 3 AD(×10−4) 34 107 187 255 311 IAD(×10−4) 12 47 113 183 246 以上分析表明:本文在Active demons算法中引入平衡系数k是非常正确的, 不仅提高配准的精确度, 收敛速度也更快.那么, 对于固定的某一均化系数α, 该如何确定平衡系数k的值才可使配准性能达到最好呢?接下来, 本文将具体分析这个问题.
2.2.2 最优平衡系数的取值分析
从理论上分析, 在式(6)中, 因为 $|k\nabla S(x, y)|^2 +\alpha^2 (M(x, y)-S(x, y))^2 \geq 2\alpha k \nabla S(x, y)(M(x, y)-S(x, y))$ , 所以驱动力的上限为 $u\leq (1/\alpha k)$ .若α不变, u与k成反比. k越大, u就越小, 形变程度小, 配准的精度高, 但收敛速度较慢. k越小, u越大, 允许的形变度越大, 收敛速度很快, 但配准精度较低.
在实验中, 设定α为0.5, 均化系数k取10个不同的值(0.05, 0.1, 0.4, 0.5, 0.6, 1.0, 1.5, 2, 2.5, 3), 考察平衡系数对配准结果的影响, 其他实验条件与上节中相同, 结果如图 5~图 7所示.图 5是取不同参数值时的配准结果, 图 6是相应的配准结果的差值图, 图 7是配准过程的收敛曲线图.从图 5和图 6的实验结果可以明显看到, 在200次迭代条件下, 当α取0.5时, k取0.5的配准结果的均方误差最小, 说明配准精度最高; 从图 7的曲线图中可以发现, 随着迭代次数的增加, 均方差基本上保持单调递减, 收敛后保持微幅振动. k越小, 配准的收敛速度越快; k越大, 收敛速度变慢, 但是配准精度更高.这也与上面的理论分析相符.
图 7 平衡系数对改进Active demons算法的影响Fig. 7 The impact of balance coefficient k on the improved active demons registration algorithm图 8为不同α和k组合情况下, 以固定迭代200次之后的配准结果的均方差关系曲面图, 其最小值点位于α=0.5, k=0.5处.从实验结果中可以看出, α和k较小时, 配准精度高平滑性较差; α和k较大时, 配准精度相对低.建议二者的取值在[0.5, 1.5]之间较适宜.
图 8 均方差与均化系数和平衡系数的关系曲线Fig. 8 The relations of the mean square error with the α and the k2.3 基于多分辨率策略实现
Active demons算法容易陷入局部极小值, 而且该算法的运行速度慢、配准精度不高.多分辨率策略可以避免陷入局部极值, 在提高算法的准确性、速度和鲁棒性方面都很有帮助, 其基本思路如下:
1) 通过降采样的方法将待配准图像分解成多种不同分辨率的图像;
2) 在低分辨率图像上快速进行粗配准;
3) 对在低分辨率图像上得到的形变场进行升采样, 将其作为更高一级分辨率的初始变换;
4) 逐层的迭代, 最终实现整幅图像的配准.
为了进一步提高配准的收敛速度、准确度以及避免陷入局部极值, 本文基于多分辨率策略实现, 具体实现步骤如下所示.
算法1. 基于多分辨率策略实现的配准过程 1:设定初始形变矩阵 $T_0$ , 最大迭代次数k, 最大分解尺度L, 当前分解级数l, 当前迭代次数n;
2:对参考图像S和浮动图像M进行L级尺度分解;
3: ${\rm while}\ (l < L)\ {\rm do}$
4: ${\rm while}\ (n < N)\ {\rm do}$
5: 根据式(6)计算第n次迭代后浮动图像M的位移向量 ${\pmb u}{_n^l}$ , 并更新形变矩阵 $T_l=T_{l-1}+u{_n^l}$ , 将 $T_l$ 应用于浮动图像M使其形变;
6: $n=n+1$
7: end while
8: $l=l+1, \, n=0$
9:对前一级图像的形变矩阵 $T_{l-1}$ 升采样, 作为当前图像配准的初始形变矩阵;
10: end while
11:将最终的形变矩阵T应用于浮动图像M, 应用双三次插值得到配准结果.
3. 实验结果分析
本文先用主观评价方法对所得配准结果进行分析, 并获得初步的质量评价结果, 再用5种客观评价方法进一步准确地分析配准效果.
为了分析和验证本文提出算法的实际效果, 分别采用测试图像、自然图像和医学图像进行实验.实验中同时采用了Demons算法、Active demons算法(AD)、单尺度实现的改进算法(Single-scale implementation active demons, SIAD)、多分辨率实现的改进算法(Multi-scale implementation active demons, MIAD)和TV-L1光流场算法(TV-L1).为了公平起见, 在下面的实验中, 各种算法的参数选择都是相对于该算法的最佳结果.
3.1 Checkboard测试图像对比分析
图 9和图 10是针对典型的Checkboard测试图像的配准实验, 浮动图像为参考图像模拟球状扭曲而得到的, 现在欲通过配准将其复原, 以此来检验各种配准算法的性能.图 9(a)为浮动图像; 图 9(b)为参考图像; 图 9(c)为Demons算法的配准结果; 图 9(d)为Active demons算法的配准结果, α为0.5; 图 9(e)为SIAD算法的配准结果, α为0.5, k为0.5; 图 9(f)为MIAD算法的配准结果, α为0.5, k为0.5, 分解尺度为3; 图 9(g)为TV-L1光流场算法的配准结果.为了能清晰地分辨配准结果的优劣, 图 10为配准后的浮动图像与参考图像的差值图, 图 9(a)~(f)分别为原始差值、Demons算法、Active demons算法、SIAD算法、MIAD算法以及TV-L1光流场算法的差值图.表 2为5种算法的配准结果的客观比较.
图 10 图像的配准结果与参考图像的差值Fig. 10 The differences between the registration results and the static image表 2 配准结果的客观分析Table 2 The objective analysis of registration results评价方法 Demons AD SIAD MIAD TV-L1 均方差(×10−4) 43 40 29 25 11 相互系数(%) 99.16 99.24 99.44 99.52 99.98 峰值信噪比 54.43 55.31 58.36 59.97 87.74 归一化互信息 1.42 1.42 1.43 1.43 0.71 结构相似度(%) 95.42 95.74 96.82 97.41 99.97 通过对比图 9和图 10, 可以发现, Demons算法由于仅依靠参考图像的梯度信息配准图像, 配准结果与参考图像的误差最大, 配准效果最差; Active demons算法的配准结果相比Demons算法有了一定的改善, 尤其是对Demons算法无法配准的大形变区域实现了较好的配准; 本文提出的MIAD算法的配准误差最小, 效果最好; SIAD算法的配准结果次之, 但也好于原来的两种配准算法.虽然TV-L1算法配准的效果最好, 但其与参考图像的差值明显有不规则的边缘效应.
从表 2的5种客观评价方法中可以看出, MIAD算法的均方差值最小, 相关系数最大, 峰值信噪比最大, 归一化互信息最大, 结构相似度也是最大的. SIAD算法的5种评价值也均优于Demons算法和Active demons算法, 这也与前面的主观分析一致.综合以上分析说明, 对于测试图像, 无论从配准图像的主观效果和客观指标, 显然目前最新的TV-L1光流场算法是最好的, 但本文提出的SIAD算法、 MIAD算法的配准性能与之非常接近, 且好于原来的同类算法.
3.2 自然图像对比分析
下面是针对自然图像的配准实验, 本文选取两种不同形变程度的图像分别进行配准实验.其中, Lena图像代表大形变图像配准, 自然彩色图像代表小形变图像配准.
图 11和图 12是针对典型的大形变Lena图像进行的配准实验.图 11(a)为浮动图像; 图 11(b)为参考图像; 图 11(c)为Demons算法的配准结果; 图 11(d)为Active demons算法的配准结果, α为1.5;图 11(e)为SIAD算法配准结果, α为1, k为1.5;图 11(f)为MIAD算法的配准结果, α为1, k为1.5, 分解尺度为3;图 11(g)为TV-L1算法配准结果.图 12为配准后的浮动图像与参考图像的差值图, 图 12(a)~(f)分别对应初始差值、Demons算法、Active demons算法、SIAD算法、MIAD算法和TV-L1算法的差值图.表 3为5种算法的配准结果的客观比较.
表 3 配准结果的客观分析Table 3 The objective analysis of registration results评价方法 Demons AD SIAD MIAD TV-L1 均方差(×10−4) 13 8.11 6.19 2.81 4.61 相互系数(%) 96.35 97.69 98.25 99.22 99.87 峰值信噪比 66.59 71.15 73.87 81.75 91.49 归一化互信息 1.37 1.29 1.39 1.41 3.96 结构相似度(%) 88.55 90.99 92.56 96.28 99.11 从图 11和图 12的配准结果差值图中可以明显看出, Demons算法完成了图像中小形变部分的配准, 在一定程度上纠正了浮动图像中面部五官的形变, 但头发部位的形变加强了, 眉毛和帽子也发生了畸变; Active demons算法相比Demons算法在处理大形变图像配准问题上优越性明显, 尤其是头发部位的矫正.但对于形变较小的四周区域, 在经过了Active demons算法配准之后, 明显形变过大; SIAD算法不仅解决了Active demons算法存在的小形变区域过配准的问题, 而且保留了大形变区域配准的准确性, 配准误差更小, MIAD算法的配准误差最小, 配准的差值图很光滑, 证明了多分辨率策略实现的优越性.从配准图像及与参考图像的差值来看, TV-L1算法配准的效果是最好的, 但其差值图边缘具有不规则的边缘效应.
图 12 图像的配准结果与参考图像的差值Fig. 12 The differences between the registration results and the static image从表 3中的5种客观评价方法可以看出, 对于TV-L1光流场算法, 除了均方误差指标不及本文提出的算法外, 其他的指标都是最好的.同理, 对于本文提出的算法, 除了TV-L1外, 其5个指标均优于其他同类算法, 特别是结构相似度, 比其他算法都高出很多.这也进一步表明了本文算法优越性, 能够较好地保持图像的拓扑结构.
图 13是针对小形变的彩色自然图像进行的配准实验, 浮动图像只出现了局部小形变.图 13(a)为浮动图像; 图 13(b)为参考图像; 图 13(c)为Demons算法的配准结果; 图 13(d)为Active demons算法的配准结果, α为2;图 13(e)为SIAD算法配准结果, α为0.5, k为1.5;图 13(f)为MIAD算法的配准结果, α为0.5, k为1.5, 分解尺度为3;图 13(g)为TV-L1算法配准结果.由于该测试图像是彩色小形变图像, 图像的信息丰富, 由5种配准算法得到的配准图像的配准精度都较高, 从图 13中的各种算法配准结果很难能够辨别出差异; 同样, 由这些算法得到的配准图像与参考图像的差值图, 依靠肉眼观察差异性较小, 整体背景将是黑色, 所以这里就不把这些差值图列出来.我们可以依靠表 4中5种算法的配准结果的客观指标比较进行评判.
表 4 配准结果的客观分析Table 4 The objective analysis of registration results评价方法 Demons AD SIAD MIAD TV-L1 R 5.09 5.03 5.12 18.69 2.50 均方差(×10−5) G 5.07 5.10 5.33 17.68 2.47 B 4.69 4.79 5.01 16.71 3.20 R 99.84 99.92 99.92 99.90 99.99 相互系数(%) G 99.86 99.86 99.85 99.83 99.98 B 99.92 99.84 99.83 99.81 99.71 R 91.07 91.12 91.04 90.19 94.16 峰值信噪比 G 91.08 91.05 90.86 90.43 94.20 B 91.42 91.33 91.14 90.67 93.08 R 3.25 3.26 3.25 3.27 3.98 归一化互信息 G 3.17 3.17 3.16 3.17 4.03 B 3.02 3.02 3.00 3.01 3.63 R 98.87 97.49 97.56 97.71 99.21 结构相似度(%) G 98.99 97.46 97.41 97.44 99.27 B 98.70 97.30 97.23 97.23 98.89 从表 4中的5种客观指标(R、G、B彩色通道分列)可以看出, Demons算法、Active demons算法、SIAD算法、TV-L1算法的配准均方误差均明显小于MIAD算法的配准均方误差.但整体来讲, 这5种算法的其余4个客观指标基本相当, 这说明对于小形变的彩色自然图像, 本文提出的SIAD算法与Demons算法、Active demons算法的效果基本等价, 且提出的MIAD算法不如Demons算法、Active demons算法和SIAD算法要好, 这主要是由于Demons算法、Active demons算法本身擅长小形变区域的配准.由于采用的测试图像是变形小的彩色图像, 且信息量丰富, MIAD算法的优势显现不出来.主要原因是由于图像形变太小, 呈现不出MIAD算法的优势.
3.3 医学图像的对比分析
为了考察5种算法在医学图像上的配准性能, 本节进行了2组试验.前1组是来源于网络检验配准算法优劣的常用实验数据, 代表小形变图像配准的测试.后1组为真实的实验数据, 为磁共振肝图像的位置矫正实验, 代表大形变图像配准的测试.
图 14和图 15是对医学图像的小形变配准实验.图 14(a)为参考图像; 图 14(b)为浮动图像; 图 14(c)为Demons算法的配准结果; 图 14(d)为Active demons算法的配准结果, α为2;图 14(e)为SIAD算法的配准结果, α为0.5, k为2;图 14(f)为MIAD算法的配准结果, α为0.5, k为2, 分解尺度为3;图 14(g)为TV-L1算法配准结果.图 15为配准后的浮动图像与参考图像的差值图, 表 6为5种算法的配准结果的客观比较.
图 15 图像的配准结果与参考图像的差值图Fig. 15 The differences between the registration results and the static image表 5 配准结果的客观分析Table 5 The objective analysis of registration results评价方法 Demons AD SIAD MIAD TV-L1 均方差(×10−5) 4.78 7.53 3.06 0.81 8.56 相互系数(%) 99.81 99.70 99.88 99.97 99.97 峰值信噪比 99.51 94.94 103.94 117.26 98.81 归一化互信息 1.43 1.31 1.44 1.54 3.49 结构相似度(%) 98.71 98.44 99.03 99.76 99.34 表 6 配准结果的客观分析Table 6 The objective analysis of registration results评价方法 Demons AD SIAD MIAD TV-L1 均方差(×10−4) 30 19 16 14 11.4 相互系数(%) 97.21 98.26 98.50 98.75 99.06 峰值信噪比 58.08 62.81 64.29 66.05 77.56 归一化互信息 1.38 1.39 1.41 1.43 2.90 结构相似度(%) 83.55 89.27 89.61 94.40 98.77 通过视觉观察图 14和图 15可以发现, Demons算法的配准效果比Active demons算法要好; Active demons算法单纯依靠均化系数调节形变程度, 会导致小形变区域矫正过大, 产生误配准的问题; SIAD算法的配准误差相比Active demons算法和Demons算法明显减少, 差值图更平滑; 基于多分辨率策略实现的MIAD算法的配准差值图总体比较平滑, 没有出现较大的配准误差区域. TV-L1算法的配准差值图是最平滑的, 但表 6中其均方误差及峰值信噪比要次于MIAD算法.
从表 5的5种客观评价方法可以看出, Active demons算法的各项客观评价指标确实不如Demons算法, 配准性能最差, 这也与前面的视觉观察相呼应.除了TV-L1的归一化互信息最高外, MIAD算法的其他4种客观评价指标最好, 所以配准效果相对最优.
肝癌患者接受分次放射治疗的过程中, 治疗部位的位置和形状可能发生变化, 靶区形状以及靶区与周围危及器官的位置关系也会发生变化, 因此需要利用图像配准技术来跟踪形变[27].图 16和图 17是采用真实的磁共振肝图像进行的配准实验.图 16(a)为参考图像; 图 16(b)为浮动图像; 图 16(c)为Demons算法的配准结果; 图 16(d)为Active demons算法的配准结果, α为1;图 16(e)为SIAD算法的配准结果, α为0.4, k为1.5;图 16(f)为MIAD算法的配准结果, α为0.5, k为2, 分解尺度为3;图 16(g)为TV-L1算法配准结果.图 17为5种算法的配准结果与参考图像的差值图, 表 7为5种算法的配准结果的客观比较.
图 17 图像的配准结果与参考图像的差值图Fig. 17 The differences between the registration results and the static image从图 16和图 17中可以看出, TV-L1算法的配准效果是最好的; Demons算法对于图像中部的大形变区域与参考图像存在很大的偏差; Active demons算法比Demons算法有了一定程度的改善, 但在轮廓的小形变区域里出现了明显的误配准现象; SIAD算法相比Active demons算法, 配准误差减小了很多; MIAD算法得到的配准效果最好, 无论是对大形变区域还是小形变区域, 均实现了准确的配准.表 6的客观分析验证了前面的主观分析, TV-L1算法是最好的, MIAD算法其次, 进一步说明了本文提出的算法相比原来的算法在医学图像配准应用上更有优势.
为了衡量某种算法在配准精度方面优于其他算法的程度, 文献[18]定义了 $P_e$ 来权衡, 计算公式如下:
$\begin{align} P_e =\frac {MSE_e -MSE_p}{MSE_e}\times 100% \end{align}$
(7) 其中, $MSE_e$ 为已有的配准算法的均方差, $MSE_p$ 为提出算法的均方差.
本文也采用上面的方法来衡量提出的算法在配准精度上的优越性.由于Demons算法的配准精度不如Active demons算法, 所以本文仅与Active demons算法相比较.关于实验结果的精确度分析如表 7所示, 其中最后两列表示本文提出的两种算法相对于Active demons算法在配准精确度上的提高程度, 经过计算分别平均提高了28.45%、54.28%.
表 7 几种算法的配准结果均方差值比较Table 7 The objective analysis of difference algorithmsCases AD SIAD MIAD Pe(s) Pe(M) (×10−4) (×10−4) (×10−4) (%) (%) 1 40 29 25 27.5 37.5 2 8.11 6.19 2.81 23.00 65.00 3 0.497 0.515 1.769 16.67 53.33 4 0.753 0.306 0.081 59.32 89.28 5 19 16 14 15.79 26.32 Average — — — 28.45 54.28 4. 结束语
针对Active demons算法存在的问题, 本文提出了一种改进算法.经过在Checkboard测试图像、自然图像和医学图像上的实验表明, 本文算法的图像配准精确度要好于Active demons算法、Demons算法, 接近于最新的TV-L1光流场配准算法, 配准的速度较Active demons算法和Demons算法大大提高, 但不如TV-L1光流场配准算法快.对于相对大形变图像的配准, 改进算法解决了Active demons算法可以配准大形变区域而无法兼顾小形变区域的配准问题.对于小形变图像的配准, 改进算法提高了Active demons算法的配准精度, 实现了小形变区域的更准确的配准.所以本文提出的算法比经典的Active demons算法更有优势, 是一种有效的处理非刚性配准的方法.
应该说无论是TV-L1光流场配准算法、Demons算法、Active demons算法, 还是本文提出的新算法都是基于灰度信息的, 针对小形变的图像配准有明显效果, 相对大形变的图像还有基于图像特征的配准算法有优势.尽管如此, 本文算法也存在如下不足.由于改进算法引入了新参数, 在进行配准实验时, 为了达到高配准精度, 针对不同的图像类型, 可能需要进行一定的参数取值试验工作, 以选择最优的配准结果.
-
图 7 平衡系数对改进Active demons算法的影响
Fig. 7 The impact of balance coefficient k on the improved active demons registration algorithm
图 8 均方差与均化系数和平衡系数的关系曲线
Fig. 8 The relations of the mean square error with the α and the k
图 10 图像的配准结果与参考图像的差值
Fig. 10 The differences between the registration results and the static image
图 12 图像的配准结果与参考图像的差值
Fig. 12 The differences between the registration results and the static image
图 15 图像的配准结果与参考图像的差值图
Fig. 15 The differences between the registration results and the static image
图 17 图像的配准结果与参考图像的差值图
Fig. 17 The differences between the registration results and the static image
表 1 两种算法的配准结果均方差对比
Table 1 The comparison of two registration algorithms on the MSE
α 0.05 0.1 0.4 0.5 0.6 AD(×10−4) 3.55 4.02 6.49 8.18 10 IAD(×10−4) 2.25 2.35 3.71 4.29 5.27 α 1.0 1.5 2 2.5 3 AD(×10−4) 34 107 187 255 311 IAD(×10−4) 12 47 113 183 246 
下载: 导出CSV
表 2 配准结果的客观分析
Table 2 The objective analysis of registration results
评价方法 Demons AD SIAD MIAD TV-L1 均方差(×10−4) 43 40 29 25 11 相互系数(%) 99.16 99.24 99.44 99.52 99.98 峰值信噪比 54.43 55.31 58.36 59.97 87.74 归一化互信息 1.42 1.42 1.43 1.43 0.71 结构相似度(%) 95.42 95.74 96.82 97.41 99.97
下载: 导出CSV
表 3 配准结果的客观分析
Table 3 The objective analysis of registration results
评价方法 Demons AD SIAD MIAD TV-L1 均方差(×10−4) 13 8.11 6.19 2.81 4.61 相互系数(%) 96.35 97.69 98.25 99.22 99.87 峰值信噪比 66.59 71.15 73.87 81.75 91.49 归一化互信息 1.37 1.29 1.39 1.41 3.96 结构相似度(%) 88.55 90.99 92.56 96.28 99.11
下载: 导出CSV
表 4 配准结果的客观分析
Table 4 The objective analysis of registration results
评价方法 Demons AD SIAD MIAD TV-L1 R 5.09 5.03 5.12 18.69 2.50 均方差(×10−5) G 5.07 5.10 5.33 17.68 2.47 B 4.69 4.79 5.01 16.71 3.20 R 99.84 99.92 99.92 99.90 99.99 相互系数(%) G 99.86 99.86 99.85 99.83 99.98 B 99.92 99.84 99.83 99.81 99.71 R 91.07 91.12 91.04 90.19 94.16 峰值信噪比 G 91.08 91.05 90.86 90.43 94.20 B 91.42 91.33 91.14 90.67 93.08 R 3.25 3.26 3.25 3.27 3.98 归一化互信息 G 3.17 3.17 3.16 3.17 4.03 B 3.02 3.02 3.00 3.01 3.63 R 98.87 97.49 97.56 97.71 99.21 结构相似度(%) G 98.99 97.46 97.41 97.44 99.27 B 98.70 97.30 97.23 97.23 98.89
下载: 导出CSV
表 5 配准结果的客观分析
Table 5 The objective analysis of registration results
评价方法 Demons AD SIAD MIAD TV-L1 均方差(×10−5) 4.78 7.53 3.06 0.81 8.56 相互系数(%) 99.81 99.70 99.88 99.97 99.97 峰值信噪比 99.51 94.94 103.94 117.26 98.81 归一化互信息 1.43 1.31 1.44 1.54 3.49 结构相似度(%) 98.71 98.44 99.03 99.76 99.34
下载: 导出CSV
表 6 配准结果的客观分析
Table 6 The objective analysis of registration results
评价方法 Demons AD SIAD MIAD TV-L1 均方差(×10−4) 30 19 16 14 11.4 相互系数(%) 97.21 98.26 98.50 98.75 99.06 峰值信噪比 58.08 62.81 64.29 66.05 77.56 归一化互信息 1.38 1.39 1.41 1.43 2.90 结构相似度(%) 83.55 89.27 89.61 94.40 98.77
下载: 导出CSV
表 7 几种算法的配准结果均方差值比较
Table 7 The objective analysis of difference algorithms
Cases AD SIAD MIAD Pe(s) Pe(M) (×10−4) (×10−4) (×10−4) (%) (%) 1 40 29 25 27.5 37.5 2 8.11 6.19 2.81 23.00 65.00 3 0.497 0.515 1.769 16.67 53.33 4 0.753 0.306 0.081 59.32 89.28 5 19 16 14 15.79 26.32 Average — — — 28.45 54.28
下载: 导出CSV
-
[1] Sotiras A, Davatzikos C, Paragios N. Deformable medical image registration: a survey. IEEE Transactions on Medical Imaging, 2013, 32(7): 1153-1190 doi: 10.1109/TMI.2013.2265603 [2] Cachier P, Pennec X, Ayache N. Fast non rigid matching by gradient descent: study and improvements of the "Demons" algorithm. Technical Report-3706, India, 1999. [3] Thirion J P. Image matching as a diffusion process: an analogy with Maxwell's demons. Medical Image Analysis, 1998, 2(3): 243-260 doi: 10.1016/S1361-8415(98)80022-4 [4] Brown L G. A survey of image registration techniques. ACM Computing Surveys, 1992, 24(4): 325-376 doi: 10.1145/146370.146374 [5] Lester H, Arridge S R. A survey of hierarchical non-linear medical image registration. Pattern Recognition, 1999, 32(1): 129-149 doi: 10.1016/S0031-3203(98)00095-8 [6] Jie T, Xue J, Dai T K, Chen J, Zheng J. A novel software platform for medical image processing and analyzing. IEEE Transactions on Information Technology in Biomedicine, 2008, 12(6): 800-812 doi: 10.1109/TITB.2008.926395 [7] Cao Z L, Dong E Q, Zheng Q, Sun W Y, Li Z Z. Accurate inverse-consistent symmetric optical flow for 4D CT lung registration. Biomedical Signal Processing and Control, 2016, 24: 25-33 doi: 10.1016/j.bspc.2015.09.005 [8] Cao Z L, Dong E Q. Distinctive local binary pattern for non-rigid registration of lung computed tomography images. Electronics Letters, 2015, 51(22): 1742-1744 doi: 10.1049/el.2015.0598 [9] Nithiananthan S, Schafer S, Mirota D J, Stayman J W, Zbijewski W, Reh D D, Gallia G L, Siewerdsena J H. Extra-dimensional demons: a method for incorporating missing tissue in deformable image registration. Medical Physics, 2012, 39(9): 5718-5731 doi: 10.1118/1.4747270 [10] Reaungamornrat S, Wang A S, Uneri A, Otake Y, Khanna A J, Siewerdsen J H. Deformable image registration with local rigidity constraints for cone-beam CT-guided spine surgery. Physics in Medicine and Biology, 2014, 59(14): 3761-3787 doi: 10.1088/0031-9155/59/14/3761 [11] Hellier P, Barillot C, Corouge I, Gibaud B, Le Goualher G, Collins D L, Evans A, Malandain G, Ayache N, Christensen G E, Johnson H J. Retrospective evaluation of intersubject brain registration. IEEE Transactions on Medical Imaging, 2003, 22(9): 1120-1130 doi: 10.1109/TMI.2003.816961 [12] Rogelj P, Kovačič S. Symmetric image registration. Medical Image Analysis, 2006, 10(3): 484-493 doi: 10.1016/j.media.2005.03.003 [13] Wang H, Dong L, O'Daniel J, Mohan, Garden A S, Ang K K, Kuban D A, Bonnen M, Chang J Y, Cheung R. Validation of an accelerated'demons' algorithm for deformable image registration in radiation therapy. Physics in Medicine and Biology, 2005, 50(12): 2887-2905 doi: 10.1088/0031-9155/50/12/011 [14] Vercauteren T, Pennec X, Perchant A, Ayache N. Non-parametric diffeomorphic image registration with the demons algorithm. In: Proceedings of 10th International Conference on Medical Image Computing and Computer-Assisted Intervention (MICCAI 2007). Brisbane, Australia, 2007. 319-326 [15] Vercauteren T, Pennec X, Perchant A, Ayache N. Diffeomorphic demons: efficient non-parametric image registration. NeuroImage, 2009, 45: S61-S72 http://cn.bing.com/academic/profile?id=2103857226&encoded=0&v=paper_preview&mkt=zh-cn [16] Lorenzi M, Ayache N, Frisoni G B, Pennec X. LCC-Demons: a robust and accurate symmetric diffeomorphic registration algorithm. NeuroImage, 2013, 81: 470-483 doi: 10.1016/j.neuroimage.2013.04.114 [17] 林相波, 邱天爽, 阮素, Nicolier F. Demons非刚性配准算法拓扑保持性的研究.自动化学报, 2010, 36(1): 179-183 doi: 10.3724/SP.J.1004.2010.00179Lin Xiang-Bo, Qiu Tian-Shuang, Ruan Su, Nicolier F. Research on the topology preservation of the demons non-rigid registration algorithm. Acta Automatica Sinica, 2010, 36(1): 179-183 doi: 10.3724/SP.J.1004.2010.00179 [18] Liu X Z, Yuan Z M, Zhu J M, Xu D R. Medical image registration by combining global and local information: a chain-type diffeomorphic demons algorithm. Physics in Medicine and Biology, 2013, 58(23): 8359-8378 doi: 10.1088/0031-9155/58/23/8359 [19] Lu C, Mandal M. Improved image registration technique based on demons and symmetric orthogonal gradient information. In: Proceeding of the 2010 International Conference on Signal Processing and Communications (SPCOM). Bangalore: IEEE, 2010. 1-5 [20] Cazoulat G, Simon A, Dumenil A, Gnep K, de Crevoisier R, Acosta O, Haigron P. Surface-constrained nonrigid registration for dose monitoring in prostate cancer radiotherapy. IEEE Transactions on Medical Imaging, 2014, 33(7): 1464-1474 doi: 10.1109/TMI.2014.2314574 [21] Lin X B, Qiu T S, Nicolier F, Ruan S. An improved method of'Demons'non-rigid image registration algorithm. In: Processing of the 9th International Conference on Signal Processing. Beijing, China: IEEE, 2008. 1091-1094 [22] Sharp G C, Kandasamy N, Singh H, Folkert M. GPU-based streaming architectures for fast cone-beam CT image reconstruction and demons deformable registration. Physics in Medicine and Biology, 2007, 52(19): 5771-5783 doi: 10.1088/0031-9155/52/19/003 [23] Guimond A, Roche A, Ayache N, Meunier J. Three-dimensional multimodal brain warping using the demons algorithm and adaptive intensity corrections. IEEE Transactions on Medical Imaging, 2001, 20(1): 58-69 doi: 10.1109/42.906425 [24] Xu Sheng-Zhou, Song En-Min, Xu Xiang-Yang. Non-rigid mammogram registration based on improved demons algorithm. Journal of Image and Graphics, 2009, 14(12): 2566-2571 http://cn.bing.com/academic/profile?id=2352702219&encoded=0&v=paper_preview&mkt=zh-cn [25] 徐胜舟, 宋恩民, 许向阳.基于改进Demons算法的乳腺X线摄片非刚性配准.中国图象图形学报, 2009, 14(12): 2566-2571 http://www.cqvip.com/qk/90287a/200912/32418755.htmlPock T, Urschler M, Zach C, Beichel R, Bischof H. A duality based algorithm for TV-L^1-optical-flow image registration. In: Proceedings of the 10th International Conference on Medical Image Computing and Computer Assisted Intervention (MICCAI). Brisbane, Australia: Springer, 2007. 511-518 http://www.cqvip.com/qk/90287a/200912/32418755.html [26] Sánchez J, Meinhardt-Llopis E, Facciolo G. TV-L1 optical flow estimation. Image Processing on Line, 2013, 3: 137-150 doi: 10.5201/ipol [27] Lin Xiang-Bo, Qiu Tian-Shuang, Nicolier F, Ruan Su. The study of active demons algorithm for deformable image registration. Chinese Journal of Biomedical Engineering, 2008, 27(4): 636-640 http://cn.bing.com/academic/profile?id=2381391515&encoded=0&v=paper_preview&mkt=zh-cn -
下载:
计量
- 文章访问数: 1820
- HTML全文浏览量: 291
- PDF下载量: 1539
- 被引次数: 0
