医学图像的配准与融合

第八章医学图像的配准与融合第一节概述-、医学图像配准与融合的应用背景随着计算机技术的飞速发展，与计算机技术密切相关的医学成像技术也是日新月异。但是，各种成像技术和检查方法都有它的优势与不足，并非一种成像技术可以适用于人体所有器官的检查和疾病诊断，也不是一种成像技术能取代另一种成像技术，它们之间是相辅相成、相互补充的。如CT和X线机对骨等密度较高的组织能提供高清晰的图像，MRI对人体软组织的成像具有较高的分辨率，而PET和SPECT则能够提供人体组织或器官的功能性代谢的图像。成像原理的不同造成了某一种成像模式所能提供的图像信息具有一定的局限性，有时单独使用某一类图像难以获得正确的诊断结论。因此，为了提高诊断正确率，需要综合利用患者的各种图像信息。图像配准与融合技术为医学图像的综合利用提供了很好的技术手段。根据医学图像所提供的信息，可将医学图像分为两大类：解剖结构图像(CT、MRkX线图像等)和功能图像(SPECT,PET等)。这两类图像各有其优缺点：解剖图像以较高的分辨率提供了脏器的解剖形态信息，但无法反映脏器的功能情况。功能图像分辨率较差，但它提供的脏器功能代谢信息是解剖图像所不能替代的，这些信息是对疾病特别是肿瘤进行早期诊断的重要依据。目前医学影像学的•个明显的发展趋势是利用信息融合技术，将多种医学图像结合起来，充分利用不同医学图像的特点，在一幅图像上同时表达来自人体的多方面信息，使人体内部的结构、功能等多方面的状况通过影像反映出来，从而更加直观地提供人体解剖、生理及病理等信息。二、医学图像配准与融合的关系医学图像配准和融合有着密切的关系，特别是对多模态图像而言，配准和融合是密不可分的。待融合的图像往往来自于不同的成像设备，它们的成像方位、角度和分辨率等因子都是不同的，所以这些图像中相应组织的位置、大小等都有差异，若事先不对融合图像进行空间上的对准，那么融合后的图像豪无意义。因此，图像配准是图像融合的先决条件，必须先进行配准变换，才能实现准确地融合。三、医学图像配准和融合在临床中的应用医学图像配准和融合具有很重要的临床应用价值。对使用各种不同或相同的成像手段所获得的医学图像进行配准和融合不仅可以用于医疗诊断，还可用于外科手术计划的制定、放射治疗计划的制定、病理变化的跟踪和治疗效果的评价等各个方面。1、在外科手术中的应用了解病变与周围组织的关系对制定手术方案，决定手术是否成功至关重要。如对脑肿瘤患者，一般是采用外科手术切除肿瘤。患者的生存时间和生活质量与病灶（如肿瘤、血肿等）的切除程度密切相关。如果对病灶过度切除，会造成对病灶周围重要功能区域的损害，而这种损害是不可逆转的，严重影响患者的生活质量；反之如果对病灶切除不够，残余病灶会严重影响患者的生存时间。最大程度地切除病灶，同时使主要的脑功能区域（如视觉、语言和感知运动皮层等）得以保留是神经外科手术的目标。为此，在手术前，•般要利用CT或MRI获取患者的脑肿瘤结构信息，利用PET或fMRI获取患者脑肿瘤周围的脑功能信息，通过对结构成像和功能成像的配准、融合，对脑肿瘤及其周围的功能区进行精确定位，在此基础上制定出外科手术计划，是对患者进行精确手术的基础。2,在放射治疗中的应用大约70%的病人在肿瘤的治疗过程中接受放疗。放疗的目的就是最大限度的把放射能量集中在靶位上，从而使周围的正常组织的损害达到最小。放射治疗中，应用CT和MR图像的配准和融合来制定放疗计划和进行评估，用CT图像精确计算放射剂量，用MR图像描述肿瘤的结构。用PET和SPECT图像对肿瘤的代谢、免疫及其他生理方面进行识别和特性化处理，整合的图像可用于改进放射治疗计划或立体定向活检或手术。此外，放射治疗后扫描的MRI图像中，坏死组织往往表现为亮区，很容易与癌症复发混淆。把MRI图像与PET或SPECT图像进行配准，可区分坏死组织（没有代谢）与肿瘤复发（通常表现为高代谢）。3、在癫痫病治疗中的应用原发癫痫病灶的准确定位一直是困扰影像界的一大难题，许多学者利用配准和融合技术对此做出了富有成效的探索。例如：Pelizzari等人对癫痫病人的MREPET图像融合处理后，可观察到病人的脑外伤、炎症、硬化症等的变化,还可看到手术及麻醉前后的区别；Lewis等研究表明，在发作期和发作期间，对癫痫患者分别进行SPECT检查，将两者的图像相减，再分别与MRI图像融合，可使功能损伤的解剖学标记更准确，以SPECT所示的局部脑血流对大脑新皮质的癫痫病灶进行准确定位，从而为手术提供重要依据。第二节医学图像配准技术-、医学图像配准的概念对几幅不同的图像作定量分析，首先要解决这几幅图像的严格对齐问题，这就是我们所说的图像配准（imageregistration）o医学图像配准是指对于一幅医学图像寻求一种（或一系列）空间变换，使它与另一幅医学图像上的对应点达到空间上的•致。这种一致是指人体上的同一解剖点在两张匹配图像上有相同的空间位置（位置一致，角度一致、大小一致）。配准的结果应使两幅图像上所有的解剖点，或至少是所有具有诊断意义的点及手术感兴趣的点都达到匹配。图8-1是配准示意图。同一个人从不同角度、不同位置拍摄的两张照片，由于拍摄条件不同，每张照片只反映某些方面的特征。要将这两张照片一起分析,就要将其中的一张人像做移动和旋转，使它与另一幅对齐。这一对齐过程就是配准过程。保持不动的叫做参考图像，做变换的称作浮动图像。将配准后的图像进行融合就可以得到反映人的全貌的融合图像。图8-1图像配准示意图二、医学图像配准方法的分类到目前为止，医学图像配准方法的分类始终没有一个统一的说法。目前比较

流行的是1993年VandenEisen等人对医学图像配准进行的分类，归纳了七种分类标准。（一）按图像维数分类按图像维数分为2D/2D,2D/3D,以及3D/3D配准。2D/2D配准通常指两个断层面间的配准；2D/3D配准通常指空间图像和投影图像（或者是单独的「个层面）间的直接配准；3D/3D配准指2幅三维空间图像间的配准。（二）根据医学图像的模态分类根据医学图像的模态分为单模态医学图像配准和多模态医学图像配准。单模态图像配准是指待配准的两幅图像是用同一种成像设备获取的。•般应用在生长监控、减影成像等。多模态图像配准是指待配准的两幅图像来源于不同的成像设备，主要应用于神经外科的诊断、手术定位及放疗计划设计等。比如将MRLCT、DSA等解剖图像与SPECT、PET和EEG等功能信息相互结合，对癫痫进行手术定位。另外，由于MR适于肿瘤组织的轮廓描述而通过CT又可精确计算剂量，因此，在放疗中常需要将二者进行配准。多模态图像配准是医学图像配准的重点研究课题。（三）根据变换性质分类根据变换的性质可分为刚性变换、仿射变换、投影变换和曲线变换四种。刚性变换只包括平移和旋转；仿射变换将平行线变换为平行线；投影变换将直线映射为直线；曲线变换则将直线映射为曲线。（四）根据用户交互性的多少分类根据用户参与的程度，分为自动配准、半自动配准和交互配准。自动配准是用户只需提供相应的算法和图像数据；半自动配准是用户需初始化算法或指导算法（如拒绝或接受配准假设）；交互配准是用户在软件的帮助下进行配准。（五）根据配准所基于的图像特征分类根据配准所基于的图像特征分为基于外部特征和基于内部特征两大类。外部特征的图像配准是指在研究对象上设置一些标志点（这些标记点可以是立体定位框架、在颅骨上固定螺栓和在表皮加上可显像的标记），使这些标记点能在不同的影像模式中显示，然后再用自动、半自动或交互式的方法用标记将图像配准。基于内部特征的配准方法主要包括三个方面：基于标记的配准方法、基于分割的

