光学分子影像手术导航系统中图像融合算法研究

徐同花, 赵建民, 朱信忠, 迟崇巍, 叶津佐

(1.浙江师范大学 数理与信息工程学院,浙江 金华 321004;2.中国科学院 自动化研究所,北京 100190)

0 引 言

据文献[1]报道,我国每年新发肿瘤病例约为312万例，且未来10年,我国的癌症发病率与死亡率仍将继续攀升.手术导航系统对于提高术中肿瘤定位精度、减少手术损伤、提高手术成功率等具有重要意义[2].光学分子影像手术导航系统具有快速安全、定位精确、分辨率高、灵敏度高、成像价格低等优点[3].图像融合是导航过程中肿瘤定位的关键技术，可以在3个不同层次上进行:像素级融合、特征级融合和决策级融合.像素级图像融合的准确度最高,能够提供其他层次融合处理所不具有的更丰富、更精确、更可靠的细节信息[4].像素级图像融合需要在精确配准的前提下进行,图像配准方法主要分为两大类:基于区域的方法和基于特征的方法[5-6].目前图像融合和图像配准方面的算法很多,但是应用于手术导航系统中时都存在运行时间长、无法满足实时性要求的问题.本文对光学分子影像手术导航系统中的图像融合算法进行了研究,采用像素级图像融合和基于区域的配准方法,仅利用图像的荧光信息和患者体表信息就可以在 50 s内完成配准,并且配准的精度和融合的效果都得到了很大的提高.

为了便于显示和理解，2幅图像的相对大小已做了处理图1 手术导航过程中采集的1组原始图像数据

1 方 法

1.1 实验数据

本文选用的图像数据是中国科学院分子影像重点实验室研制的光学分子影像手术导航系统在实施预临床及临床实验过程中采集的,目前该系统已在汕头大学医学院附属肿瘤医院成功实施了22例人体乳腺癌手术的实时导航.

每组图像数据包括大小为2 448×2 056的彩色图像(见图1(a))及大小为1 024×1 024的荧光图像(见图1(b)).图中荧光部分是病灶区域.每组图像数据都是通过2台CCD相机在术中实时采集的,其中彩色图像是在系统提供的可见光源下由彩色CCD相机拍摄的,荧光图像是采用激发荧光成像技术,利用近红外LED灯进行激发,由近红外CCD相机采集.彩色图像能够清晰显示病灶的解剖结构信息,包括病灶区域和手术器械以及它们的相对位置等.荧光图像用于显示病灶的特征信息,包括大小、边界等.本研究的目的是将荧光图像中荧光显示的病灶区域融合到彩色图像上(以伪绿的形式),让医生可以在彩色图像上看到病灶的准确位置及大小.

图2 算法流程图

手术过程中,导航系统的软件控制图像采集设备采集图像,并对采集到的每一组图像数据进行预处理、配准、融合,如图2所示.最后输出融合后的图像,引导手术的精确实施.

本算法运行的硬件平台为 Lenovo Win7 PC, Intel Core i3 2.53GHz双核CPU,2 GB RAM,软件平台为Matlab.R2010b.

1.2 图像预处理

为了减小手术过程中系统提供的光源强弱对采集到的图像光强的影响,以尽量消除图像亮度过亮或过暗对配准精度和融合效果的影响,首先对彩色图像和荧光图像的光强进行调节.为了提高算法的运行速度,减少图像融合所花费的时间,将彩色图像的像素(2 448×2 056)在水平和垂直方向以4为间隔采样,得到大小为612×514的彩色图像.将重采样后的彩色图像做灰度处理,得到其灰度图像.

1.3 图像配准

本文采用了基于区域的配准方法,以经过重采样的彩色图像的灰度图像gray_image为模板,选用穷尽搜索策略[5],在荧光图像的子图像中搜索匹配图像,分别以平均绝对差、均方差和归一化积相关[7]作为配准算法中的相似性度量.

1.3.1 匹配图像搜索原理[8]

设模板图像为X,大小为m×n;搜索图像为Y,大小为u×v,如图3所示.在图像Y中以某点S(i,j)为基点,截取一个大小与图像X一样的子图像Yk,利用相似性测度比较图像X和Yk的相似程度.这样的基点有(u-m+1)×(v-n+1)个,相应的子图像就有(u-m+1)×(v-n+1)个.在(u-m+1)×(v-n+1)个子图像中找1个与模板图像X最相似的图像Ykfinal,这样得到的基点S(i,j)就是最佳配准点O(ifinal,jfinal).

