改进Demons算法的非刚性医学图像配准

基于改进Demons算法的三维肺部医学影像配准的研究发表时间：2017-06-20T10:48:44.433Z 来源：《健康世界》2017年第7期作者：侯立川王凯君[导读] 本研究认为基于改进Demons算法能够对三维肺部医学影像进行有效配准，明显提高肺部影像资料质量，值得进行深入推广应用。

山东省威海市文登整骨医院放射科 264400摘要：目的：探讨基于改进Demons算法的三维肺部医学影像配准。

方法：获取1例肺部疾病患者一个呼吸运动周期内十个时相的肺部CT影像资料，根据均方误差确定实验所需时相影像，最终确定第10时相为最大吸气相，第6时相作为最大呼气相，应用混合层次配准模型进行影像配准，以最大吸气相作为参考影像，以最大呼气相作为浮动影像，完成全局图像间非刚性配准，然后应用改进Demons算法进行局部优化配准，采用的算法策略为多分辨率策略，对比配准前后影像间的均方误差值。

结果：配准后全局影像间的均方误差值为10.7418，较配准前的21.9306明显降低，局部影像间的均方误差值为3.8935，较配准前的13.0452明显降低。

结论：基于改进Demons算法能够对三维肺部医学影像进行有效配准。

关键词：改进Demons算法；三维肺部医学影像；配准；应用价值影像学技术在辅助疾病诊断和指导疾病治疗中发挥着十分重要的作用。

影像学资料的质量直接影响着临床医师对诊断信息的收集，如何能够提高影像学资料质量是我国医疗领域长期研究的一项课题[1]。

配准（Registration）是指同一区域内通过不同成像方法获得的不同图形的的地理坐标的匹配，是对复杂图像信息进行综合分析的有效手段[2]。

本研究主要探讨基于改进Demons算法的三维肺部医学影像配准，旨在提高肺部影像学诊断资料的质量，为肺部疾病的诊断、治疗及预后评估提供客观参考依据。

材料与方法1.1一般材料获取1例肺部疾病患者一个呼吸运动周期内10个时相的肺部CT影像资料进行本次研究。

非刚性医学图像配准的必要性
由于人体内部器官和组织的形态和位置会随 时间、姿势、呼吸等因素而发生变化，因此 非刚性医学图像配准在临床应用中具有重要 意义。
非刚性医学图像配准能够更好地适应人体内 部器官和组织的形态和位置变化，提高图像 配准的精度和可靠性，从而更好地支持疾病
诊断和治疗方案制定。
非刚性医学图像配准的关键技术
基于深度学习的配准方法利用深度神经网络对图 像进行学习和配准，具有强大的自适应能力和较 高的配准精度，但需要大量的训练数据和计算资 源。
03 非刚性医学图像 配准算法
基于特征的方法
特征提取
提取医学图像中的局部特征，如边缘、角点等，作为 配准的依据。
特征匹配
通过匹配图像中的特征点，确定图像间的对应关系。
研究内容和方法
本研究的主要内容是针对非刚性医学图像配准技术进行研究，提出一种基于特征提取和优化算法的非刚性医学图像配准方法 。
首先，对医学图像进行预处理，包括去噪、分割、提取特征等步骤；然后，利用特征匹配算法将两幅图像进行初步配准；最 后，采用优化算法对初步配准的结果进行优化，得到最终的配准结果。
02 非刚性医学图像 配准概述
医学图像配准的定义
医学图像配准是一种基于图像处理技 术的过程，旨在将不同时间点、不同 设备或不同个体之间的医学图像进行 对齐和比较，以便准确地进行疾病诊 断和治疗方案制定。
VS
医学图像配准通常涉及对图像进行平 移、旋转、缩放、剪切等变换，以使 不同图像之间的对应关系达到最佳匹 配。
随着精准医疗的发展，非刚性医学图像配准技术将更加注重个体差异，针对不同患者和疾病状态开展精细化配准。同时，随着可穿戴设备和移动设备的普及 ，如何在这些设备上实现非刚性医学图像配准也是一个值得研究的方向。