配准方法、基于像素特性的配准。基于标记的配准方法分为解剖知识的标记（如：利用人体特殊的解剖结构，i般由人工直接描述）和基于儿何知识的标记（如：运用数学知识得到大量的点、线、面的曲率、角落特征等）；基于分割的配准指通过图像分割获得一些配准标志；基于像素特性的配准方法是把图像内部的灰度信息值作为配准的依据。又可分为两种。一是把图像灰度信息简约成具有一定的尺度和方向的集合（如力矩主轴法）二是在配准过程中始终使用整幅图像的灰度信息（如：互相关法、最大互信息法）。（六）根据配准过程中变换参数确定的方式分类根据配准过程中变换参数确定的方式可以分为两种：一是通过直接计算公式得到变换参数的配准，二是通过在参数空间中寻求某个函数的最优解得到变换参数的配准。前者完全限制在基于特征信息（例如小数目的特征点集、二维曲线、三维表面）的配准应用中。在后者中，所有的配准都变成一-个能量函数的极值求解问题。（七）根据主体分类根据主体可分为：（1）同一患者（Intrasubject）的配准。是指将来自同一个病人的待配准图像，用于任何种类的诊断中；（2）不同患者（Intersubject）的配准。指待配准图像来自不同病人，主要用在三维头部图像（MR、CT）的配准中，既可以基于分割也可以基于灰度。变换方式多为非线性的曲线变换，有时也采用刚性变换。（3）患者与图谱的（Atlas）图像配准。是指待配准图像一幅来自病人，一幅来自图谱；主要用于收集某些特定结构、大小和形状的统计信息。目前典型的数字化医学图谱是法国Talairach和Toumoux制作的Talairach-Tournoux图谱（TTAtlas）o图谱和实际图像配准后，能更直观和方便地应用图谱中的信息。三、医学图像配准的基本过程医学图像配准一般由以下三个步骤组成：（1）根据待配准图像（浮动图像）b与参考图像（基准图像）L,提取出图像的特征信息组成特征空间；（2）根据提取出的特征空间确定出一种空间变换（T）,使待配准图像L经过该变换后与参考图像L能够达到所定义的相似性测度，即I尸T（L）；

(3)在确定变换的过程中，还需采取一定的搜索策略也就是优化措施以使相似性测度更快更好地达到最优值。当然，配准过程并不绝对按上述步骤进行，比如一些自动配准方法，一般不包括特征提取步骤。此外，(2)和(3)的过程在实际计算过程中是彼此交叉进行的。图8-2给出了医学图像配准的流程图：图8-2医学图像配准流程图第三节医学图像配准的理论基础一、图像配准原理由图8-2,可以看出：对于在不同时间或/和不同条件下获取的两幅图像4(X)和B(X)的配准，就是要定义一个相似性测度并寻找一个空间变换关系，使得经过该空间变换后两幅图像间的相似性达到最大(或者差异性最小)。即使图像A上的每一个点在图像B上都有唯的点与之对应，并且这两点应对应同一解剖位置。用公式表示如下：S(T)=S(4(X),B(T(X))) zon式中，S是相似性测度，配准的过程可归结为寻求最佳空间变换的过程。

(8.2)T=argmaxS(T)(8.2)由于空间变换包含多个参数，是一个多参数最优化问题，所以一般由迭代过程实现：T'=T+AT (8.3)增量AT可以用计算梯度的方法计算：g包 (8.4)dT二、空间变换图像A和B的配准就是寻找一种映射关系T：Xa~Xb，使得Xa上的每一点在Xb上都有唯一的点与之对应。这种映射关系表现为一组连续的空间变换，如整幅图像应用相同的空间变换，则称之为全局变换(globaltransformation),否则，称之为局部变换(localtransformation)o根据图像变换形式的不同，有线性变换(lineartransfbrmation)和非线性变换(nonlineartransformation,也称作弯曲变换(curvedtransformation))两种。线性变换包括刚体变换(rigidbodytransformation),仿射变换(affinetransformation)和投影变换(projectivetransformation)。如图8-3所示。图8-3图像配准的基本变换(一)刚体变换所谓刚体，是指物体内部任意两点间的距离保持不变。刚体变换使得一幅图像中任意两点间的距离在变换前后保持不变。例如人体的头部由坚硬的颅骨支

撑，在处理时通常忽略头部皮肤的微小变形，将整个人脑看作是一个刚体。两幅图像之间的刚体变换可由一个刚体模型描述：V=sRU+T （8.5）其中，S是比例变换因子。7=（%仆山），是图像之间沿为，y,Z方向上的平移量。R是3X3的旋转矩阵，满足约束条件：R'R=RR'=1,det（/?）=1 （8.6）相对笛卡尔坐标系的三个坐标轴，/?有三种不同的形式：’10 o、&=0cos0xsin0x(8.7)-sin4cos0K)'cos40-sing'R、=0 1 0(8.8)sin0 0cos0\y yz'cos夕-sin0.0、R：=sin0_cos6_ 0(8.9)Looi,其中8，ev,d分别表示围绕x,r,z坐标轴的旋转角度。（二）仿射变换二维仿射变换的旋转矩阵叫3、/«23 (8.10)1,当式（8.6）的约束条件不满足时，方程式（二维仿射变换的旋转矩阵叫3、/«23 (8.10)1,'叫1叫2R'=加21m22

、00三维仿射变换的旋转矩阵R表示为：叫41团叫41团24加341JTOC\o"1-5"\h\zr,=W21加22 m23m3\叫2 m33、0 0 0仿射变换的具体表现可以是各个方向尺度变换系数••致的均匀尺度变换或变换系数不一致的非均匀尺度变换及剪切变换等。均匀尺度变换多用于使用透镜系统的照相图像，在这种情况下，物体的图像和该物体与成像的光学仪器间的距离有直接的关系，一般的仿射变换可用于校正CT台架倾斜引起的剪切或MR梯度线圈不完善产生的畸变。（三）投影变换与仿射变换相似，投影变换将直线映射为直线，但不保持平行性质。投影变换主要用于二维投影图像与三维体积图像的配准。二维投影变换按照下式将图像A（x，m）映射至图像8（%2，%）：(8.12)_/丙+。|2弘+卬3