图3 搜索原理图 图4 确定基点搜索区域R2的示意图

1.3.2 基点搜索区域

确定基点搜索区域R2:

1)分别求出灰度图像的梯度矩阵Ag-g和荧光图像的梯度矩阵Ag-f,并以矩阵Ag-f和Ag-g的左上角点作为各自坐标空间的(0,0)点;

2)以矩阵Ag-f的一个最大值点I(i,j)作为兴趣点(ig-f,jg-f)=max(max(Ag-f,1)),如图4(a)所示.由于兴趣点偏离荧光图像中心太远会导致该配准算法的失败,本文对其进行了域法修正,当兴趣点I的坐标超出一定范围时,对其重新赋值,使兴趣点向图像中心靠拢;

3)将兴趣点I按矩阵Ag-f和矩阵Ag-g的维数比例对应到矩阵Ag-g中,如图4(b),得到兴趣点I在Ag-g中的对应点I′的坐标(i′,j′),即 i′=i*514/1 024, j′=j*612/1 024.由于荧光图像与灰度图像在视野范围上并不相同,因此,此对应点只是假设;

4)根据点I′的坐标得到Ag-g的左上角点P′ 在矩阵Ag-f中的对应点P的坐标为(i-i′,j-j′).将P点邻域内100×100个点组成的区域作为配准时基点的搜索区域R2,如图4(a)中虚线框所示.

1.3.3 图像相似性测度

图像的相似性度量方法有很多种,本文选用的平均绝对差度量、均方差度量和归一化积相关度量是其中应用比较广泛且比较简单的3种,前2种属于距离度量,后1种属于相关度量.

1)平均绝对差度量[7].计算公式如下:

(1)

基于平均绝对差度量的匹配算法的基本原理是以模板图像X与子图像Yk的灰度值作差,将差值的绝对值取平均得到d值.d值越小,说明子图像与模板图像的差异越小,即越相似.使得d值取得最小时的子图像Yk,即为目标图像Ykfinal,此时得到的基点(i,j)即为最佳配准点(ifinal,jfinal).

2)均方差度量[7].计算公式如下:

基于均方差度量的匹配算法的基本原理是将模板图像X与子图像Yk的灰度值作差,对每一项差值取平方后再取均值,即求得d值.d值越小说明子图像与模板图像的差异越小,即越相似.使得d值取得最小时的子图像Yk,即为目标图像Ykfinal,此时得到的基点(i,j)即为最佳配准点(ifinal,jfinal).

以上2种相似性度量是最小度量,采用这2种度量的匹配图像搜索过程如图5所示,搜索完成后输出最佳配准点(ifinal,jfinal)和度量值dmin.

图5 匹配图像搜索过程(采用最小度量) 图6 匹配图像搜索过程(采用最大度量)

3)归一化积相关度量[7].计算公式如下:

(3)

归一化积相关度量是最大度量,采用这种度量的匹配图像搜索过程如图6所示,搜索完成后输出最佳配准点(ifinal,jfinal)和度量值ρmax.

1.4 原始算法

本文所研究的光学分子影像手术导航系统中的图像融合原始算法利用相关系数(CorrelationCoefficient,CC)度量方法对灰度图像和荧光图像的子图像进行模板匹配搜索匹配图像.相关系数度量计算公式如下：

(4)

式(4)中：Cov(X,Yk) 表示X与Yk的协方差;D(X)与D(Yk)分别为X与Yk的方差.

图像的相关系数是两幅图像之间近似程度的一种线性描述.ρ值越接近于1,说明模板图像X和子图像Yk越近似于线性关系,即两幅图像的相似性越大.使ρ值取得最大时的子图像Yk即为目标图像Ykfinal,此时得到的基点即为最佳配准点(ifinal,jfinal).

1.5 图像融合