第36卷第1期自动化学报Vol.36,No.1 2010年1月ACTA AUTOMATICA SINICA January,2010Demons非刚性配准算法拓扑保持性的研究林相波1,2邱天爽1阮素2NICOLIER Frédéric2摘要在基于配准的图像分割应用中,拓扑保持性是非刚性图像配准算法的一个重要约束.本文从矢量场特性出发,分析了Demons非刚性图像配准算法导致目标拓扑改变时变形场的特点.根据变形场特点与其雅可比行列式之间的关系,给出了校正该算法拓扑保持性的方法.实验表明,改进后的变形场具有了拓扑保持性.关键词非刚性配准,拓扑,雅可比行列式,矢量场分析DOI10.3724/SP.J.1004.2010.00179Research on the Topology Preservation of the Demons Non-rigid RegistrationAlgorithmLIN Xiang-Bo1,2QIU Tian-Shuang1RUAN Su2NICOLIER Frédéric2Abstract Topology preservation is an important constraint of image non-rigid registration algorithms in the applications of registration based image segmentation.In this paper,the defor-mationfield of the Demons non-rigid registration that make the object s topology change is analyzed according to the features of the vectorfield.An method that can enforce the algorithm s topology preservation is provided on the basis of relation between the deformationfield and its Jacobian.The experiment results show that the modified deformationfield is topology preserving. Key words Non-rigid registration,topology,Jacobian deter-minant,vectorfield analysis非刚性医学图像配准是近年来的一个研究热点,实现不同个体间的配准是其中的一个主题,它的一个重要应用是利用图谱提供的信息,分割不同个体的组织器官结构[1].图像配准是通过将模板图像变形为目标图像,实现两图像中相同结构间的良好对应,这种对应性用配准得到的空间变形场T 表达.设可以变形的模板图像为R,目标图像为S,则T表示从模板图像的支撑域ΩR向目标图像的支撑域ΩS变形的一个空间变换,即T:ΩR→ΩS.通常可以通过迭代方式求解能量最小化问题得到最优解T,其配准能量的一般表达式如下E(S,R◦T)=E dis(S,R◦T)+E reg(T)(1)其中E dis(S,R◦T)表示图像特征之间的匹配程度,可以直接采用图像的灰度信息[2],也可以采用图像的其他特征[3];收稿日期2008-10-06录用日期2009-04-13Manuscript received October6,2008;accepted April13,2009国家自然科学基金(30170259,30570475,60372081),教育部博士点基金项目(20050141025)资助Supported by National Natural Science Foundation of China (30170259,30570475,60372081),Ph.D.Programs Foundation of Ministry of Education of China(20050141025)1.大连理工大学电子与信息工程学院电子工程系大连1160242.兰斯大学特鲁瓦技术学院通用电气与计算机工程系100261.Department of Electronic Engineering,School of Electronic and Information Engineering,Dalian University of Technology,Dalian 116024,P.R.China2.Department of Electrical Engineering and Computer Engineering,University of Reims,10026,France E reg(T)是惩罚项,旨在保证得到的变形场满足某些必要的属性,如平滑性、正反变换的一致性、拓扑保持性等,其中要求变形场具有拓扑保持性在完全基于配准的图像分割中是非常重要的,因为从医学上讲,不同正常个体的同源结构虽大小不同,但具有同样的拓扑.变形场具有拓扑保持性意味着各结构在变形后的图像中仍保持内部的连接性和结构之间的邻接关系不变,不会发生撕裂或粘贴,不会有新的结构出现,原有的结构也不会消失.通过将这一约束加在变形场上,可以限制式(1)的最优解空间,使它符合真实材料的变形特性.判断给定的变形场是否具有拓扑保持性,一般是考察该变形场的雅可比行列式的值是否总是大于零[4].在3D情况下,变形场向量可以表示成T=(X,Y,Z),X,Y,Z分别表示点p(x,y,z)变形后的位置,则变形场在p点处的雅可比行列式可以表示成J T(p)=∂X∂x∂X∂y∂X∂z∂Y∂x∂Y∂y∂Y∂z∂Z∂x∂Z∂y∂Z∂z(2)也就是说,只有保证在任意点处J T(p)>0,所得到的变形场才具有拓扑保持性[5−8].由于在变形场迭代寻优过程中约束雅可比行列式的值大于零的计算代价太大,文献[9]通过对正反变换进行处理以保证变形场的双射属性,然后再进行平滑处理,力求其具有拓扑保持性.然而,研究表明,很难从理论上证明利用双射和平滑技术可以保证变形场具有拓扑不变性[7,10],实践中,我们也通过实验发现该算法得到的变形场并不总是具有拓扑不变属性.为此,本文从矢量场的特点出发,详细分析了该算法得到的变形场拓扑性质,并在此基础上给出了强化变形场拓扑保持性的方法.1算法简介文献[9]提出的Demons算法应用于医学图像领域时,常基于光流场理论求解其变形场.光流场理论的假设前提是图像在运动的过程中保持灰度不变,即I(x(t),y(t),z(t),t)=const(3) I表示图像的灰度,是像素点坐标的函数.对式(3)求微分得∂I∂x∂x∂t+∂I∂y∂y∂t+∂I∂z∂z∂t=−∂I∂t(4)若将待配准的两幅图像看成是单位时间间隔内同一目标的成像结果,则上式变为v·∇s=r−s(5)其中v=(∂x∂t,∂y∂t,∂z∂t),r和s分别表示图像R和S的灰度.经过近似处理,得到运动速度的表达式为v=(r−s)∇s[||∇s||2+(r−s)2](6)这个速度代表了变形网格点受到的作用力.其寻优过程如下:1)初始化变换矩阵T;180自动化学报36卷2)相似度函数寻优,得到源图像中每个点p 的位移矢量:u (p )=−(r ◦T (p )−s (p ))(r ◦T (p )−s (p ))2+||∇s (p )||2∇s (p )(7)3)对得到的位移矢量进行高斯平滑处理:u ←G σ∗u 4)T ←T +u .重复迭代,直到满足停止条件.2变形场分析图像配准得到的变形场是一个矢量场,它反映了两图像中相对应结构间的变化.为了更好地理解变形场所模拟的物理过程和其中隐含的拓扑信息,有必要对得到的矢量变形场进行特性分析.矢量场的特性可以用它的特征点来表征[11],其中最重要的几个特征点是源点、汇点和漩涡中心点,如图1所示.图1矢量场的主要特征点Fig.1The important attribute points of a vector field根据式(7)可知,位移矢量是对标量场求梯度得到的,而梯度场是无旋场,所以本文所用配准算法的特征点主要是源点和汇点.矢量场的特征点代表了与应用相关的不同的物理属性.以不同个体间脑部核磁共振图像的配准应用为例,源点表征了对应脑部器官的增大,而汇点则表示对应脑部器官的缩小,它与变形场雅可比行列式的值所表征的物理特性是相吻合的.J T >1表征了对应脑部器官的增大,J T <1表征了对应脑部器官的缩小,J T=1表征了物质既没有增加也没有减少,J T ≤0表示存在拓扑改变[7−8,12].如前所述,图像变形带来的拓扑改变主要表现在图像中的目标结构发生粘贴或撕裂,参考矢量场特征点的图示(图1),可见汇点处容易发生挤压粘贴和反折,源点处容易因膨胀过度而撕裂.2D 情况下常见的几种情况如图2所示.(a)反折粘贴(a)Fold(b)挤压粘贴(b)Overlap(c)撕裂(c)Tear图22D 空间中变形导致的拓扑改变示例Fig.2Examples of the changed topology in 2D space due tothe deformations图2中的编号1、2、3、4表示初始状态下一个矩形网格的四个顶点,粗实线表示变形后的网格,其顶点为1 、2 、3 、4 ,分别与初始状态下的四个顶点相对应,带箭头的虚线表示各顶点的变化路径.图2(a)表示网格变形后出现了反折粘贴,就像是一个平面不仅形状发生了变化,而且由正面翻扣过来变成了背面;图2(b)表示网格变形后出现了挤压粘贴,2、3两点在变形过程中有交叉现象,导致网格扭曲;图2(c)表示网格变形后向两侧膨胀的幅度过大,超过了两点间连线允许的拉伸强度而撕裂.5 和6 表示网格断裂后插入的新顶点.研究表明,Demons 算法得到的变形场有拓扑改变现象,在某些区域得到的雅可比行列式的值不为正.观察这些区域的变形场特点,可见到与图2示例相近的情况,如图3所示.图3(a)为可变的模板图像,图3(b)为目标图像,图3(c)为模板图像在Demons 变形场作用下变形得到的图像,可以看出,二者具有较高的相似度.上面的白斑表示变形场雅可比行列式的值不大于零的位置,意味着得到该结果的变形场在某些位置不具有拓扑保持性.我们选择了三个拓扑发生改变的区域1,2,3,在图3(c)中用白色矩形框标出.图3(d)∼(f)给出了这三个区域的变形场图示,可以看到在这些矢量场的特征点附近有交叉或反折现象.(a)模板图像(a)Template(b)目标图像(b)Target(c)变形后的模板图像(c)Deformed template(d)1区放大(d)Enlarged region 1(e)2区放大(e)Enlarged region 2(f)3区放大(f)Enlarged region 3图3原Demons 算法得到的变形场特点Fig.3The deformation field obtained from the originalDemons algorithm1期林相波等:Demons非刚性配准算法拓扑保持性的研究181从图2可以看出,变形时发生反折主要表现在沿变形方向上相邻点的位移幅度不匹配,后面的点位移后的位置超过前面的点位移后的位置.变形时发生交叉主要是在矢量场汇点附近,相邻的点位移方向相对,各点的位移幅度相对较大导致相互交叉.因此,既直观又简单的解决方案应是在尽量保持原有变形场几何特性的前提下合理地改变这些点处的位移,达到拓扑校正的目的.3增强变形场拓扑保持属性的方法若连续可微的变形场具有拓扑保持性,则其雅可比行列式的值在其定义域内应处处为正.数字图像受分辨率的限制,是定义在离散网格点上的,不是一个连续函数,因而得到的变形场也是离散的.文献[10]指出,连续变形场是离散变形场经过插值得到的,它的特性决定于离散场的特性.如果离散变形场在各个离散网格点处雅可比行列式的值为正,则经过线性插值后得到的连续变形场在定义域内具有拓扑保持性.设连续变形场T =(X,Y,Z )可以表示成X (x,y,z )=x +u x (x,y,z )(8)Y (x,y,z )=y +u y (x,y,z )(9)Z (x,y,z )=z +u z (x,y,z )(10)其中u x (x,y,z )、u y (x,y,z )和u z (x,y,z )分别表示点p (x,y,z )在x,y,z 坐标方向上的位移,该点处雅可比行列式的值为J (x,y,z ).定义另外一个变形场T k ,其表达式为X k (x,y,z )=x +ku x (x,y,z )(11)Y k (x,y,z )=y +ku y (x,y,z )(12)Z k (x,y,z )=z +ku z (x,y,z )(13)该变形场在p (x,y,z )点处的雅可比行列式的值为J k (x,y,z ).显然,当k =1时,J k (x,y,z )=J (x,y,z );若p (x,y,z )点处J (x,y,z )<0,意味着其不具有拓扑保持性.J k (x,y,z )是k 的连续函数,根据函数的连续性可知,至少存在一个k ∗∈[0,1],使得J k ∗(x,y,z )>0.离散场是对连续场的采样,则对于离散场,通过在原有位移场上乘以一个0∼1之间的系数,可以校正原离散变形场的拓扑性.变形场的雅可比行列式建立在对变形场求梯度的基础上.在离散条件下,我们采用中间差分方式求梯度,任一点的梯度在不同坐标方向上的分量由该坐标方向上待求点两侧的点(可以称这些点为关联点)确定.具体校正过程如下:步骤1.根据得到的变形场T 计算各点的雅可比行列式值J (x,y,z );步骤2.若在点(x,y,z )处J (x,y,z )<0,则从k =1开始,等步长减小k 值,得到点(x,y,z )的关联点的新位移.对每一个k ,根据这些新位移计算J k (x,y,z ).重复这个过程,直到找到k ∗,使所有的J k ∗(x,y,z )>0;步骤3.根据得到的T k ∗,检验变形后的模板图像与目标图像的匹配程度,若不满足要求,则以新得到的变形后的模板图像为可变形模板,重新配准,对得到的T 重复这一过程,直到满足条件为止.4实验本文实验采用真实的脑部医学核磁共振图像.其中模板图像为SPL [13]标准图谱,T1加权,包含256×256×160个体素.目标图像的采集仪器为1.5T GE 核磁共振仪,T1加权,256×256×124个体素.实验目的主要是验证本文提出的增强变形矢量场拓扑保持性方法的性能,主要进行了两方面的对比,分别比较了变形场拓扑校正前后其矢量空间特性的变化和不同的目标结构分割结果.图4所示为原Demons 变形场经过拓扑校正处理后的结果.实验中改变k 值时步长设为0.01,迭代终止条件为变形后的模板图像与目标图像的互相关系数达到给定的阈值或者以新得到的图像为可变形模板重新配准次数到达上限(本实验中设为3次).本文实验的互相关系数阈值设置方式为:设模板图像与目标图像的初始互相关系数为CC 0,则阈值CC t =(1−CC 0)/α+CC 0,参数α的值设为1.2.分别比较图4(a)∼(d)与图3(c)∼(f),可以看出校正后的变形场雅可比行列式的值满足拓扑保持性的要求,而且校正后的变形场更加规则、平滑.从图3(c)给出的实验结果可以看出,原变形场雅可比行列式的值不大于零的点所占比例较小,也就是说需要校正的点所占比例较小,因此相对于原算法其计算量只大约增加了5%.(a)变形后的模板图像(a)Deformed template(b)1区放大(b)Enlarged region 1(c)2区放大(c)Enlarged region 2(d)3区放大(d)Enlarged region 2图4拓扑校正后的变形场Fig.4The deformation field with corrected topologies图5所示为变形场拓扑校正后,各点位移量的变化值,颜色越浅表示该位置处的位移矢量的幅度改变得越多,最大改变量为1.86个像素,说明校正后的变形场基本保持了原有182自动化学报36卷场的几何特性.同时,图5表明位移幅值发生明显变化的像素点的个数相对较少,也证明了采用本文的拓扑校正方法不会使计算量大幅增加.图6比较了采用拓扑校正前和校正后的变形场分割丘脑结构的效果.图6(a)为可变形模板图像,上面叠加了丘脑的专家分割结果.图6(b)为待分割的目标图像,图6(c)中的白斑给出了校正前的变形场拓扑发生改变的区域,采用该变形场得到的丘脑分割结果见图6(d),显然得到的结果有撕裂现象,丘脑的拓扑结构发生了改变,相应位置的变形场有交叉现象,见图6(e).经过本文提出的拓扑校正处理后,表征雅可比行列式的值不为正的白斑消失,见图6(f),说明此时的变形场具有拓扑保持性,用该变形场进行同样的分割,避免了丘脑结构的撕裂,如图6(g)所示.对比校正前后该位置的变形场(图6(e)和图6(h)),可以看出本文提出的方法消除了其中的矢量交叉,使变形场更平滑、更合理.不过,从图6(g)可以看出,虽然校正后的变形场具有了拓扑保持性,但其分割结果仍存在误分割处.这是因为原配准算法是单纯基于图像灰度信息的,当邻近的不同目标的灰度值相似时,得到的变形场可能是局部最优解,而本算法在对变形场进行拓扑校正时,遵循尽可能少地改变原有变形场几何特性这一原则,一旦满足拓扑保持条件,即停止校正,所以不能完全改变这一局部最优变形场的情况.可以考虑增加配准时所用的信息(如位置,形状等),尽量避免初始变形场陷入局部最优.图5变形场校正前后的幅度差(以像素值度量)Fig.5The magnitude differences of the deformation fields before and after topology correction (measured in pixels)5结论在非刚性图像配准的研究中,保证变形场的拓扑不变性是一个重要问题.研究表明,经典的基于灰度的Demons 配准算法虽然采用了双射技术和平滑处理,所得变形场仍然存在拓扑改变问题.本文从矢量场几何特性的角度,分析了该变形场导致拓扑改变的原因,提出在尽量保持原有变形场几何特性的前提下直接调整点位移的方法校正原变形场的拓扑性,并给出了具体的理论支持.实验表明,本文给出的方法能够校正原变形场在拓扑保持方面的缺陷,得到更平滑、合理且具有拓扑保持性的变形场.(a)模板图像(a)Template(b)目标图像(b)Target(c)校正前拓扑状况(c)Topology beforecorrection(d)校正前丘脑分割(d)Thalams segmentationbefore correction(e)校正前的场放大图(e)Enlarged field beforecorrection(f)校正后拓扑状况(f)Topology aftercorrection(g)校正后丘脑分割(g)Thalams segmentationafter correction(h)校正后的场放大图(h)Enlarged field beforecorrection图6变形场校正前后脑结构分割对比Fig.6Comparison of brain structures segmentation based on the deformation fields with and without topology corrections1期林相波等:Demons非刚性配准算法拓扑保持性的研究183References1Gholipou A,Kehtarnavaz N,Briggs R,Devous M,Gopinath K.Brain functional localization:a survey of image regis-tration techniques.IEEE Transactions on Medical Imaging, 2007,26(4):427−4512Cachier P,Bardinet E,Dormont D,Pennec X,Ayache N.Iconic feature based nonrigid registration:the PASHA puter Vision and Image Understanding,2003, 89(2-3):272−2983Shen D G,Davatzikos C.HAMMER:hierarchical attribute matching mechanism for elastic registration.IEEE Transac-tions on Medical Imaging,2002,21(11):1421−14394Holden M.A review of geometric transformations for non-rigid body registration.IEEE Transactions on Medical Imag-ing,2008,27(1):111−1285Noblet V,Heinrich C,Heitz F,Armspach J P.Retrospec-tive evaluation of a topology preserving non-rigid registra-tion method.Medical Image Analysis,2006,10(3):366−3846Musse O,Heitz F,Armspach J P.Topology preserving deformable image matching using constrained hierarchical parametric models.IEEE Transactions on Image Processing, 2001,10(7):1081−10937Noblet V,Heinrich C,Heitz F,Armspach J P.3-D de-formable image registration:a topology preservation scheme based on hierarchical deformation models and interval anal-ysis optimization.IEEE Transactions on Image Processing, 2005,14(5):553−5668Leow A D,Yanovsky I,Chiang M C,Lee A D,Klunder A D, Lu A.Statistical properties of jacobian maps and the realiza-tion of unbiased large-deformation nonlinear image registra-tion.IEEE Transactions on Medical Imaging,2007,26(6): 822−8329Thirion J P.Image matching as a diffusion process:an anal-ogy with Maxwell s demons.Medical Image Analysis,1998, 2(3):243−26010Karacah B,Davatzikos C.Estimating topology preserving and smooth displacementfields.IEEE Transactions on Med-ical Imaging,2004,23(7):868−88011Tittgemeyer M,Wollny G,Kruggel F.Visualising deforma-tionfields computed by non-linear image -puting and Visualization in Science,2002,5(1):45−5112Rohlfing T,Maurer C R J,Bluemke D A,Jacobs M A.Volume-preserving nonrigid registration of mr breast images using free-form deformation with an incompressibility con-straint.IEEE Transactions on Medical Imaging,2003,22(6): 730−74113Kikinis R,Shenton M E,Iosifescu D V,McCarley R W, Saiviroonporn P,Hokama H H.A digital brain atlas for surgical planning,model-driven segmentation,and teaching.IEEE Transactions on Visualization and Computer Graph-ics,1996,2(3):232−241林相波大连理工大学电子与信息工程学院博士研究生.主要研究方向为生物医学图像与信号处理,模式识别.E-mail:linxbo@ (LIN Xiang-Bo Ph.D.candidate in the School of Electronic and Information Engineering,Dalian University of Technology. Her research interest covers biomedical image and signal pro-cessing,pattern recognition.)邱天爽大连理工大学电子与信息工程学院教授.主要研究方向为非高斯、非平稳信号分析与处理理论,生物医学信号与图像处理,通信信号处理理论.本文通信作者.E-mail:qiutsh@ (QIU Tian-Shuang Professor in the School of Electronic and Information Engineering,Dalian University of Technology. His research interest covers non-Gaussian,non-stable signal pro-cessing theory,biomedical signal and image processing,commu-nication signal processing theory.Corresponding author of this paper.)阮素法国兰斯大学通用电气与计算机工程系教授.主要研究方向为图像分割与分类,数据融合,模式识别,医学图像处理.E-mail:su.ruan@univ-reims.fr(RUAN Su Professor in the Department of Electrical Engi-neering and Computer Engineering,University of Reims,France. Her research interest covers image segmentation and classifica-tion,data fusion and pattern recognition,medical image pro-cessing.)NICOLIER Frédéric法国兰斯大学通用电气与计算机工程系副教授.主要研究方向为图像比较,分割评价.E-mail:frederic.nicolier@univ-reims.fr(NICOLIER Frédéric Associate professor in the Depart-ment of Electrical Engineering and Computer Engineering,Uni-versity of Reims,France.His research interest covers image com-parison by similarity measurements,segmentation evaluation.)。

