医学图像配准中的数据抽样方法研究

光电子激光第19卷第2期2008年2月JOURNALOFOPTOELECTRONICSLASERVOL19NO2FEB2008医学图像配准中的数据抽样方法研究3胡顺波1,2,刘常春1331山东大学控制科学与工程学院,山东济南2500612临沂师范学院物理系,山东临沂276005摘要针对在基于灰度的医学图像配准中,传统数据抽样方法在过抽样时产生较多的局部极值点问题,提出了基于浮动图像灰度概率分布和其梯度信息的抽样方法。通过对人体脑部的刚体配准实验,从函数曲线和收敛性能方面,对比分析了五种数据抽样方法。实验结果表明,新抽样方法可以有效地减少局部极值点,提高归一化互信息测度的收敛性能。关键词图像配准归一化互信息抽样中图分类号TN91173文献标识码A文章编号10052008620080220237205DATASUBSAMPLINGINMEDICALIMAGEREGISTRATIONHUSHUN2BO1,2,LIUCHANG2CHUN1331SCHOOLOFCONTROLSCIENCEANDENGINEERING,SHANDONGUNIVERSITY,JINAN250061,CHINA2DEPARTMENTOFPHYS2ICS,LINYINORMALUNIVERSITY,LINYI276005,CHINAABSTRACTINORDERTODECREASETHENUMBEROFLOCALEXTREMACAUSEDBYOVERSUBSAMPLINGOFCONVENTIONALSUBSAMPLINGMETHODSININTENSITY2BASEDMEDICALIMAGEREGISTRATION,ANEWSUBSAMPLINGMETHODBASEDONINTENSITYPROBABILITYDISTRIBUTIONANDGRA2DIENTINFORMATIONOFFLOATINGIMAGEISPROPOSEDTHEFUNCTIONCURVESANDCONVERGENCEARECOMPAREDBYAPPLYINGFIVESUBSAM2PLINGMETHODSTOTHERIGIDREGISTRATIONOFBRAINIMAGESTHERESULTSOFTESTSSHOWTHATTHENEWSUBSAMPLINGMETHODOUTPER2FORMSOTHERSUBSAMPLINGMETHODSINLESSEXTREMAANDCONVERGENCEPERFORMANCEOFTHENORMALIZEDMUTUALINFORMATIONKEYWORDSIMAGEREGISTRATIONNORMALIZEDMUTUALINFORMATIONSUBSAMPLING1引言最大互信息配准是目前研究较多,应用较广的一种基于像素灰度的配准方法。在1995年被COLLIGNON等1和VIOLA等2应用到医学图像配准之后,它的有效性得到了广泛的认可36。这种方法配准精度高且不需要图像预处理,可以用于任何不同模态医学图像的配准。但是这种方法运算量大,配准时间长,难以满足临床上实时处理的要求。在基于像素灰度的配准方法中7,数据抽样是提高配准速度的一种有效手段。MAES等3使用了均匀抽样法,在浮动图像的坐标轴方向上进行均匀抽样。但是这种方法容易丢失一些重要信息,如图像的边缘等。SABUNCU等8提出了两种不同上述几种数据抽样方法都没有利用到灰度概率分布信息,而这是基于像素灰度配准中最重要的信息。为此本文提出一种新的抽样方法,这种方法综合了梯度信息和灰度概率分布信息,能保证浮动图像的灰度概率分布在抽样前后基本不变,对相同灰度级的像素点,优先抽样高梯度的像素点。通过不同抽样方法的对比实验表明,该方法可以提高归一化互信息测度的收敛性能,减少局部极值点。2方法21互信息和归一化互信息已知随即变量A和B,它们之间的联合概率分布为PAB,的非均匀抽样法确定性非均匀抽样DETERMINISTICNONUNIFORMSUBSAMPLING和随机性非均匀抽样STOCHASTICNONUNIFORMSUB2相应的边缘概率分布为A定义为A和PB,则A和B间的互信息PABA,BSAMPLING,这两种方法都是根据图像的梯度大小进行抽样,梯IA,BPABA,BLOGPAPB度大的像素抽样间隔小,梯度小的像素抽样间隔大。SABUNCUA,BAB的方法考虑了图像的边缘等重要信息,但是抽样前后图像数据其中的灰度概率分布变化太大,造成测度局部极值点过多,影响了配准。HAHBHA,B1HAPAALOGPAA2A3收稿日期20072062163基金项目国家“863”计划资助项目2006AA02Z4D9山东省自然科学基金资助项目Z2006C0533E2MAILHSBTIGER7748MAILSDUEDUCN3PB光电子激光2008年第19卷HA,BPABA,BLOGPABA,B3A,BA和B间的归一化互信息定义为4如果函数FX满足以下两点FX0,XR。FN1,NN。则两幅图像间的联合直方图可NNMIA,BHAHBHA,B很明显NMIA,B1IA,BHA,B45以通过以下方式更新HFX,Y,RIP,JQHFX,Y,RIP,JQFPIFQJ6B样条函数满足F函数的条件。当F是一阶B样条函归一化互信息可以有效地解决“互信息测度对交叠区域减少敏感”的问题4,本文用归一化互信息测度进行配准。22数据抽样数据抽样都是在浮动图像F上进行,抽样后的像素点的集合记为SF。浮动图像经过几何交换T后,SF在参考图像R上数时,GPVE插值法就是PV插值法。CHEN还验证了当F是二阶B样条函数或三阶B样条函数时,GPVE插值法可以有效地消除局部极值现象。正整数阶B样条函数可以由下面的公式7和8导出。的对应点的集合记为SR,SRTSF。集合SR中的点不一定在参考图像的整数网格上,这就需要插值,详见23节。本文提出的数据抽样方法的详细步骤如下NXN130XN1T0TXDT71输入浮动图F和抽样率SAMPLERATE。0X11/2X1/20其它82计算F的灰度概率分布PFA,A为F的任一灰度级。确定抽样后每个灰度级对应的像素个数NANFPFA其中,一阶形式B样条函数和二阶B样条函数写成如下的SAMPLERATE,NF为F中的像素总数。初始化DA0。DA统计集合SF中灰度值为A的像素个数。1X1X1X091X0X13计算F的梯度图像AF,对AF中的值按从小到大的顺序排序,得到Q。1X3223X1224按顺序扫描Q的所有像素,如果第I个像素的灰度值为A,且该像素不是噪声,则DADA1。如果DAN2X3X21422X1210A,则把像素I并入到集合SF中。5输出集合SF。1X3221X322步骤2和4保证了抽样前后数据集合的灰度概率分布稳定性。因为高梯度可能是噪声引起的,所以步骤4不统计这些噪声引起的梯度点。步骤4按梯度的大小顺序进行抽样,可以采集到重要的边缘信息。我们称这种抽样方法为基于概率分布和梯度的抽样方法SUBSAMPLINGBASEDONPROBABILITYDISTRIBUTIONANDGRADIENT,简称为PDGS抽样法。23插值把公式9代入公式6得1阶B样条GPVE插值法,即PV插值法。把公式10代入公式6得二阶B样条GPVE插值法。根据公式6,一一对应扫描SF和SR中的所有像素,得两幅图像抽样后数据间的联合直方图HTA,B,A为SF中的灰度级,B为R中的灰度级。24测度和最佳几何变换计算根据两幅图像抽样后数据的联合直方图HTA,B,计算联合概率分布和边缘概率分布。HTA,BHUA2MEICHEN等9提出了GPVEGENERALIZEDPARTIALVOL2PTA,BHTA,B11UMEESTIMATION插值方法,解决了传统插值方法导致的局部极值问题。下面简述GPVE方法。A,BPTAPTA,B12B设X,Y是SF中一个像素点的坐标,经过几何变换T,X,Y对应R中点的坐标为II,JJ。图1是几何变换过程PBTPTA,B13A的示意图。把公式11、12和13代入到公式4,可得抽样后图像数据的归一化互信息测度。NMITPTALOGPTAPTBLOGPTBABPTA,BLOGPTA,BA,B14238图1几何变换过程示意图FIG1ILLUSTRATIONOFTHETRANSFORMATIONT在图像完全配准时,几何变换T称为最佳几何变换T3,此时的测度NMI的值达到最大T3ARGMAXNMIT153实验结果与分析我们通过实验对比分析了5种数据抽样方法,分别为一2第2期胡顺波等医学图像配准中的数据抽样方法研究维数据均匀抽样法UNL、MAES使用的均匀抽样法UN2、SABUNCU提出的确定性非均匀抽样法DNS、SABUNCU提出的随机性非均匀抽样法SNS、本文新提出的抽样法PDGS。UNL抽样法是把二维数据变成一维数据后,再进行均匀抽样。实验数据来源于网站10。实验平台配置是28GHZCPU,1G内存,操作系统为WINDOWSXP。31抽样结果对比图2描述了5种不同抽样方法的抽样结果,交叉号的中心表示抽样点。由图2可见,B图和C图是均匀抽样,容易丢失重要的边缘信息D图和E图的非均匀抽样结果和SABUN2CU的分析一致,D图的确定性抽样把梯度值分成10类进行抽样,抽样间隔与梯度类别值成反比,E图对应随机性非均匀抽样,这种方法不抽样梯度值为0的点,故在梯度值为0的背景区域无抽样点F图对应本文提出的方法,即能抽样重要的边缘轮廓信息,又能抽样到较高梯度值的背景区域边缘轮廓附近,从而保证了抽样后灰度概率分布的稳定性,又保证不丢失图2不同抽样方法的结果对比FIG2THECONTRASTOFDIFFERENTSUBSAMPLINGMETHODS重要信息。32插值和抽样对测度的影响把图2A作为参考图像,以1/16抽样率抽样后的SF作为浮动图像。令浮动图像沿Y轴保持不变,在X轴方向的平移范围是55,平移的步长为01像素,对每种抽样方法,计算归一化互信息测度值,得到图3BF。图3A对应浮动图像没有抽样的结果。图3比较了两种插值方法PV法和二阶B样条GPVE法。从图3可见二阶B样条GPVE法得到的测度曲线黑实线比PV插值法得到的曲线光滑蓝虚线,PV插值法在整数平移处有局部极大值二阶B样条GPVE法得到的测度幅值比PV插值法小,这是由于“支撑域增大”对联合直方图的平均效应造成的在使用二阶B样条GPVE法时,无抽样浮动图像对应的测度曲线最光滑,PDGS抽样法的曲线最接近A图,优于UN1、UN2、DNS、SNS插值法得到的曲线,前者有一个唯一的极大值点,后四者的极大值数目分别为2、3、4、6。33抽样率对测度的影响分别以1/4、1/16、1/64的抽样率得到SF,浮动图像的几何变换和32节中的相同,插值方法采用二阶B样条GPVE法,得到测度值,并画出曲线图4。图4比较了三种抽样方法UN1图3不同抽样方法时归一化互信息测度对X轴平移的函数曲线FIG3FUNCTIONSOFX2SHIFTSOFNORMALIZEDMUTUALINFORMATIONMEASUREUNDERDIFFERENTSUBSAMPLINGMETHODS239光电子激光2008年第19卷图4抽样率对归一化互信息测度的影响FIG4THEINFLUENCEOFSUBSAMPLINGRATEONNORMALIZEDMUTUALINFORMATIONMEASURE蓝虚线、SNS红点划线、PDGS黑实线。UN2的结果类似UN1,DNS的结果类似SNS。从图4可见,抽样率越低,抽样后用于配准的数据越少,测度曲线越不光滑,抽样率越高,测度曲线越光滑在同样抽样率的情况下,PDGS法得到的曲线最光滑,其次是UN1法,最后是SNS法在抽样率为1/4,各抽样方法都有唯一的极值点,这说明在抽样率很高且使用二阶B样条GPVE插值法时,各种抽样方法都不会产生伪局部极值点。34抽样方法对测度收敛性能的影响本文通过两个指标来表示收敛性能,分别是AFAAREAOFFUNCTIONATTRACTION和MAX局部极大值数。CAPEK等人11首先使用AFA来描述测度函数的收敛性能。AFA计算从某初始位置通过爬山法能收敛到全局极大值的体素数。本文爬山法的步长是01像素。AFA的值越大,测度的收敛范围越广。MAX统计某一区域范围内极大值点的数目。MAX的值越大,测度曲线越不光滑。参考图像和SF的选择类似32节,抽样率为1/4、1/16、1/64,插值方法为二阶B样条GPVE法。选定以下平移范围X55、Y55,平移步长为01像素,计算不同抽样方法对应的归一化互信息测度的AFA和MAX。结果见表1。表1中数据表明,在较低抽样率的情况下,对于归一化互信息测度,PDGS方法的收敛性能优于其它四种抽样方法,因为它有最小的MAX,最大的AFA随着抽样率的增大,各抽样方法的收敛性能都提高当抽样率增大到1/4时,PDGS和UN2抽样法在整个平移范围内都收敛,UN1和DNS的收敛区域较大,SNS的收敛性能最差。表1不同抽样方法下测度收敛性能比较TAB1CONVERGENCEUSINGDIFFERENTSUBSAMPLINGMETHODSSUBSAMPLINGRATE列出抽样率为1/16时PDGS和SNS插值法对应的曲面。由图5可见,PDGS插值法得到的曲面比SNS插值法得到的曲面光滑,故其MAX小于SNS的MAXPDGS插值法得到的曲面中间的“帽子”粗于SNS,故其AFA大于SNS的AFA。1/41/161/64AFA1161132710905412136图5是表1一些配准数据对应的曲面,由于篇幅所限,仅定性,又不丢失重要的边缘轮廓信息。在抽样率较低的情况下,使用二阶B样条GPVE插值法,利用抽样后的图像结果进240MEASURESUN1UN2DNSSNSPDGS图5不同抽样方法对应的引力盆曲线MAX214591FIG5ATTRACTIONBASINOFDIFFERENTSUBSAMPLINGMETHODSAFA100821020199283822102014结论MAX4661549038AFA26422477192810882799根据浮动图像的灰度概率分布和梯度信息,按照本文提出MAX7669779550的PDGS抽样方法进行数据抽样,即保证了灰度概率分布的稳第2期胡顺波等医学图像配准中的数据抽样方法研究行归一化互信息配准,相对于两种均匀抽样法和两种基于梯度的非均匀抽样法,PDGS抽样方法减少了局部极值现象,提高了6LIUCHANG2CHUN,YANGJIN2BAO,HUSHUN2BOMEDICALIMAGEREGISTRA2TIONBASEDONIMPROVEDLA2INFORMATIONJJOURNALOFOPTOELECTRONICS测度的收敛性能。LASER光电子激光,2007,1810123121234INCHINESE参考文献7LIUGUI2XI,WANGLEIANIMAGEREGISTRATIONMETHODBASEDONREGIONSELECTINGANDFEATUTEPOINTSMATCHINGJJOURNALOFOPTOELECTRONICS1COLLIGNONA,MAESF,DELAERED,ETALAUTOMATEDMULTIMODALITYMEDI2LASER光电子激光,2007NCHINESECALIMAGEREGISTRATIONUSINGINFORMATIONTHEORYCINFORMATIONPRO28SABUNCUMR,RAMADGEPJGRADIENTBASEDNONUNIFORMSUBSAMPLINGCESSINGINMEDICALIMAGING,DORDRECHTKLUWERACAEMICPUBLISHERS,FORINFORMATION2THEORETICALIGNMENTMETHODSCPROCEEDINGOFTHE1995,263227426THANNUALINTERNATIONALCONFERENCEOFTHEIEEEEMBS,SANFRANCIS22VIOLAP,WELLSWMALIGNMENTBYMAXIMIZATIONOFMUTUALINFORMATIONCPROCEEDINGSOFTHE5THINTERNATINALCONFERENCEONCOMPUTERVI2SION,WASHINGTONDCIEEECOMPUTERSOCIETY,1995,162233MAESF,COLLIGNONA,VANDERMEULEND,ETALMULTIMODALITYIMAGEREG2ISTRATIONBYMAXIMIZATIONOFMUTUALINFORMATIONJIEEETRANSACTIONSONMEDICALIMAGI