本文选用了像素级图像融合中的灰度加权平均法[9-11]进行图像的融合.首先利用灰度级-彩色变换技术[12]对荧光图像的荧光区域加伪绿,然后根据图像配准得到的最佳配准点将彩色图像与伪绿图像进行融合.本算法中只在荧光图像的荧光区域加伪绿,荧光图像中除了有荧光的区域外,其他区域的灰度值很小(或者为0),因此,当灰度值小于一定值时(此时肉眼看不到有荧光),将3个颜色分量的值都设为0.由于伪绿图像的颜色值除了RGB中的G分量外其他都为0,因此,融合算法中灰度值加权平均的方法简化为像素值直接相加,将伪绿直接融合到彩色图像上.

2 结 果

本文选用了32组图像数据,其中包括预临床实验中对小鼠、兔子乳腺癌手术的成像和临床实验中对人体乳腺癌手术的成像.对彩色图像重采样使运算时间减少到原来的1/16.另外,通过选取兴趣点将基点搜索区域从200×200(原始算法)缩小到100×100,进一步减少了计算量.

2.1 融合效果

图7(a)、图8(a)是荧光图像;图7(b)、图8(b)是彩色图像;图7(c)、图8(c)是荧光图像加伪绿后得到的以伪绿显示病灶区域的图像.图7、图8中的(d),(e),(f)分别是采用MAD,MSD,NPROD3种度量方法配准后对(a),(b),(c)3幅图像进行融合得到的融合图像,(g)是原始算法对(a),(b),(c) 3幅图像进行融合得到的融合图像,图中荧光部分显示了肿瘤的位置、大小等详细信息.分别比较图7、图8中的(d),(e),(f),(g)可以看出:图7中,(d),(e)的融合效果相同,都比较好;(f)的融合效果较差;(g)的融合效果最差;图8中,(e),(f)的融合效果相同,都比(d),(g)的融合效果好,其中(g)的融合效果最差.

图7 切开前的成像及其融合图像 图8 切开后的成像及其融合图像

表1列出了32组实验数据分别采用3种相似性度量的算法和原始算法得到的配准时间和最佳配准点，其中原始算法在处理第 18 组数据时配准失败在32组图像数据的处理结果中：

表1 采用3种相似性度量的算法与原始算法所用配准时间及得到的最佳配准点

注：表中每组图像数据得到的4对最佳配准点根据融合效果分别用实下划线、点式下划线、波浪线标出,融合效果按下划线的顺序依次递减,未标线的点融合效果最差.

1)采用平均绝对差度量与采用均方差度量相比,有11组融合效果相同或相近,有6组采用平均绝对差度量的算法融合效果好,有15组采用均方差度量的算法融合效果好;

2)采用均方差度量与采用归一化积相关度量相比,有26组融合效果相同或相近,有5组采用均方差度量的算法融合效果好,有1组采用归一化积相关度量的算法融合效果好;

3)采用归一化积相关度量与原始算法相比,有12组融合效果相同或相近,有20组采用归一化积相关度量的算法融合效果好;

4)采用平均绝对差度量与采用归一化积相关度量相比,有9组融合效果相同或相近,有10组采用平均绝对差度量的算法融合效果好,有13组采用归一化积相关度量的算法融合效果好;

5)采用平均绝对差度量与原始算法相比,有3组融合效果相同或相近,有21组采用平均绝对差度量的算法融合效果好,有8组原始算法融合效果好;

综合以上的分析比较可以得出:采用3种相似性度量的算法的融合效果都比原始算法的融合效果好,其中采用均方差度量的算法融合效果最好；采用归一化积相关度量的算法融合效果次之；采用平均绝对差度量的算法融合效果稍差.

2.2 配准时间

由表2可知,对本文中32组图像数据进行配准时,采用平均绝对差度量、均方差度量、归一化积相关度量的算法和原始算法的平均配准时间分别为43.200 6，46.685 3，85.978 6 和539.384 s.从图9中可以非常直观地看到采用3种相似性度量的算法与原始算法配准时间长短的对比.从算法的稳定性来看,采用平均绝对差度量的算法最稳定；采用均方差度量的算法次之；原始算法的稳定性最差.

表2 采用3种相似性度量的算法与原始算法所用配准时间的最小值、最大值、平均值以及标准差

便于作图,图中数据均已标准化