第42卷第9期自动化学报Vol.42,No.9 2016年9月ACTA AUTOMATICA SINICA September,2016基于平衡系数的Active Demons非刚性配准算法薛鹏1杨佩1曹祝楼1贾大宇1董恩清1摘要经典的Active demons算法利用参考图像和浮动图像的梯度信息作为驱动力,并使用均化系数调节两种驱动力之间的强度.该算法克服了Demons算法单一使用参考图像的梯度信息作为驱动力的缺点,但是Active demons算法中的均化系数无法同时兼顾大形变和小形变区域的准确配准,还会导致配准的收敛速度和精确度相互制约的问题.为此,本文提出一种新的Active demons非刚性配准算法.提出的算法在Active demons扩散方程中引入一个称为平衡系数的新参数,与均化系数联合调整驱动力,不仅可以兼顾图像中同时具有的大形变和小形变区域的准确配准,而且在一定程度上缓和了收敛速度和精确度相互制约的问题.为了进一步提高配准的收敛速度和精确度,避免陷入局部极值,在新的配准算法的实现中引入由粗到细的多分辨率策略.在Checkboard测试图像、自然图像和医学图像上的实验结果表明,提出的算法较经典的Active demons算法收敛速度更快,配准精度平均提高了54.28%,接近最新的TV-L1光流场图像配准算法的配准精度,解决了Active demons 算法存在的问题.关键词非刚性配准,Active demons算法,光流场图像配准,驱动力,多分辨率策略引用格式薛鹏,杨佩,曹祝楼,贾大宇,董恩清.基于平衡系数的Active demons非刚性配准算法.自动化学报,2016,42(9): 1389−1400DOI10.16383/j.aas.2016.c150186Active Demons Non-rigid Registration Algorithm Based on Balance Coefficient XUE Peng1YANG Pei1CAO Zhu-Lou1JIA Da-Yu1DONG En-Qing1Abstract 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 by54.28% on average,and is close to that of the latest TV-L1opticalflow image registration algorithm.Key words Non-rigid image registration,active demons algorithm,opticalflow image registration,driving force,multi-resolution strategyCitation Xue Peng,Yang Pei,Cao Zhu-Lou,Jia Da-Yu,Dong En-Qing.Active demons non-rigid registration algorithm based on balance coefficient.Acta Automatica Sinica,2016,42(9):1389−1400图像配准是图像处理和计算机视觉等领域的热收稿日期2015-04-22录用日期2016-03-20Manuscript received April22,2015;accepted March20,2016国家自然科学基金(81371635),高等学校博士学科点专项科研基金(20120131110062),山东省科技发展计划项目(2013GGX10104) Supported by National Natural Science Foundation of China (81371635),Specialized Research Fund for the Doctoral Pro-gram of Higher Education of China(20120131110062),Science and Technology Development Project of Shandong Province (2013GGX10104)本文责任编委杨健Recommended by Associate Editor YANG Jian1.山东大学(威海)机电与信息工程学院威海2642091.School of Mechanical,Electrical&Information Engineering, Shandong University,Weihai264209点研究问题,广泛应用于遥感图像拼接、目标定位和疾病诊断治疗等方面.一般我们将配准过程中固定不变的图像称为参考图像,将施加变换的图像称为浮动图像.图像配准通过进行参考图像和浮动图像的空间变换,最终使得参考图像和浮动图像对齐.一般来说,图像配准在不同领域存在一定差异,比如:将多幅遥感图像拼接成一幅更大的图像,多幅图像之间有重叠部分,重叠部分的特征提取是配准的关键;而在医学图像配准中,两幅图像之间几乎是完全重叠的,研究关注的核心问题是图像的非均匀形变.1390自动化学报42卷按照配准过程中图像是否存在形变来划分,图像配准分为刚性配准和非刚性配准.其中,非刚性图像配准是近年来的研究热点之一[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根据李群理论提出了Diffeomor-phic demons配准算法.基于该理论,后人提出了很多的改进算法.例如,Lorenzi等[16]提出在形变场中增加SVFs(Stationary velocityfields)约束和使用新的相似性测度函数LCC(Local correlation co-efficient).另外,林相波等[17]从矢量场的特点出发,详细分析了该算法得到的形变场拓扑性质,并在此基础上给出了强化形变场拓扑保持性的方法.策略2是在驱动力计算公式中增加新的信息项,弥补仅依靠灰度梯度信息配准图像的不足;常见的信息项有灰度梯度场的相似性[18]、对称正交梯度信息[19]、几何形状约束[20]等.策略3从分力的影响程度上开展研究,如Rogelj等[12]提出对称梯度平均化的Symmetric demons算法;Lin等[21]对两个分力进行了取舍,如果某个分力指向图像匹配方向,则对其进行加强,否则将其设为零,取消该分力的影响.文献[21]中的算法有利于提高配准的精度,但是分力的取舍需要较长的时间.为了加快配准的收敛速度,Sharp等[22]研究了Demons算法基于硬件的加速方法,利用GPU中的并行计算能力,大幅减少算法运行时间.也有一些学者采用多分辨率策略[21,23],将低分辨率层的空间变换作为高分辨率层的初始形变矩阵,按照由粗到细的方式进行逐层配准.文献[24]提出一种在配准迭代过程中采用自适应地调节均化系数的算法.Pock等[25]及S´a nchez等[26]在传统光流模型的基础上,利用全变差正则项(Total variation reg-ularization,TV正则项)和鲁邦性数据项(L1范数项)构建能量函数得到的TV-L1光流场图像配准算法的配准精度较高,与Demons类算法相比运算时间缩短不小,但是同样不能解决大形变的问题.本文在Active demons算法的基础上,提出了一种能够同时提高配准精度和收敛速度的非刚性配准算法.主要思路如下:1)在Active demons算法的驱动力计算公式中引入平衡系数k,该系数联合均化系数α精细地调整驱动力的强度,使得Active demons算法适合形变范围更大的图像配准;2)通过实验分析,给出平衡系数k的一般取值范围,增强算法的普适性;3)引入了多分辨率策略.在实验测试分析中,采用配准精度和收敛速度指标,将提出的算9期薛鹏等:基于平衡系数的Active Demons非刚性配准算法1391法分别与经典的Active demons配准算法和TV-L1光流场图像配准算法进行对比.本文剩余部分结构安排如下:在第1节,本文回顾经典的Demons算法以及Active demons算法;在第2节,介绍本文提出的改进Active demons算法;在第3节,通过对Checkboard测试图像、自然图像和医学图像的实验验证本文算法的有效性;第4节对本文工作进行总结.1方法与原理1.1Demons算法Demons算法[3]将参考图像和浮动图像看成是连续运动图像序列中的两帧,并估计从浮动图像到参考图像的形变向量u.通常该向量理解为施加在浮动图像上的形变驱动力.假设图像在运动的过程中保持灰度不变,对于空间中任意一点(x,y),如果参考图像S和浮动图像M上的灰度分别为S(x,y)和M(x,y),∇S(x,y)为参考图像在(x,y)处的梯度值,根据光流场方程可以得到下式:u∇S(x,y)=M(x,y)−S(x,y)(1)则u可以用下式估计:u=M(x,y)−S(x,y)|∇S(x,y)|2∇S(x,y)(2)显然,当|∇S(x,y)|2很小时,u可能无限大,与实际情况不符.为了解决该问题,在式(2)分母上增加一新项,变为u=M(x,y)−S(x,y)|∇S(x,y)|2+(M(x,y)−S(x,y)2)∇S(x,y)(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.2Active demons算法在Demons算法中,使浮动图像发生形变的力单纯取自参考图像的梯度信息.这样的模型只适合于处理小形变问题,不能满足较大形变图像配准应用的需要.为了扩大该算法的应用范围,Wang等[13]根据牛顿第三定律的作用力与反作用力的原理,提出了Active demons算法,将浮动图像的梯度信息也作为一种正内力引入到驱动力计算公式中,得到下式:u=M(x,y)−S(x,y)|∇S(x,y)|2+(M(x,y)−S(x,y))2∇S(x,y)+M(x,y)−S(x,y)|∇M(x,y)|2+(M(x,y)−S(x,y))2∇M(x,y)(4)为了能够调整驱动力的强度,在上式中引入了一个均化系数α,式(4)变为u=M(x,y)−S(x,y)|∇S(x,y)|2+α2(M(x,y)−S(x,y))2∇S(x,y)+M(x,y)−S(x,y)|∇M(x,y)|2+α2(M(x,y)−S(x,y))2∇M(x,y)(5)文献[27]详细分析了α的取值对配准结果的影响,该分析结果表明:较小的α可以加快算法的收敛速度,减少配准时间,但配准精度较差;较大的α会增加配准时间,但得到的配准结果较好.2基于Active demons算法的非刚性配准方法2.1一种改进的Active demons算法本文针对Active demons算法在图像配准中存在的问题,提出了一种能够同时提高配准精度和收敛速度的改进算法.在Active demons扩散方程中,为了能够调节驱动力,Wang在分母中灰度差的平方项上引入了一个均化系数α.受此启发,将该思想应用在另一分母项上,即在梯度的模值的平方项上引入一个新的系数k,我们称之为平衡系数.然后通过调整均化系数α和平衡系数k的取值,联合调节驱动力的强度,较依靠单一系数能够更加精细地控制配准过程中的的形变程度.改进后的Active demons扩散方程如下所示:u=M(x,y)−S(x,y)k2|∇S(x,y)|2+α2(M(x,y)−S(x,y))2∇S(x,y)+M(x,y)−S(x,y)k2|∇M(x,y)|2+α2(M(x,y)−S(x,y))2∇M(x,y)(6)式(6)中,如果均化系数α不变,平衡系数k对驱动力u的大致影响可以分为三种情况:当k<1时,分母项减小,u变大,形变程度增大;当k>1时,分母项增大,u变小,形变程度减小;当k=1时,式(6)等价于Active demons扩散方程的形式.下一节将具体分析平衡系数k的取值对配准结果的影响.1392自动化学报42卷2.2平衡系数的引入及取值分析2.2.1引入平衡系数的有效性分析为了清晰地了解参数k对配准结果的影响,本节采用简单的二值图像进行配准.如图1所示,正方形为浮动图像,圆为参考图像.实验参数设置如下: 1)图像分辨率为256×256;2)Demons估计浮动图像中所有点的形变大小;3)初始形变矩阵为零,空间变换为自由项变换;4)最大迭代次数为200,单尺度实现;5)图像插值使用三线性插值法.图1二值图像Fig.1The binary images为了进行详细的分析,我们先使用经典的Ac-tive demons算法进行配准实验,均化系数α取10个不同的值,依次为0.05,0.1,0.4,0.5,0.6,1,1.5, 2,2.5,3,得到迭代过程的均方差曲线图,如图2所示.图2表明:随着迭代次数的增加,均方差基本上保持单调递减,收敛后保持微幅振动.α越小,配准的收敛速度越快,α较大时收敛变慢,但是配准精度更高,这也与文献[27]的分析相符.图2均化系数α对Active demons算法的影响Fig.2The impact ofαon the active demons algorithm实验条件不变,采用本文提出的改进Active demons算法进行配准实验.平衡系数k取0.5,实验结果的均方差与迭代次数的关系曲线图如图3所示.为了比较两种算法的配准精度,表1为两种算法迭代200次时得到的均方差(Mean square error,MSE)值.AD为Active demons算法,IAD(Im-proved active demons)为引入平衡系数k之后的Active demons算法.可以很明显地看到,IAD算法的均方差值均小于AD算法,这说明改进算法配准精度更高.为了考察收敛速度,分别从两个实验中取出4组具有代表性的实验数据进行对比分析.图4是两种算法的配准曲线对比图,其中,α取为0.05, 0.4,1,2,图4中带标志的曲线表示Active demons 算法,未带标志的曲线表示改进算法,相同的线型代表相同的均化系数.从曲线图中可以明显看到,改进后的Active demons算法的配准曲线与经典的Active demons算法相比更贴近y轴,这说明本文算法的收敛速度更快.表1两种算法的配准结果均方差对比Table1The comparison of two registration algorithmson the MSEα0.050.10.40.50.6 AD(×10−4) 3.55 4.02 6.498.1810IAD(×10−4) 2.25 2.35 3.71 4.29 5.27α 1.0 1.52 2.53 AD(×10−4)34107187255311IAD(×10−4)1247113183246图3均方误差与迭代次数的关系曲线Fig.3The relations between the mean square error anditerations以上分析表明:本文在Active demons算法中引入平衡系数k是非常正确的,不仅提高配准的精确度,收敛速度也更快.那么,对于固定的某一均化系数α,该如何确定平衡系数k的值才可使配准性能达到最好呢?接下来,本文将具体分析这个问题.2.2.2最优平衡系数的取值分析从理论上分析,在式(6)中,因为|k∇S(x,y)|2+α2(M(x,y)−S(x,y))2≥2αk∇S(x,y)(M(x,y)−S(x,y)),所以驱动力的上限为u≤(1/αk).若α不9期薛鹏等:基于平衡系数的Active Demons非刚性配准算法1393变,u与k成反比.k越大,u就越小,形变程度小,配准的精度高,但收敛速度较慢.k越小,u越大,允许的形变度越大,收敛速度很快,但配准精度较低.图4两种算法对比曲线Fig.4The comparison of two registration algorithms在实验中,设定α为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越大,收敛速度变慢,但是配准精度更高.这也与上面的理论分析相符.图8为不同α和k组合情况下,以固定迭代200次之后的配准结果的均方差关系曲面图,其最小值点位于α=0.5,k=0.5处.从实验结果中可以看出,α和k较小时,配准精度高平滑性较差;α和k较大时,配准精度相对低.建议二者的取值在[0.5,1.5]之间较适宜.图7平衡系数对改进Active demons算法的影响Fig.7The impact of balance coefficient k on theimproved active demons registration algorithm2.3基于多分辨率策略实现Active demons算法容易陷入局部极小值,而且该算法的运行速度慢、配准精度不高.多分辨率策略可以避免陷入局部极值,在提高算法的准确性、速度和鲁棒性方面都很有帮助,其基本思路如下:1)通过降采样的方法将待配准图像分解成多种不同分辨率的图像;2)在低分辨率图像上快速进行粗配准;图5不同平衡系数k的配准结果Fig.5The registration results with different k variants图6不同参数k的配准结果差值Fig.6The registration results error with different k variants1394自动化学报42卷3)对在低分辨率图像上得到的形变场进行升采样,将其作为更高一级分辨率的初始变换;4)逐层的迭代,最终实现整幅图像的配准.为了进一步提高配准的收敛速度、准确度以及避免陷入局部极值,本文基于多分辨率策略实现,具体实现步骤如下所示.图8均方差与均化系数和平衡系数的关系曲线Fig.8The relations of the mean square error withtheαand the k算法1.基于多分辨率策略实现的配准过程1:设定初始形变矩阵T0,最大迭代次数k,最大分解尺度L,当前分解级数l,当前迭代次数n;2:对参考图像S和浮动图像M进行L级尺度分解;3:while(l<L)do4:while(n<N)do5:根据式(6)计算第n次迭代后浮动图像M的位移向量u l n,并更新形变矩阵T l=T l−1+u l n,将T l应用于浮动图像M使其形变;6:n=n+17:end while8:l=l+1,n=09:对前一级图像的形变矩阵T l−1升采样,作为当前图像配准的初始形变矩阵;10:end while11:将最终的形变矩阵T应用于浮动图像M,应用双三次插值得到配准结果.3实验结果分析本文先用主观评价方法对所得配准结果进行分析,并获得初步的质量评价结果,再用5种客观评价方法进一步准确地分析配准效果.为了分析和验证本文提出算法的实际效果,分别采用测试图像、自然图像和医学图像进行实验.实验中同时采用了Demons算法、Active demons算法(AD)、单尺度实现的改进算法(Single-scale im-plementation active demons,SIAD)、多分辨率实现的改进算法(Multi-scale implementation active demons,MIAD)和TV-L1光流场算法(TV-L1).为了公平起见,在下面的实验中,各种算法的参数选择都是相对于该算法的最佳结果.3.1Checkboard测试图像对比分析图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种算法的配准结果的客观比较.通过对比图9和图10,可以发现,Demons算法由于仅依靠参考图像的梯度信息配准图像,配准结果与参考图像的误差最大,配准效果最差;Active demons算法的配准结果相比Demons算法有了一定的改善,尤其是对Demons算法无法配准的大形变区域实现了较好的配准;本文提出的MIAD算法的配准误差最小,效果最好;SIAD算法的配准结果次之,但也好于原来的两种配准算法.虽然TV-L1算法配准的效果最好,但其与参考图像的差值明显有不规则的边缘效应.9期薛鹏等:基于平衡系数的Active Demons非刚性配准算法1395图9图像的配准结果图Fig.9The image registration results图10图像的配准结果与参考图像的差值Fig.10The differences between the registration results and the static image从表2的5种客观评价方法中可以看出,MIAD 算法的均方差值最小,相关系数最大,峰值信噪比最大,归一化互信息最大,结构相似度也是最大的. SIAD算法的5种评价值也均优于Demons算法和Active demons算法,这也与前面的主观分析一致.综合以上分析说明,对于测试图像,无论从配准图像的主观效果和客观指标,显然目前最新的TV-L1光流场算法是最好的,但本文提出的SIAD算法、MIAD算法的配准性能与之非常接近,且好于原来的同类算法.表2配准结果的客观分析Table2The objective analysis of registration results 评价方法Demons AD SIAD MIAD TV-L1均方差(×10−4)4340292511相互系数(%)99.1699.2499.4499.5299.98峰值信噪比54.4355.3158.3659.9787.74归一化互信息 1.42 1.42 1.43 1.430.71结构相似度(%)95.4295.7496.8297.4199.97 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种算法的配准结果的客观比较.从图11和图12的配准结果差值图中可以明显看出,Demons算法完成了图像中小形变部分的配准,在一定程度上纠正了浮动图像中面部五官的形变,但头发部位的形变加强了,眉毛和帽子也发生了畸变;Active demons算法相比Demons算法在处理大形变图像配准问题上优越性明显,尤其是头发部位的矫正.但对于形变较小的四周区域,在经过了Active demons算法配准之后,明显形变过大; SIAD算法不仅解决了Active demons算法存在的小形变区域过配准的问题,而且保留了大形变区域配准的准确性,配准误差更小,MIAD算法的配准误差最小,配准的差值图很光滑,证明了多分辨率策略实现的优越性.从配准图像及与参考图像的差值来看,TV-L1算法配准的效果是最好的,但其差值图边缘具有不规则的边缘效应.从表3中的5种客观评价方法可以看出,对于TV-L1光流场算法,除了均方误差指标不及本文提出的算法外,其他的指标都是最好的.同理,对于本文提出的算法,除了TV-L1外,其5个指标均优于其他同类算法,特别是结构相似度,比其他算法都高出很多.这也进一步表明了本文算法优越性,能够较好地保持图像的拓扑结构.图13是针对小形变的彩色自然图像进行1396自动化学报42卷的配准实验,浮动图像只出现了局部小形变.图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种算法的配准结果的客观指标比较进行评判.表3配准结果的客观分析Table3The objective analysis of registration results 评价方法Demons AD SIAD MIAD TV-L1均方差(×10−4)138.11 6.19 2.81 4.61相互系数(%)96.3597.6998.2599.2299.87峰值信噪比66.5971.1573.8781.7591.49归一化互信息 1.37 1.29 1.39 1.41 3.96结构相似度(%)88.5590.9992.5696.2899.11从表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算法的优势.表4配准结果的客观分析Table4The objective analysis of registration results 评价方法Demons AD SIAD MIAD TV-L1R 5.09 5.03 5.1218.69 2.50均方差(×10−5)G 5.07 5.10 5.3317.68 2.47B 4.69 4.79 5.0116.71 3.20R99.8499.9299.9299.9099.99相互系数(%)G99.8699.8699.8599.8399.98B99.9299.8499.8399.8199.71R91.0791.1291.0490.1994.16峰值信噪比G91.0891.0590.8690.4394.20B91.4291.3391.1490.6793.08R 3.25 3.26 3.25 3.27 3.98归一化互信息G 3.17 3.17 3.16 3.17 4.03B 3.02 3.02 3.00 3.01 3.63R98.8797.4997.5697.7199.21结构相似度(%)G98.9997.4697.4197.4499.27B98.7097.3097.2397.2398.89图11图像的配准结果图Fig.11The image registration results图12图像的配准结果与参考图像的差值Fig.12The differences between the registration results and the static image。