。3111+432M+(8.12)(8.13)_〃2内+〃22yl+。23/2― (8.13)a3\X\+-32yl+。33其中的变换参数与是依赖于图像本身的常数。另外一种类似的变换方式称为透视变换。透视变换是投影变换的子集。某些医疗成像设备，如内窥镜、显微镜等，获取的图像都是通过将三维物体投影到二维平面，由此产生的儿何变换称为透视变换。（四）非线性变换非线性变换是把直线变换为曲线。它反映的是图像中组织或器官的严重变形或位移。典型的非线性变换是多项式函数，如二次、三次函数及薄板样条函数。有时也使用指数函数。非线性变换多用于使解剖图谱变形来拟合图像数据或对有全局性形变的胸、腹部脏器图像的配准。.二阶多项式变换基于二阶多项式变换的公式如下：x2=«00+«01-^1+a02yl+fl03^l+aO4xl+«05-9，+«06^+»07>,|2+a08K+«0922丫2=fll0+a\\x\+a12yl+0132!+a14xl+al5-xy+al6x^+anyi+al8yz+a1922（8.14）Z2=a20+a2lxl+a22yl+023Z1+。24》：+。25孙+。26》2+。27N：+a2sy^+a29^这里涉及30个变换参数.薄板样条变换薄板样条函数变换可以表示为仿射变换与径向基函数的线性组合:f(X)=AX+B+y1¥,(/(|P,-X\) (8.15)i=l其中，X是坐标向量，A与8定义一个仿射变换，。是径向基函数，在二维图像配准中TOC\o"1-5"\h\zU(r)=r2logr2 (8. 16)r=^Jx2+y2 (8. 17)对三维图像配准t/(r)=|r| (8. 18)r=ylx2+y2+z2 (8. 19)三、参数的优化搜索配准的儿何变换参数根据求解方式可分成两类，-是根据获得的数据用联立方程组直接计算得到的，二是根据参数空间的能量函数最优化搜索得到。前者完全限制在基于特征信息的配准应用中。在后者中，所有的配准都变成一个能量函数的极值求解问题。因此图像配准问题本质上是多参数优化问题，所以优化算法的选择至关重要。常用的优化算法有:Powell法、梯度下降法、遗传算法、模拟退火法、下山单纯形法、Levenberg-Marquadrt法等。下面简要介绍最为常用的搜索方法:Powell法、梯度下降法和遗传算法。(―)Powell法Powell法是一种传统的确定性优化方法，又称为方向加速法，由M.J.D.Powell于1964年首先提出。基本含义是:对于n维极值问题，首先沿着n个坐标方向求极小，经多n次之后得到n个共挽方向，然后沿n个共挽方向求极小，经过多次迭代后便可求得极小值。原理如下：对于某一问题，将其归结为求取某一目标函数J=j(y)的极小值。其中y为一个向量：y=[yi，y2,-yn]T-设置一个满秩的步长矩阵:

Dl(8.20)u\ndDl(8.20)其中，O,=[dn,di2,--din]o对于某-初始值丫=%=[yoi，y1）2,…，Powell算法的迭代过程如下:（1）首先在方向上搜索,即求5，使1/=1/*（）+5。1）为极小，并令=丫0+。1。1（2）依次求。2，使•/=（匕+。2。2）为极小，并令彩=匕+。2。2，如此下去（3）最后求Q"，使•/=•/（匕i+Q“O“）为极小，并令匕==+£]“£>“（4）令i=l,2,…,”-1,=人（匕-%）在新的。“方向上在搜索一次,即求a使J=J（匕+Q。,,）为极小，并令新的Yq为丫°=Y,,+OD„o至此，完成了第一轮n+1次的搜索。接下去进行下一轮的搜索，直至性能指标满意或满足某种停止条件为止。该法不需要对目标函数进行求导计算，具有直接法的优点，且具有二次收敛性，收敛速度快、精度高、可靠性好，被公认是目前解无约束最优化问题十分有效的直接法。（二）梯度下降法该算法在求最小化过程中直接利用梯度信息，沿着起始点梯度方向的反方向，求出最小值点，然后移动到最小值点，再重复上面的过程，直到前后点的函数值的差小于给定的误差值，则结束迭代过程。梯度下降法不是一个好算法，梯度方向并非直接指向优化的最终方向，优化过程中前后两步之间方向相互垂直,步长很小，需要多次迭代才能达到最优结果。实际上，当目标函数F的一阶导数连续时，梯度算法效果比较好。（三）遗传算法遗传算法（GeneticAlgorithm）是模拟达尔文的遗传选择和自然淘汰的生物进化过程的计算模型，是一种通过模拟自然进化过程搜索最优解的方法，它是由美国Michigan大学J.Holland教授于1975年首先提出来的。在求解优化问题时，遗传算法将优化问题当作一个生存环境，问题的一个解当作生存环境中的•个个体，以目标函数值或其变化形式来评价个体对环境的适应能力，模拟由一定数量个体所组成的群体的进化过程，优胜劣汰，最终获得最好的个体，即问题的最优解。它呈现出的是一种通用算法框架，该框架不依赖于问题的种类，因而具有较强的鲁棒性，特别是对于一些大型复杂非线性系统，表现出比其它传统优化方法更加独特和优越的性能。其隐含并行性和全局搜索特性，保证算法能够在大区域中作快速搜索，有较大把握寻找到全局最优解。目前国内外研究中，最受关注的是Powell法，因为Powell法与遗传算法都是无需求导数的直接优化法，因此可以适用于搜索中的任何空间限制。遗传算法中的杂交和变异操作可以避免使算法陷入局部最优，从而有很强的优化能力，但是速度较慢，而Powell法的优化速度较快，但容易陷入局部最优。遗传算法中实现了并行计算。若以增加时间为代价来找到更多的命中参数，则遗传算法较为理想，尤其是在有能量约束时。对于参数相对较少的配准来说，一般还是选择Powell算法，以减少配准所需的时间。在实际应用中，经常使用附加的多分辨率和多尺度方法加速收敛，降低需要求解的变换参数树木、避免局部最小值，并且多种优化算法混合使用，即开始时使用粗略的快速算法，然后使用精确的慢速算法。四、插值方法在图像配准中，空间坐标变换后得到的像素坐标位置可能不在整数像素上，因此需要用灰度插值的方法对像素值进行估计。常用的插值方法有:最近邻插值(NN,NearestNeighbor)法、双线性插值(BLBilinearInterpolation)法和部分体积分布(PV,PartialVolumnDistribution)法等。(一)最近邻插值(NN)该方法是一种简单的插值算法，也称为零阶插值。设需要插值的点为〃，在二维图像中，临近该点的落在坐标网格上的像素点分别为nt,小，小，最近邻法直接计算〃和邻近四个点之间的距离，并将与该点距离最小的点的灰度值赋给〃，如图8-4所示。计算公式如下：(8.21)/(«)=/(v),v=argmin(J(«,«,))

(8.21)这种方法简单快捷，但当邻近点之间的像素灰度差别很大时，会产生较大的误差。图8-4最近邻插值示意图(二)双线性插值(B1)双线性插值法又称为一阶插值算法。它是使用线性插值来求像素灰度的一种方法。具体计算方法为先沿着一个坐标轴方向使用线性插值方法求出两点的插值灰度，然后沿另一个坐标轴，利用这两个点对目标点进行线性插值来求灰度。计算方法如图8-5所示：计算公式为：/(")=工3,/(",) (8.22)上式中/(〃,)为它们的灰度值，3,为各相邻点的权重，与它们到〃的距离成反比，表达式如下:(8.23)31=(1-必(1-dy)

0)2=cbe(}-dy)

(1)3=dx^dy

(x)4(8.23)其中：dx,dy分别是〃与力之间沿x,y方向的距离。双线性插值方法由于考虑到直接邻近点对待插值点的灰度的影响，因此一般能得到令人满意的插值效果。但这种方法具有低通滤波性质，使高频分量受到损失。此外，由插值所得到的灰度值是经过数字计算出来的，一般不会是整数值，而且也有可能产生原始图像中所没有的灰度值，因此可能会改变图像中的灰度分布，特别是当图像中有很多需要进行插值的像素点时。图8-5线性插值和部分体积插值示意图(三)部分体积插值法(PV)部分体积分布法是F.Maes等人提出来的，是对双线性插值方法的一个改进。主要是为了克服双线性插值方法在图像中会产生新的灰度值而引起图像灰度分布发生变化的缺点，以便得到比较光滑的目标函数，有利于优化搜索。PV是根据线性插值的权重分配原则，将每对像素对联合直方图的贡献分散到联合直方图中与之相邻的各个像素对上，这样联合直方图上各个像素对的频度值以小数增加，因此不会出现新的灰度值而破坏目标函数值分布的光滑性。PV方法具体的计算公式为：献〃“),/(，))=W(“)J(Vj))+3, (8.24)3,.为权重，其取值同BI法。实际上，PV方法只是用灰度统计来代替插值，这反应了医学图像配准问题与其它医学图像处理问题中的不同之处。特别地，对基于灰度的配准方法来说，需要进行统计的是图像中的灰度信息而不是每点的灰度值，这意味着在处理的过程中不一定要得到每点的灰度值，因此可以采用PV方法来得到分布较好的目标函数。五、相似性测度配准过程在得到几何变换后，进一步的工作就是要找到一种合适、最优的描述量，用以表征相似或者差异，称这种描述量为相似性测度。下面介绍儿种常用的相似性测度。(一)灰度均方差设力⑴和方⑴分别表示参考图像和浮动图像中的数据，两幅图像像素值的均方差可以表示为:/—J(/r⑴-。(力•⑴))成 (8.25)II"llxeV其中，V表示参与计算的图像区域，|比||表示参与计算的像素总量，QCD表示对图像数据的变换。灰度均方差作为相似性测度适用于单模医学图像的配准问题，它的优点是易于理解且优化过程相对简单。(二)归一化互相关归一化互相关公式如下：R= "3 _ __ (8.26)X(/何(iJ)-%)2^(/A,(i,j)-7A,)其中：几/和心，分别是参考图和浮动图在(i,j)eT区域内的像素灰度平均值，R为相关系数。(三)互信息相似测度形式多样，在基于全图像信息的图像配准中，以互信息量作为相似测度的方法以其计算复杂度低、鲁棒性好等特性逐渐成为当前研究的热点。本章下一节将对基于互信息量的图像配准方法做详细的介绍。第四节医学图像配准的主要方法医学图像配准的方法有很多种，目前主要的配准方法大体上可以分为两类:基于特征的配准方法和基于灰度的配准方法。下面重点介绍这两种常用的图像配准方法。一、基于特征的配准方法基于特征的配准方法首先要对待配准图像进行预处理，也就是特征提取的过程，然后利用提取到的特征完成两幅图像特征之间的匹配。由于图像中有很多种可以利用的特征，因而产生了多种基于特征的方法。常用到的图像特征有:点、直线段、边缘、轮廓、闭合区域、特征结构以及统计特征如矩不变量、重心等等。(一)基于点特征的配准点特征是图像配准中最为常用的图像特征之一，分为外部特征点与内部特征点两种。外部特征点是成像时固定在患者身体上的标记物，不同成像时灌入不同的显影物质使得标记物在所有图像模态中均能清楚可视和精确检测。这种方法所求参数可用联立方程组直接计算得到，但标记物的固定对人体是侵入性的。这种方法的配准变换被限制为刚性变换，临床上常使用的基于立体框架的配准方法精度最高，它是用螺丝旋入头骨将其固定在患者的外颅表面，可作为其他配准算法评估的金标准，这种方法主要应用在神经外科手术的定位和导航（精度1mm之内）。目前已出现了多种对患者友好的非侵入性标记物，或是为个体定制的泡沫面具，或是用定位栓将特制的面具固定在患者头颅上，或是用特制的牙套，或是使用个体定制的鼻部支撑物和两耳的插件形成一种头部固定架，这些方法的配准误差均不超过2mm。内部特征点是一些有限的可明显识别的点集，标记点可以是解剖点（一般由用户识别出），也可以是儿何点（包括边缘点、角点、灰度的极值点、曲率的极值点、两个线性结构的交点或某一封闭区域的质心等）。这种方法主要求解刚体或仿射变换，如果标记点数目足够多，也能用来更复杂的非刚体变换。识别出来的标志点集与原始图像信息量相比是稀疏的，这样参数优化相对比较快。在基于点特征的图像配准方法中，Besl等首先将ICP策略引入到图像配准算法中，通过迭代过程使两点集间距离不断减小，最终实现3D点特征的配准Chui等提出一种更为通用的配准框架TPS-RPM（ThinPlateSplines-RobustPointMatching）,能够确定特征点的对应关系并弹性地配准点特征。在国内，张煜等提出了一种通过离散轮廓半自动地提取特征点的方法，然后采用平滑的薄板样条函数对特征点进行插值，该方法能够有效地减弱特征点对应位置误差对配准结果产生的影响。张二虎等将互信息相似性测度引入到点配准算法中，首先建立起两特征点集间匹配对应关系的联合概率分布匹配矩阵，通过最大化嫡和互信息最大化建立一个包含匹配矩阵和空间变换参数的能量函数，最后采用退火算法获得最优解。（二）基于直线特征的配准线段是图像中另一个易于提取的特征。Hough变换是提取图像中直线的有效方法。Hough变换可以将原始图像中给定形状的曲线或直线变换到变换域空间的•个点位置。它使得原始图像中给定形状的曲线或直线上所有的点都集中到变换域上的某一个点位置从而形成峰值。这样，原图像中的直线或曲线的检测问题就变成寻找变换空间中的峰点问题。正确地建立两幅图像中分别提取的直线段的对应关系依然是该方法的重点和难点。综合考虑直线段的斜率和端点的位置关系，可以构造一个这些信息指标的直方图，并通过寻找直方图的聚集束达到直线段的匹配。（三）基于轮廓与曲线特征的配准近年来，随着图像分割、边缘检测等技术的发展，基于边缘、轮廓的图像配准方法逐渐成为配准领域的研究热点。分割和边缘检测技术是这类方法的基础，目前已报道的有很多图像分割方法可以用来做图像配准需要的边缘轮廓和区域的检测，比如Canny边缘提取算子，拉普拉斯-高斯算子（LoG）,动态阈值技术,区域增长等。在特征提取的基础上，很多学者针对轮廓、边缘等进行了配准研究。Govindu等采用轮廓上点的切线斜率来表示物体轮廓，通过比较轮廓边缘的分布确定变换参数。Davatzikos等提出了…种二阶段大脑图像配准算法，在第•阶段使用活动轮廓算法建立一一影射，第二阶段采用弹性变换函数确定轮廓的最佳变换。李登高等提出了一种对部分重叠的图像进行快速配准的方法，该方法是基于轮廓特征的随机匹配算法。通过提取轮廓上的“关键点”作为特征点，随机选择若干特征点对得到候选变换，随后的投票阶段对其变换参数进行检验和求精。赵训坡等提出一种基于证据积累的图像曲线粗匹配方法，比较有效地解决了将图像中提取的一条曲线（较短）与一条参考曲线（较长）相匹配的问题。（四）基于面特征的配准基于面的配准方法中最典型的算法是由Pelizzari和Chen提出的“头帽法”（Head-HatMethod）o即从一幅图像中提取一个表面模型称为"头"（Head）,从另外一幅图像轮廓上提取的点集称为“帽子”（Hat）。用刚体变换或选择性的仿射变换将“帽子”的点集变换到“头”上，然后采用优化算法使得“帽子”的各点到“头”表面的均方根距离最小。头帽法最初用于头部的SPECT和CT（或MRI）配准，参考特征是头部的皮肤表面；然后用于头部的SPECT图像之间的配准，参考特征是头卢页骨表面和大脑表面。优化算法目前一般用Powell法。均方距离是六个待求刚体变换参数的函数，其最小时可得刚体变换参数。许多学者对

该算法作了重要改进，例如用多分辨金字塔技术克服局部极值问题；用距离变换拟合两幅图像的边缘点(EdgePoints),斜面匹配技术(ChamferingMethod)可有效地计算距离变换。另外比较常用的配准方法还有迭代最近点算法(ICP)。迭代最近点(ICP)配准算法:由Besl和Mckay提出的，它将•般的非线性最小化问题归结为基于点的迭代配准问题。迭代最近点算法(ICP)中必须先采样出图像结构上的特征点，然后用迭代的方法不断求出一幅图中相对于另一幅图中所有采样点的最近点，直到两个点集的均方差低于设定阈值，这时可得到匹配变换参数。除了采用分割的方法提取两幅图像中脑外表面轮廓特征外，还有用多尺度算子提取脑内部儿何特征，然后采用相关方法在多尺度空间结合外表面特征和内部特征进行自动配准的方法。也有采用平面变形轮廓和样条插值提取手术前CT图像的表面轮廓点集，通过最小化从二维轮廓到三维表面的投影线的能量而达到与手术中所获得的脊椎点集配准的目的。(五)实用算法举例下面给大家介绍一种基于特征点的刚体变换配准算法。设在待配准的两幅图像上选择的特征点集分别为：丫,=｛为"=(4,仇),《也为直角坐标值｝,i=0,l,2…N(N为点的个数)Xi={Xj|Xj=(a^bj),i=0,l,2…N。由于X,和匕上的点是一一对应的，所以它们点的个数是相同的。设刚性变换F为待求的最佳变换，F可以表示为一个旋转变换R和一个平移变换T的组合。R和T可以表示为：(cos6-sinfrR=.八八，7=t (8.27)^sin0cos6J可得下式：%=F(Xj)+j,j=0,1,2,-^ (8.28)其中E,.为误差矢量项。于是求取最佳变换F可表示为最小化下面的均方误差:minEminE=min—

ffNXu1=0(8.29)我们称上式为代价函数。为了求得F变换中的刚性变换，传统的方法是应用迭代法，这种方法的时间开销比较大，并且需要较多的配准点。为此我们选用基于奇异值分解(SVD)的最小二乘算法。此算法只需较少的配准点就能快速计算出旋转变换矩阵并同时算出平移矢量。将最小化均方误差该为对下式最小化：।N<8.30)Ni=Q算法过程如下：①计算X,和y,的坐标中心点N _NX=—>X,,F=—>V, (8.31)②计算每个特征点相对与中心点的位移x\=Xj-X,y\=y,-V (8.32)③计算矩阵H：NT(8.33)i=0④求H的奇异值分解：//=1/aV7 (8.34)其中 =/\=diag(入i，入2,入3)，且人।士入2 >0⑤求R:R=Vdiag(\,\,det(VU))UT (8.35)⑥求T:t=Y-RX (8.36)二、基于灰度的配准方法基于灰度的配准方法是目前研究得较多的一种方法，它直接利用图像的灰度信息进行配准，从而避免了因分割带来的误差，因而具有精度较高、稳健性强、不需要预处理而能实现自动配准的特点。基于灰度的配准有两类主要的方法，…类是通过图像灰度直接计算出代表性的比例和方向等要素；另一类是配准过程中使用全部的灰度信息。第一种方法以力矩和主轴法为代表，第二种方法一般称为体素相似性。(一)力矩和主轴法力矩和主轴法是指先用经典力学物体质量分布的原理计算出两幅图像的质心和主轴，再通过平移和旋转变换使两幅图像达到配准。利用此方法，图像可以模型化为椭圆形区域的点分布。这样的分布可以用这些点的位置的一阶和二阶矩描述。该方法对数据的缺失较敏感，要求整个物体必须完整地出现在两幅图像中。从整体上来说，配准精度较差，所以目前它更多地用来进行粗配准，使两幅图像初步对齐，以减少后续主要配准方法的搜索步骤。矩和主轴法的一个应用是将多个三维图像合成到标准脑的坐标系下，从而得到统计学平均意义上的脑模型，我们以这个例子重点讲述矩和主轴法的配准过程。首先使用阈值法和数学形态学等方法逐个层面地进行处理，半自动地得到脑的二值化模板，并通过基于形状的插值方法得到三维脑二值化模板B。然后按照下面的步骤将各个脑配准：①定义模板B：D,、_U 如果(x,y,z)在脑结构内[0如果(x,y,Z)不在脑结构内 (8.37)②计算模板中心，用(Xg,%,Zg)表示：y,z)xty,zXg=干 y,z)x,ytzZyB(x,y,z)ye (8.38)2^8(x,y,z)x,y,zZzB(x,y,z)"Zs(x,y,z)③计算旋转角度:B(x,y,z)的主导轴为下列矩阵的特征向量(8.39)(8.39)其中；&=>(丫-工)2+(z-Zg)2]B(x,y,z)x.y.zlyy=Z1(x-Xg)2+(Z-ZJ]fi(x,y,z)x,y,z心=Z[(x-Xg)2+(y-yg)2]B(x,y,z)x,y,zIxy=^(x-xg)(y-yK)B(x,y,z) (8.40)x,y,zIyz=2L(y- )(z-zg)B(x,y,z)x,y,z%=X(x-4)(z-Zg)8(x,y,z)x,y,zl的归一化特征向量矩阵E等于三个方向的旋转矩阵的乘积。设E为(8.41)TOC\o"1-5"\h\ze\\ e\2 e\3(8.41)E= ^21 e22 e23e31 e32 e33设x,y,z三个轴方向的旋转角度为aPy,对应矩阵为此,阵乘积cosy/?Y/?p/?Qcosy/?Y/?p/?Q=-siny0sinY0cosy00 1cosP00 1sin30-sinpFl0 0cosP00 0cosasina

-sinacosa(8.42)(8.43)(8.43)(8.44)求得：3=arcsin(^31)Y=arcsin(-e2(/cos0)a=arcsin(-e32/cos0)④计算缩放比例：比例因子：八松%为标准脑体积，彩为实际脑模板体积。最后，对每个实际脑执行上述步骤，并利用移位、旋转和缩放变换，将实际脑转换到标准脑的坐标空间，利用简单的加权平均就可以得到脑的空间概率分-Ar仲o(―)体素相似性法体素相似性法是目前研究较多的一类方法。由于它利用了图像中的所有灰度信息，因此这种方法一般都较为稳定，并能获得相当准确的结果。该方法还有一个优点是它是完全自动的，且不需要特殊的预处理，但这种方法由于需要大量的复杂计算，因此最近几年才转入实际应用。常见的基于体素相似性的配准方法有：①互相关法②基于傅立叶域的互相关法和相位相关法③灰度比的方差垠小化法④直方图的互信息最大化法等。.互相关法(Cross-CorrelationMethod) 互相关函数配准是应用广泛的一种配准方法。灰度互相关配准就是在参考图像上选择一临时窗口卬亿m),在目标图像SG,/)上寻找与其对应的最相似的窗口57(/,机)(移动窗口)，用互相关相似性测度为：MMR(i,j)= («*，/) (8.45)/=!m=\或定义归一化互相关函数：MMr (」，/*) l<i,j<L-M+\ (8.46)MM MM..xxs2‘4")1_1=1m=1 _in=\其中W(i,j)是大小为MxM的窗口图像，S7(/,⑼是参考图像中的临时窗口。窗口图像W(i,/)在搜索图像S中以扫描方式搜索，对该过程中每一位置相关函数值以亿力中最大值点的位置(门广)，并认为此处就是配准点位置。在基于傅立叶分析的测度中主要有两种，••是将空域中的互相关在频域中进行计算，另一个是相位相关法。由傅立叶的相关定理可知，两个函数在定义域中的卷积等于它们在频域中的乘积，而相关是卷积的一种特定形式。因此，可以利用快速傅立叶变换来求解相似度。相位相关法可以用于计算两幅图像的平移。基于傅立叶变换的性质，两幅图像的平移可以认为是傅立叶变换的角度差别，计算两幅图像的互功率谱就可以得到两幅图像的角度差别。相关法主要限于单模图像配准，特别是对一系列图像进行比较，从中发现由疾病引起的微小改变。.最大互信息法最大互信息法以互信息作为相似性测度。互信息(MutualInformation,MI)是信息论中的一个基本概念，用于描述两个系统间的统计相关性，或者是在一个系统中包含的另一个系统的信息的多少，一般用场来表示，表达的是•个系统的复杂性或不确定性。1995年分别被Viola和Collignon等首次用于医学图像配准中。对于概率分布函数为p(a)的随机变量集A,其病，(4)定义如下：H(A)= p(a)logp(a)aeA (8.47)对两个离散的随机变量A和B,它们的边缘概率分布函数分别是p(a)和p(b),联合概率分布函数是P33,则随机变量A和B的联合燧定义如下：H(A,B)=p(a,b)\ogp(a,h)aeAbeB (8.48)如果”(A|8)表示已知系统B时A的条件嫡，那么”(A)与H(A|8)的差值，就代表了在系统B中包含的A的信息，即互信息。因此两个系统间的互信息定义为：=(849)在医学图像配准中，虽然两幅图像来源于不同的成像设备，但是它们基于共同的人体解剖信息，所以当两幅图像的空间位置达到完全•致时，其中--幅图像表达的关于另一幅图像的信息，也就是对应像素灰度的互信息应为最大。通常用联合概率分布和完全独立时的概率分布间的广义距离来估计互信息：/(4,B)=Zp(aMlog要2、 (8.50)J P(a)P(b)对于离散的数字图像，联合概率分布外8(。向可以用归一化的联合直方图表

75：pAB（i,j）=J，（iJ） （8.51）边缘概率分布勿⑷表示为：PA（i）=ZPA8（i，/） （8.52）边缘概率分布/⑶表示为：PB（j）=£pAB（i，j） （8.53）则有：(8.(8.54)这就是用互信息表示的相似性测度。接下来的任务是寻找一个变换使得一幅图像经过此变换后和另-一幅图像的互信息最大。一般采用刚体变换，即在三维空间中寻找三个方向上的平移值和旋转角度。对于大规模断层扫描医学图像来说，三维体积数据集包含的数据量极大，无法满足临床上实时处理的要求，因此必须采取优化措施。常用无需计算梯度的Powell多参数优化算法。最大互信息法是目前应用较多的一种方法，其配准精度一般高于基于特征的方法。由于该方法不需要对图像做分割、特征提取等预处理，儿乎可以用于任何不同模式图像的配准，并具有较强的鲁棒性，特别是当其中一幅图像的数据部分缺损时，也能得到很好的配准效果。因此，从它一开始出现，就得到了学者的普遍重视和广泛应用。但是，研究者也发现该方法并不是尽善尽美的。从1998年开始就有相关文献指出基于最大互信息量的图像配准方法并不像Maes在文章中描述的那样，没有误配，精确度可以达到亚像素水平。Studholme通过研究发现，互信息本身的大小与待配准两图像间的重叠度具有一定的关联性。为了消除这种关联关系，他提出标准化互信息的方法，实验证明，它比标准的互信息方法更具有鲁棒性。其公式表示为(8.55)i(a,b)(8.55)其中”（A）,"（B）分别是图像A,B的边缘嫡,H（A,8）是它们的联合嫡，则配准过程就是寻找最优变换”的过程：7。=argmaxI(A,TB)=argmax叱);黑)T TH(A,TB) (8.56)它将B变换为t08,并且最可能多的包含参考图像A的信息。F.Maes等人也提出了类似的标准化互信息配准准则。第五节医学图像配准的评估医学图像配准，特别是多模医学图像配准结果的评估一直是件很困难的事情。由于待配准的多幅图像基本上都是在不同时间或/和条件下获取的，所以没有绝对的配准问题，即不存在什么金标准(goldstandard)。只有相对的最优(某种准则下的)配准。在此意义上，最优配准与配准的目的有关。常用的评估方法有以下几种：一、体模(Phantom)体模又有硬件体模和软件体模之分。后者是计算机图像合成的结果。体模法用已知的图像信息验证新配准算法的精度。由于体模都比较简单，与实际临床图像差异较大，因此只能对配准方法作初步的评估。二、准标(FiducialMarks)立体定向框架系统(StereotacticFrameSystems)包括立体定向参考框架、立体定向图像获取、探针或手术器械导向儿部分。优点是定位准确，不易产生图像畸变。使用立体定向框架系统的体积图像数据可以用来评估其它配准方法的精度。使用人工记号作准标的方法很多。一种准标是使用9根棍棒组成的3个方向的N字型结构。在作CT测试时，棒内充以硫酸铜溶液；作PET测试则填充氟18o这样，在两组图像中都可见此N字型准标，从而可对图像准确空间定位。例如用在人脑表面嵌螺丝作标记(每人8个)的方法对多个病人做CT、MRCTUT2及PD)和PET实测，得到多组数据。这些数据专门用于多模医学图像配准算法评估使用。三、图谱(Atlas)Thompson用随机向量场变换构造一个可变形的概率脑图谱。包括从多个受试者到单一-解剖模板的功能、血管、组织诸方面映射，三维图谱到新受试者的扫描图像的映射。VisibleHumanCD的CT骨窗图像、MR图像及彩绘的冷冻切片照片由于具有清晰的解剖结构和高度的分辨(1毫米/每层片)，近来也被用来做新配准方法精度的评估。四、目测检验(VisualInspection)对多模医学图像配准的结果请领域专家用目测方法检验，听起来有些主观，但在一定程度上的确是一种相当可信的方法。第六节医学图像融合技术一、医学图像融合概念医学图像融合(MedicalImageFusion)是指将两幅(或两幅以上)来自不同成像设备或不同时刻获取的已配准的图像，采用某种算法，把各个图像的优点或互补性有机地结合起来，获得信息量更丰富的新图像的技术。在图像融合处理中，图像配准是图像融合的第一步，也是实现图像融合的先决条件，只有实现了待融合图像的配准，才能实现相应组织之间的融合，如果对应组织的位置有较大的偏差，那么融合的图像是不准确的。只有两幅图像中同一空间位置的像素都对应相同的解剖结构，融合起来的图像才有意义。二、医学图像融合分类因为研究者的研究对象和研究目的不一样，所以发展和形成的图像融合分类也就多种多样。.按照被融合图像的成像方式不同，可以把融合分为单模融合(mono-modality)和多模融合(multi-modality)。所谓单模融合是指待融合的图像由同一设备获取的。简单地说，就是CT-CT或者MRI-MRI这种类似形式的融合处理。多模融合是指待融合的两幅或多幅图像来源于不同的成像设备，研究较多的CT与MRI图像融合和CT与核医学图像的融合都属于此类。.按照融合对象的不同，可分为单样本时间融合、单样本空间融合和模板融合。单样本时间融合是指跟踪某个病人，将其一段时间内对同一脏器所做的同种检查图像进行融合，以助于跟踪病理发展和研究该检查对该疾病诊断的特异性。单样本空间融合是指将某个病人在同一时期内（临床上视1-2周内的时间为同时）对同一脏器所做的儿种检查的图像进行融合，以便综合利用这儿种检查提供的信息（如MRI/CT可以提供脏器的结构信息，SPECT可以提供脏器的功能信息），对病情做出更准确的诊断。模板融合：从许多健康人的研究中建立一系列模板，将病人的图像与模板图像融合，有助于研究某种疾病和确立诊断标准。.按照图像处理方法的不同，又可分为数值融合法和智能融合法。数值融合法将不同来源的图像做空间归一化处理后直接融合。智能融合法将不同来源的图像做归--化处理后，根据需要选择不同图像中的所需信息再进行融合。.按图像类型不同，可以分为断层图像间相互融合、断层图像与投影图像融合以及结构图像与功能图像融合。断层图像间相互融合主要指CT与MRI图像融合；断层图像与投影图像融合主要指CT、MRI图像与DSA图像通过三维重建后进行融合；而结构图像融合与功能图像融合主要指CT、MRI图像与PET、SPECT图像进行融合。另外，还可以将图像融合分为前瞻性融合和回溯性融合。两者的区别在于：前瞻性融合在图像采集时使用特别措施（如加外部标志等），而回溯性融合在图像采集时则不采取特别措施。综上所述，依据不同的分类原则，图像融合有多种分类方式，应该指出，以上分类不是绝对的，孤立的，在实际应用中，一个融合系统的设计过程往往是综合各种分类概念来实现的。第七节常用的医学图像融合方法一、基于空域的图像融合基于空域的图像融合是指直接在空间域中对图像的像素点进行操作，该类方法简单直观，易于理解，但常常融合效果有限，只适用于有限的场合。（一）图像像素灰度值极大（小）融合法设gC,力和g2（i,及为待融合图像，/亿力为融合后的图像，其中八j为图像中某一像素的坐标，图像大小为M*N,贝hG[O,M-l],je[0,NT],g】(i,j),g2(i，j)e[0,255]=Max{gl（i,j）,g2（i,j）} （857）F（"= （i,力,g2（，,/）} （&58）此方法只需对两幅待配准图像取对应像素点灰度值较大（小）者即可。这种方法计算简单，效果有限，只适用于对融合效果要求不高的场合。（二）图像像素灰度值加权融合法加权法是将两幅输入图像g《,j）和g2。，力各自乘上一个权系数，融合而成新的图像EG,j）。=ag[j）+（\-a）g2（i,j） （859）其中：a为权重因子，且OWaWl,可以根据需要调节a的大小。该算法实现简单，其困难在于如何选择权重系数，才能达到最佳的视觉效果。（三）TOET图像融合方法①首先求输入图像由«,_/）和g2（i，，）的共同成分。gieg?=Mi〃{g”g2} （8.60）②从图像©上扣除共同成分得到图像g2的特征成分g；：g「=gi-giCg2 （8.61）同理得到g2的特征成分g2*：82=82-8x^82 （8.62）③从图像团中扣除图像g2的特征成分g2，，得到：gl-2'=（gl-g2）+giCg2 （8.63）同理，从图像g2中扣除图像的特征成分g；，得到：g2-gj=（g|-g2）+giCg2 （8.64）这项操作是为了改善图像的融合效果。（4）确定图像g2（i，j）和©力力的不同成分,g2-Si=82-81 （8.65）当|g2*|<|g:|时'定义g2*-g「=0。此操作的目的是将两幅图像的不同部分作为背景，突出图像g2（i"）的特征，以便准确判断©Gj）的位置；反之也行。该成分在融合图像中的比重由权重系

数决定，这里到底突出哪个图像的特征以及判断哪个图像的位置要根据实际情况确定。⑤将步骤3和步骤4中得到的结果按不同权重计算融合图像的灰度值:F，(hj)=a(gi-g2)+b(g2-g*)+c(g2'-gl，) 记46)其中：a,b,c为权重系数，且a+b+c=l,具体可根据需要选取。图8-6是分别采用以上方法对已配准的CT、MR图像进行的融合。CT图像MR图像像素灰度极小值法像素灰度极大值法TOET法CT图像MR图像像素灰度极小值法像素灰度极大值法TOET法图8-6几种融合方法比较二、基于变换域的图像融合变换域法，顾名思义，就是将变换后的两个或多个图像进行融合，再通过反变换得到融合后图像的方法。下面介绍儿种常用的基于变换域的图像融合方法。(一)、基于多分辨率的金字塔融合法这是最早的一种基于变换域的方法。在这种方法中，原图像不断地被滤波,形成一个塔状结构。在塔的每一层都用一种算法对这•层的数据进行融合，从而得到一个合成的塔式结构，然后对合成的塔式结构进行重构，最后得到合成的图像。合成图像包含了原图像的所有重要信息。但这类方法产生的数据有冗余，且不同级的数据之间相关。（二）、基于傅里叶变换的图像融合法傅立叶变换是数字图像处理技术的基础，其通过在时空域和频率域来回切换图像，对图像的信息特征进行提取和分析，简化了计算工作量，被喻为描述图像信息的第二种语言。基于傅里叶变换的图像融合法包括以下三个步骤：1、对每一源图像分别进行图像的二维傅里叶变换；2、对变换系数通过加权得到融合图像的傅里叶变换；3、对融合后的系数进行傅里叶反变换，得到融合图像。（三）基于小波变换的图像融合小波变换本质是一种高通滤波，当采用不同的小波基，就会产生不同的滤波效果。小波变换可将原始图像分解成系列具有不同空间分辨率和频域特性的子图像，可以针对不同频带子图像的小波系数进行组合，形成融合图像的小波系数。1、图像的二维小波分解及融合Mallat于1989年提出了图像的二维小波分解的Mallat快速算法，公式如下：CjrHCjH*Dhj+i=GCjH*'D'j+i=HCjG”Ddi+\=GCjG'j= （8.67）式中：6,分别表示水平、垂直和对角分量；H（低通）和G（高通）为两个•维滤波算子，/和G“分别是”和G的共胡转置矩阵；1/为分解层数。相应的小波重构算法为：Cj-fCjH+G*%H+H*D'jG+G*D：G （8.68）图像经二维小波变换分解后，可得到四个不同的频带LL、LH、HL、HH。其中低频带LL保留了原图的轮廓信息。HL、LH、HH分别保留了原图水平、垂直和对角方向的高频信息，代表图像的细节部分。然后在对子图像分解得到LL2,HL2,LH2及HH2,依次进行多层分解。N层小波分解后可得到（3N+1）个频带。基于小波分解的图像融合的本质是采用不同的滤波器，将源图像分解到一系列的频率通道中，然后针对系数特性采用不同的融合规则和融合策略。基于小波变换的图像融合具体步骤为：①分解：对每一源图像分别进行小波变换，得到每幅图像在不同分辨率下不同频带上的小波系数；②融合：针对小波分解系数的特性，对各个不同分辨率上的小波分解得到的频率分量采用不同的融合方案和融合算子分别进行融合处理；③逆变换：对融合后系数进行小波逆变换，得到融合图像。整个过程如图8-7所示：从图8-7可以看出设计合理的融合规则是获得高品质融合的关键。小波变换应用于图像融合的优势在于它可以将图像分解到不同的频率域，在不同的频率域运用不同的融合规则，得到合成图像的多分辨率分析，从而在合成图像中保留原图像在不同频率域的显著特征。小波分解图图8-7小波分解融合图2、基于小波变换的融合规则（1）低频系数融合规则通过小波分解得到的低频系数都是正的变换值，反映的是源图像在该分辨率上的概貌。低频小波系数的融合规则可有多种方法：既可以取源图像对应系数的均值，也可以取较大值，这要根据具体的图像和目的来定。（2）高频系数融合规则通过小波分解得到的三个高频子带都包含了一些在零附近的变换值，在这些子带中，较大的变换值对应着亮度急剧变化的点，也就是图像中的显著特征点，如边缘、亮线及区域轮廓。这些细节信息，也反映了局部的视觉敏感对比度，应该进行特殊的选择。高频子带常用的融合规则有三大类，即基于像素点的融合规则、基于窗口的融合规则和基于区域的融合规则。如图8-8所示。融合处理融合规则:图像A分解层1)2)3)4)5)基于窗口的融合规则图像F融合处理融合规则:图像A分解层1)2)3)4)5)基于窗口的融合规则图像F分解层基于区域的融合规则基于像素的融合规则系数加权绝对值选大区域能量最大系数模值极大图像B分解层图8-8小波融合规则第•类方法是逐个考虑源图像相应位置的小波系数，要求源图是经过严格对准处理的。因为基于像素的选择方法具有其片面性，其融合效果有待改善。第二类是基于窗口的融合规则，是对第一类方法的改进。由于相邻像素往往有相关性,该方法以像素点为中心，取一个MXN的窗口，综合考虑区域特征来确定融合图像相应位置的小波系数。该类方法的融合效果好，但是也相应的增加了运算量和运算时间。由于窗口是一个矩形，是规则的；而实际上，图像中相似的像素点往往具有不规则性，因此，近年来又提出了基于区域的融合规则。该类方法常常利用模糊聚类来寻找具有相似性的像素点集。下面介绍几种小波分解系数的融合规则。①小波系数加权法，如下式所示:(8.69)Cj(F,p)=aCj(A,p)+(l-a)Cj(B,p),(8.69)其中：Cj(A,p),Cj(B,p),g(尸,p)分别表示源图像A,B和融合图像F在J层小波分解时，在P点的系数，下同。②小波分解系数绝对值极大法:Cj(F,p)=Cj(A,PCj(F,p)=Cj(A,P)\CJ(A,p)\>\CJ(B,p)\"B,p)|Cj(4,p)」g(B,p)|(8.70)③小波分解系数绝对值极小法:Cj(F,p)=Cj(A,p)|CCj(F,p)=Cj(A,p)|C7(A,p)|<|Cj(B,p)\\CJ(A,p)\<\CJ(B,p)\(8.71)④区域能量最大法:在1/层小波分解的情况下，局部区域Q的能量定义为:E(4,p)=E(4,p)=Zgeo3(q)C；(4,g)(8.72)其中：0(q)表示权值，g点离p点越近极辘尢且；°是0的个邻或iwawE(B'p)Cj(F,p)=Cj(A,p)E(A,p)>E(B,p)Cj(B,p)Cj(F,p)=Cj(A,p)E(A,p)>E(B,p)Cj(B,p)E(A,p)<E(B,p)(8.73)以CT、MR两幅图像为例，进行小波融合。①选择小波db2,对所需融合图像进行小波两层分解。②对小波分解后的小波低频系数采用均值法，高频系数分别采用系数加权、绝对值极大极小法、区域能量法、直接进行融合。③以融合后的小波系数进行图像重构，得到融合后的图像如图8-9„均值、绝对值极小均值、区域能量最大均值、绝对值极小均值、区域能量最大图8-9原始图像及小波融合结果我们还对原始图像进行了不同层的小波分解,然后融合重构，实验结果显示,小波分解层数等于3、4的时候，效果最好。我们又选取了不同的小波函数sym5,sym6,db4,coif3,coifs,haar等小波函数对两原始图像进行融合，实验结果表明重构图像变换不大，基本相同。第八节医学图像融合效果的评价目前，图像融合效果的评价主要有主观评价和客观评价两种。主观评价以人作为观察者，对图像的优劣做出主观定性评价。人对图像的识别或理解不仅和图像的内容有关，而且还与观察者的心理状态有关。由于人的视觉系统很复杂，受环境条件、视觉性能、人的情绪爱好以及知识状况影响很大，因此主观评价具有主观性和不全面性，所以有必要把主观评价与客观的定量评价标准相结合，这样既便于人的观察，也便于将融合结果交于计算机进行处理。下面是一些常用的客观评价指标。一、.月(Entropy)图像的熠值是衡量图像信息丰富程度的•个重要指标，端值的大小表示图像所包含的平均信息量的多少。根据香农信息论的原理，一幅图像的信息嫡为：L-1H=PjInPj (8.74)i=0其中Pi为图像的直方图，即灰度值等于，•的像素数与图像总像素数之比。如果融合图像的嫡增大，表示融合图像的信息量增加，融合图像所包含的信息就越丰富，融合质量越好。二、交叉炳CE(CrossEntropy)交叉燧也称相对燧，直接反映了两幅图像灰度分布信息的差异。设源图像和融合图像的直方图分别为P，和“，则交叉嫡定义为：L-\CE=Yp,log2-^- (8.75)M/交叉炳越小，说明融合图像从源图像提取的信息量越多，融合效果越好。在实际应用中，可以选择平均值来描述融合图像与源图像的综合差异： C4-CCfab=、fb （8.76）三、交互信息量A//（MutualInformation）交互信息量为两个变量之间相关性的量度，或一个变量包含另一个变量的信息量的量度。假设两幅源图像A和8,将它们融合得到融合图像F,E与A、B的交互信息量分别表示为M/fa和MIfb：（&77）k-Qi=0 Pf*）Pa⑴M2骂/"加弗特 （8-78）式中Pa、Pb和Pf分别是a、B、尸的灰度直方图；PFA（k,i）和。万伏，/）分别代表两组图像的归一化联合灰度直方图。综合考虑这两个值量，用M/必和MIfb之和来表示图像融合后包含源图像A、B的交互信息量的总和。M，=MIfa+M%b （8.79）交互信息量的值越大，表示融合图像从源图像中获取的信息越丰富，融合效果越好。四、图像均值％图像均值是图像像素的灰度平均值，对人眼反映为平均亮度。图像均值的定义为：_ 1MN〃= VVG（x,y） （8.80）MxN金勺其中G（x,y）表示图像中第（x,y）个像素的灰度，图像尺寸为MxN。如果均值适中，则目视效果良好。五、灰度标准差3g图像的灰度标准差定义为为：忖（g>xp（g） （8.81）其中L为图像的总灰度级，g表示图像第（x,y）个像素的灰度，〃表示图

像均值，p(g)表示灰度值为g的像素出现的概率。标准差反映了图像灰度相对于灰度平均值的离散情况，若标准差大，则图像灰度级分布分散，图像的反差大，可以看出更多的信息。六、均方误差MSE(MeanSquareError)均方误差表示融合图像与标准参考图像之间的差异，定义为：MMZZW(M)-MSE= (8.82)M*N其中/G,j)为融合图像，R(i,j)为标准参考图像。均方误差越小说明融合图像与标准参考图像越接近。七、信噪比SNR(SignalNoiseRatio)与峰值信噪比PSNR(PeakSignalNoiseRatio)将融合图像与标准参考图像的差异看作噪声，而标准参考图像看作信息。融合图像信噪比SNR定义为：MN2£固，疥SNR=W\Sm；=，＞=1 (8.83)1=1j=l融合图像峰值信噪比PSNR为：(8.84)255(8.84)PSNR=\0\q MSE信噪比、峰值信噪比越高，说明融合效果就越好。八、平均梯度MG(MeanGradient)图像的平均梯度定义为：1MN (8.85)MG=FT汨£Z^xF(x,yy+^yF(x,y)(8.85)1V1XIVx=|y=l式中：△.vF(x,y)分别为尸(x,y)沿x方向和y方向的差分，定义如下:(8.86)(8.87)尸(x,y+1)一尸(x,y)+产(x+l,y+l)-/(x+l,y)

(8.86)(8.87)^yF-(x,y)=F(x+1,y)—F(x,y)+F(x+1,y+1)—P(x,y+^yF-(x,y)=平均梯度用来表示图像的清晰度，反映图像融合质量的改进及图像中的微小细节反差和纹理变换特征，平均梯度越大，则图像的清晰度越高，微小细节及纹理反映越好。在以上8种指标中，其中MSE、SNR、PSNR均是通过比较融合图像与标准参考图像之间的关系来评价图像融合的实际效果。在图像融合的一些应用中很难获得标准参考图像，所以这几种方法的使用受到一定限制。(巩萍，聂生东)作业8：1、简述图像配准的原理及实现过程。2、试举例说明医学图像配准在临床医学中的应用。3、常用的图像配准方法有哪些，比较各种方法的优缺点。4、试用MATLAB语言编程实现基于最大互信息的图像配准。5、常用的图像融合方法有哪些，并比较各种方法的优缺点。6、试述基于小波变换的图像融合步骤及常用的小波融合规则。7、图像配准与图像融合的评价标准有哪些。参考文献：1、罗述谦，周果宏，医学图像处理与分析，北京:科学出版社，2003。2、PelizzariCA.ImageProcessinginstereotacticplanning:volumevisualizationandimageregistration[J].MedDosim,1996,23(3):137-144.3、LewisPJ,SiegelA,SiegelAM,etal.DoesperformingimageregistrationandsubtractioninictalbrainSPECThelplocalizeneocorticalseizures[J].NuclMed,2000,41(10):216-229.J・B・Maintz,M.A.Viergever,ASurveyofMedicalImageRegistration,MedicalImageAnalysis,1998,1,1-36.5、B.Zitova,,J.Flusser,ImageRegistrationMethods:ASurvey,ImageAndVisionComputing,2003,21,977-1000.6、D.Hill,P.Batchelor,M.Holden,etahMedicalImageRegistration,PhysicsinMedicineandBiology,46,2001,1-45.7、J.P.W.Pluim,J.B.A.Maintz,M.A.Viergever,MutualInformationbasedRegistrationofMedicalImages:ASurvey,MedicalImaging,2003,986-1004.8、BrownLG,Asurveyofimageregistrationtechniques[J],ACMComputingSurveys,1992,24(4):325-376.9、RitterN,OwensR,CopperJ,etahRegistrationofstereoandtemporalimagesoftheretina[J],IEEETransMedimaging,1999,18(5):404-418.10>ZukTD,AtkinsMS,Acomparisonofmanualandautomaticmethodsforregisteringscansofthehe