综合比较表1中3种相似性度量的算法与原始算法的融合效果,以及表2和图9中对配准时间的统计分析可以得出,本算法中,采用3种相似性度量的算法的匹配性能均比原始算法的匹配性能好,其中,均方差度量的匹配性能最好,平均绝对差次之,归一化积相关度量稍差.

3 讨论分析

本文对中国科学院分子影像重点实验室研制的光学分子影像手术导航系统中的图像融合算法进行了改进和研究,在提高算法运行效率、缩短图像融合时间的同时提高了图像的融合效果.配准技术仅利用图像的荧光信息及病人体表信息,就可以快速完成图像的配准融合过程.算法的平均配准时间从原始算法的539.384 s缩短到最短43.183 6 s,并且采用3种相似性度量的算法融合效果均比原始算法好.平均绝对差、均方差、归一化积和相关系数(原始算法)的时间复杂度均为O(m×n),其中TMAD=3m×n,TMSD=3m×n,TNPROD=6m×n,TCC=10m×n.计算平均绝对差只有加减运算,计算均方差有乘法运算,采用均方差度量的算法比采用平均绝对差的算法配准时间稍长.采用归一化积相关度量的算法配准时间是采用平均绝对差度量算法的2倍左右.本文在改进原始算法时利用图像的梯度信息将配准点的搜索范围从原始算法200×200缩小到100×100的区域中.因此,采用3种相似性度量的算法比原始算法采用相关系数度量时配准时间大幅缩减.另外,本文在缩短算法运行时间的同时,还提高了图像配准的精度和图像融合的效果.平均绝对差和均方差属于距离度量,归一化积相关和原始算法中的相关系数属于相关度量,特别是相关系数是度量2幅图像间的线性相关程度,而本文中分子影像手术导航系统采集的荧光图像和彩色图像的灰度相关性很小.因此，归一化积相关和相关系数在配准本文中的图像时匹配性能不好.本文还通过域法修正避免了因兴趣点偏离荧光图像中心太远而导致的配准失败.该算法在运行效率上还有进一步改进的可能,希望在未来的研究中有更大的突破.

4 结 论

针对术中肿瘤的精确定位这一国际上的研究热点及挑战性问题,本文对中国科学院分子影像重点实验室研制的光学分子影像手术导航系统中的图像融合算法进行了研究,该算法在图像配准的精度和运行时间以及图像融合的效果方面都得到了较好的结果.通过图像配准的时间和图像融合的效果,对分别采用平均绝对差、均方差和归一化积相关3种相似性度量的算法与原始算法的匹配性能进行比较得出,采用3种相似性度量的算法的匹配性能均比原始算法的匹配性能好,其中,采用均方差度量的算法匹配性能最好.

[1]赫捷,陈万青.2012中国肿瘤登记年[M].北京:军事医学科学出版社,2012.

[2]王志刚.导航术中手术器械的光学跟踪技术研究[D].广州：华南理工大学生物科学与工程学院,2012.

[3]石立兴,张继武.光学分子影像学及其应用[J].中国医学影像技术,2008,24(12):2024-2026.

[4]李伟.像素级图像融合方法及应用研究[D].广州:华南理工大学自动化科学与工程学院,2006.

[5]Zitova B,Flusser J.Image registration methods:A survey[J].Image and Vision Computing,2003,21(11):977-1000.

[6]廉蔺,李国辉,张军,等.基于边缘最优映射的红外和可见光图像自动配准算法[J].自动化学报,2012,38(4):570-581.

[7]Giachetti A.Matching techniques to compute image motion[J].Image and Vision Computing,2000,18(3):247-260.

[8]李庆鹏.多模态医学图像配准及基于小波变换的图像融合算法的研究[D].武汉：华中科技大学生命科学与技术学院,2007.

[9]Stokking R,Zubal I G,Viergever M A.Display of fused images:Methods,interpretation,and diagnostic improvements[J].Seminars in Nuclear Medicine,2003,33(3):219-227.

[10]李添捷.生物医学图像融合显示方法的研究[D].上海:复旦大学信息科学与工程学院,2012.

[11]胡钢,刘哲,徐小平,等.像素级图像融合技术的研究与进展[J].计算机应用研究,2008,25(3):650-655.

[12]李全越,王芳.伪彩色处理在医学图像中的应用[J].微计算机信息,2008,24(3):299-300.