本发明属于图像处理与分析中的图像配准领域,更具体地,涉及一种非刚性多模医学图像的配准方法及系统。
背景技术:
非刚性多模医学图像配准对于医学图像分析和临床研究非常重要。由于各种成像技术的原理不同,在反映人体的信息方面各有其优势。计算机断层扫描(ct),磁共振成像(mri)和超声成像可以显示器官的解剖信息。正电子发射断层扫描(positronemissiontomography,pet)作为一种功能成像方式,可以显示代谢信息,但不能清楚地提供器官的解剖信息。多模态图像融合技术可以结合不同模态图像的信息,从而得到更加精准的诊断和更好治疗。
图像配准的目的是为了寻找图像中相应结构之间的正确的空间对应关系,是有效图像融合的前提。目前已经得到广泛的研究。例如,将经直肠超声图像与术前mr图像进行配准,用于引导前列腺穿刺活检手术的进行。在癫痫治疗中,mri和pet图像配准,以帮助识别脑组织的功能性区域,并指导电极的放置。
传统的解决非刚性多模图像配准问题的方法大体上被分为两类:第一类是基于互信息测度的配准方法,然而这类方法通常都没有考虑图像的局部特征结构,计算耗时,而且容易陷入局部极值导致配准结果不准确。第二类方法通过图像结构表征方法将多模图像配准简化为单模图像配准,如通过熵图、韦伯局部特征描述子(weberlocaldescriptor,wld)及基于模态独立邻域描述子(modalityindependentneighborhooddescriptor,mind)等特征对图像结构进行表征,然后利用表征结果的差值平方和(sumofsquareddifference,ssd)作为配准测度来实现图像配准。其中,基于熵图的方法是通过估计图像块的灰度概率密度函数计算得到每个图像块的熵值,由此获得整幅图像的熵图,而基于wld的方法则是利用图像的拉普拉斯运算来描述其局部结构特征。这两种方法能有效克服多模图像间灰度差异的不利影响,但是基于熵图和wld的特征都对图像中的噪声敏感,因而在图像中存在噪声的情况下,难以产生精确的配准结果。基于mind的配准方法通过图像块间的欧式距离来评价图像自相似性,利用不同图像块间的相似性之和来刻画图像局部结构特征,但该方法在评价图像自相似性时仅考虑了图像块之间的平移不变性,忽略了图像块间可能存在的旋转特性,因此在图像块间存在旋转的情况下,该方法难以提供精确的配准结果。
最近,深度学习被用来实现两种模态的图像配准。第一种方法是使用深度学习来提取图像特征,并将特征用传统的配准方法来得到配准图像。这些方法主要是基于cnn和sae等复杂的网络结构,这些方法的学习速度很慢,并且cnn和sae中的参数选择不当可能会导致局部最优。第二种方法是利用深度学习,通过直接从输入图像中学习变形场来实现端对端的图像配准。例如,利用标签驱动的弱监督学习来实现多模非刚性图像配准。但是,医学图像中的金标准数据是稀缺的。hessam等人提出的基于cnn的配准网络(regnet)可以直接从一对输入图像中估计位移矢量,这些方法多采用有监督的学习方法,但是由于使用传统的配准方法获得的形变场不精确,使用这些形变场进行学习本身就有一定的误差,其结果比传统的b样条配准方法差。
技术实现要素:
针对现有技术的以上缺陷或改进需求,本发明提供了一种基于深度学习的多模医学图像非刚性配准方法及系统,由此解决现有非刚性多模医学图像的配准精确性较低的技术问题。
为实现上述目的,按照本发明的一个方面,提供了一种基于深度学习的多模医学图像非刚性配准方法,包括:
(1)将参考图像输入到已训练的两层pcanet网络中,获得所述参考图像各级的图像特征,并将所述参考图像各级的图像特征进行合成得到参考图像的结构表征图,以及,将浮动图像输入到所述已训练的两层pcanet网络中,获得所述浮动图像各级的图像特征,并将所述浮动图像各级的图像特征进行合成得到所述浮动图像的结构表征图;
(2)根据所述参考图像的结构表征图以及所述浮动图像的结构表征图建立目标函数,根据所述目标函数获得变换参数,基于所述变换参数变换所述浮动图像,并对变换后的浮动图像进行插值处理,获得配准图像。
优选地,在步骤(1)之前,所述方法还包括:
对于n幅医学图像中的每幅图像的每个像素,无间隔的取k1×k2的块,将得到的所有块向量化,并将得到的所有向量进行组合得到目标矩阵,计算所述目标矩阵的特征向量,并将所述目标矩阵的特征值按从大到小进行排序,将前l1个特征值对应的特征向量进行矩阵化,得到第一层pcanet的l1个卷积模板;
将各卷积模板分别与输入图像进行卷积得到nl1幅图像,将所述nl1幅图像输入到第二层pcanet中,得到第二层pcanet的l2个卷积模板,并得到nl1l2幅图像,其中,k1×k2表示块的大小,l1为第一层pcanet选取的特征个数,l2为第二层pcanet选取的特征个数。
优选地,步骤(1)包括:
(1.1)将所述参考图像输入到所述已训练的两层pcanet网络中,得到所述参考图像的第一级图像特征
(1.2)由
优选地,步骤(1.1)包括:
(1.1.1)由
(1.1.2)由
(1.1.3)由
优选地,步骤(2)包括:
(2.1)由g(tτ)=ssd+αr(tτ)建立目标函数,其中,ssd表示结构表征图
(2.2)根据初始变换参数tτ变换所述浮动图像的结构表征图
(2.3)若迭代次数τ大于等于迭代次数阈值,且g(tτ)≤g(tτ-1),则根据最终的变换参数变换所述浮动图像,并对变换后的浮动图像进行插值处理,获得所述配准图像,否则返回步骤(2.2)。
优选地,相似性测度ssd为:
按照本发明的另一方面,提供了一种基于深度学习的多模医学图像非刚性配准系统,包括:
结构表征图构建模块,用于将参考图像输入到已训练的两层pcanet网络中,获得所述参考图像各级的图像特征,并将所述参考图像各级的图像特征进行合成得到参考图像的结构表征图,以及,将浮动图像输入到所述已训练的两层pcanet网络中,获得所述浮动图像各级的图像特征,并将所述浮动图像各级的图像特征进行合成得到所述浮动图像的结构表征图;
配准迭代模块,用于根据所述参考图像的结构表征图以及所述浮动图像的结构表征图建立目标函数,根据所述目标函数获得变换参数,基于所述变换参数变换所述浮动图像,并对变换后的浮动图像进行插值处理,获得配准图像。
优选地,所述系统还包括:pcanet训练模块;
所述pcanet训练模块,用于对于n幅医学图像中的每幅图像的每个像素,无间隔的取k1×k2的块,将得到的所有块向量化,并将得到的所有向量进行组合得到目标矩阵,计算所述目标矩阵的特征向量,并将所述目标矩阵的特征值按从大到小进行排序,将前l1个特征值对应的特征向量进行矩阵化,得到第一层pcanet的l1个卷积模板;将各卷积模板分别与输入图像进行卷积得到nl1幅图像,将所述nl1幅图像输入到第二层pcanet中,得到第二层pcanet的l2个卷积模板,并得到nl1l2幅图像,其中,k1×k2表示块的大小,l1为第一层pcanet选取的特征个数,l2为第二层pcanet选取的特征个数。
优选地,所述结构表征图构建模块包括:第一层有效特征构建模块以及第二层有效特征构建模块;
所述第一层有效特征构建模块,用于基于所述已训练的两层pcanet网络中,得到所述参考图像的第一级图像特征,以及,基于所述已训练的两层pcanet网络中,得到所述浮动图像的第一级图像特征;
所述第二层有效特征构建模块,用于基于所述已训练的两层pcanet网络中,得到所述参考图像的第二级图像特征,以及,基于所述已训练的两层pcanet网络中,得到所述浮动图像的第二级图像特征。
优选地,所述配准模块包括求解模块和判定模块;
所述求解模块,用于得到参考图像的结构表征图和浮动图像的结构表征图之间的相似性度量,并加上正则化项构建目标函数,获得变换参数;
所述判定模块,用于判断所述目标函数是否符合迭代停止标准,如果不符合迭代停止标准,则根据变换参数变换浮动图像的结构表征图,对变换后的浮动图像的结构表征图进行插值处理,并更新浮动图像的原始结构表征图,直至满足迭代停止标准,最终将所得空间变换作用于浮动图像,获得配准图像。
总体而言,通过本发明所构思的以上技术方案与现有技术相比,能够取得下列有益效果:
首先,与传统的配准方法相比,传统的配准方法往往只利用图像的灰度信息和一阶信息,而该方法利用深度学习,从大量的数据中提取图像的多级特征,有效表征复杂医学图像,为多模图像相似性的准确评价提供了有效依据。其次,与现有的用于配准的深度学习方法相比,该方法网络结构简单,易于训练,并且属于无监督学习,不需要标签数据,大大克服了医学图像缺少标签的问题。
附图说明
图1为本发明实施例提供的一种非刚性多模医学图像的配准方法的流程示意图;
图2为本发明实施例提供的另一种非刚性多模医学图像的配准方法的流程示意图;
图3(a)为本发明实施例以及对比例2-5所用的参考图像t1;
图3(b)为本发明实施例以及对比例2-5所用的浮动图像t2;
图3(c)为本发明实施例以及对比例2-5所用的浮动图像pd;
图3(d)为本发明实施例方法获得的配准图像t1-t2;
图3(e)为本发明对比例1方法获得的配准图像t1-t2;
图3(f)为本发明对比例2方法获得的配准图像t1-t2;
图3(g)为本发明对比例3方法获得的配准图像t1-t2;
图3(h)为本发明对比例4方法获得的配准图像t1-t2;
图3(i)为本发明实施例方法获得的配准图像gad-t2;
图3(j)为本发明对比例1方法获得的配准图像gad-t2;
图3(k)为本发明对比例2方法获得的配准图像gad-t2;
图3(l)为本发明对比例3方法获得的配准图像gad-t2;
图3(m)为本发明对比例4方法获得的配准图像gad-t2;
图4(a)为本发明实施例以及对比例1-4所用的参考图像pd加权mr;
图4(b)为本发明实施例以及对比例1-4所用的浮动图像ct;
图4(c)为本发明实施例方法获得的配准图像;
图4(d)为本发明对比例1方法获得的配准图像;
图4(e)为本发明对比例2方法获得的配准图像;
图4(f)为本发明对比例3方法获得的配准图像。
图4(g)为本发明对比例4方法获得的配准图像。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白,以下结合附图及实施例,对本发明进行进一步详细说明。应当理解,此处所描述的具体实施例仅仅用以解释本发明,并不用于限定本发明。此外,下面所描述的本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。
本发明提供了一种基于深度学习的多模医学图像非刚性配准方法及系统,通过大量的医学数据训练pcanet,并将参考图像以及浮动图像输入到训练好的pcanet中,分别获得参考图像以及浮动图像的结构表征图,从而实现非刚性多模医学图像的配准。
如图1所示为本发明实施例提供的一种非刚性多模医学图像的配准方法的流程示意图,包括:
(1)将参考图像输入到已训练的两层pcanet网络中,获得参考图像各级的图像特征,并将参考图像各级的图像特征进行合成得到参考图像的结构表征图,以及,将浮动图像输入到已训练的两层pcanet网络中,获得浮动图像各级的图像特征,并将浮动图像各级的图像特征进行合成得到浮动图像的结构表征图;
(2)根据参考图像的结构表征图以及浮动图像的结构表征图建立目标函数,根据目标函数获得变换参数,基于变换参数变换浮动图像,并对变换后的浮动图像进行插值处理,获得配准图像。
在本发明实施例中,如图2所示,在步骤(1)之前,该方法还包括对两层pcanet网络的训练过程:
对于n幅医学图像
将第一层得到的各卷积模板分别与输入的n幅图像进行卷积得到nl1幅图像,将nl1幅图像输入到第二层pcanet网络中,按照第一层的处理方式,通过对nl1幅图像中的每幅图像取块、向量化、计算特征向量得到第二层pcanet网络的l2个卷积模板,并将第二层得到的各卷积模板分别与输入的nl1幅图像进行卷积得到nl1l2幅图像,其中,k1×k2表示块的大小,一般的,当图像的复杂度较高时,优选取3×3,5×5等较小的值,配准效果会更好。l1为第一层pcanet网络选取的特征个数,l2为第二层pcanet网络选取的特征个数,l1与l2的值可以根据实际需要进行确定,优选地,取l1=l2=8。
在本发明实施例中,步骤(1)包括:
(1.1)将参考图像输入到已训练的两层pcanet网络中,得到参考图像的第一级图像特征
其中,步骤(1.1)的具体实现方式为:
(1.1.1)由
(1.1.2)由
(1.1.3)由
(1.2)根据参考图像以及浮动图像的第一层的有效特征和第二层的有效特征获得参考图像和浮动图像的结构表征图;
具体地,由
其中,
在本发明实施例中,以变换模型为基于b样条自由形变模型(free-formdeformation,ffd)构建目标函数为例,说明该配准过程,步骤(2)具体包括:
(2.1)由g(tτ)=ssd+αr(tτ)建立目标函数,g(tτ)的值越小,则浮动图像与参考图像的相似性越大,当g(tτ)最小时,完成配准,其中,ssd表示结构表征图
其中,相似性测度ssd为:
正则化项表示为:
(2.2)根据初始变换参数tτ变换浮动图像的结构表征图
(2.3)若迭代次数τ大于等于迭代次数阈值,且g(tτ)≤g(tτ-1),则根据最终的变换参数变换所述浮动图像,并对变换后的浮动图像进行插值处理,获得所述配准图像,否则返回步骤(2.2)。
其中,迭代次数阈值可以根据实际需要进行确定,在本发明实施例中不做唯一性限定。
在本发明实施例中,迭代求解的方法可以为限域拟牛顿法或梯度下降法,根据变换参数变换浮动图像采用的变换模型可以为基于b样条的自由形变模型,插值处理的方法可以为双线性插值法或b样条插值法,变换参数可以为三阶b样条函数,或者其它可以实现本发明的方法,具体采用何种方式,本发明实施例不做唯一性限定。
以下结合具体实施例对本发明方法进行详细说明。
实施例1
步骤1训练pcanet网络。输入n幅医学图像
步骤2根据pcanet深度学习网络得到一种基于pcanet的结构表征图(pcanetbasedstructuralrepresentation,简称psr)。
步骤2-1将参考图像和浮动图像输入到训练好的两层pcanet网络结构中,将得到参考图像和浮动图像的第一级特征信息
对于第一层信息:
其中,
对于第二层信息,共分两步处理,第一步,利用sigmoid函数和加权操作将l1l2信息合成l1信息
其中,
第二步,按照构建第一层有效特征的方式,构建第二层的有效特征:
步骤2-2计算图像i结构表征图psr,公式如下:
其中,
由于在本实施例中,分别为:
其中,mean(·)为均值运算符,c1和c2为调节系数,在本实施例中,c1=c2=0.8。
根据上述公式,可分别获得参考图像的结构表征图
步骤3-1采用基于b样条自由形变模型(free-formdeformation,ffd)作为变换模型,建立目标函数g(t)=ssd+αr(t)。其中,ssd表示结构表征图
步骤3-2以g(t)的值最小为目标,使用限域拟牛顿法(l-bfgs)进行迭代求解,当g(t)获得最小值时迭代停止,获得变换参数t;
步骤3-3通过步骤3-2中获得的变换参数t,对浮动图像变形,并通过双线性插值法进行插值,获得浮动图像对应的配准图像,并最终完成图像配准。
对比例1
按照(med.imageanal.16(7)(2012)1423-1435.)里的mind方法实现配准。具体参数为:图像块大小选择3×3。
对比例2
按照(med.imageanal.16(1)(2012)1-17.)里的essd方法实现配准。其中,具体参数为:选择7×7的图像块,并利用高斯权重、局部归一化方法和parzen窗估计来计算图像块对应的熵,由此获得整个图像对应的essd。
对比例3
按照(patternrecognit.32(1)(1999)71-86.)里的nmi方法实现配准。
对比例4
按照(sensors13(6)(2013)7599-7613.)里的wld方法实现配准。具体参数为:wld的计算选择半径r=1和r=2,构建相似性测度的块大小为7×7,权重项γ=0.01。
结果分析
为了进一步体现本发明的优点,我们将实施例1与对比例1-4的配准精度进行比较。配准精度采用目标配准误差tre进行评价,这里tre定义为:
其中ts表示随机形变,也是评价的金标准,tr表示由配准算法得到的形变,n表示用于图像配准性能评价的像素个数。
采用仿真mr图像进行配准精度测试,实施例1所用的仿真t1、t2和pd加权mr图像均取自brainweb数据库,表1列出了各算法所得tre的标准差和均值。从表1可看出,对不同mr图像进行配准时,实施例1提供的tre其均值和标准差皆低于其它方法,这说明本发明提出的方法在所有比较的方法中具有最高的配准精度。
表1各方法在t1-t2,pd-t2和t1-pd图像配准时的tre(mm)对比
为更直观地显示本发明相对于其余方法的优越性,我们提供了实施例与对比例1-4对应配准图像的视觉效果图,如图3所示。图3(a)为参考图像t1,图3(b)为浮动图像t2,图3(c)为浮动图像pd,图3(d)为实施例1方法获得的配准图像t2-t1,图3(e)为对比例1方法获得的配准图像t2-t1,图3(f)为对比例2方法获得的配准图像t2-t1,图3(g)为对比例3方法获得的配准图像t2-t1,图3(h)为对比例4方法获得的配准图像t2-t1,图3(i)为实施例方法获得的配准图像pd-t1,图3(j)为对比例1方法获得的配准图像pd-t1,图3(k)为对比例2方法获得的配准图像pd-t1,图3(l)为对比例3方法获得的配准图像pd-t1,图3(m)为对比例4方法获得的配准图像pd-t1,由配准结果图可以看出,我们的配准结果比其它方法更好。
针对ct与mr图像进行了配准精度的比较,我们对浮动ct图像进行了五种随机的变形处理。图4是实施例以及对比例1-4的方法对真实ct和mr图像进行配准的结果。其中,图4(a)为参考图像pd加权mr,图4(b)为浮动图像ct,图4(c)为实施例1方法获得的配准图像,图4(d)为对比例1方法获得的配准图像,图4(e)为对比例2方法获得的配准图像,图4(f)为对比例3方法获得的配准图像,图4(g)为对比例4方法获得的配准图像。可看出,对比例3-4的方法难以有效更正轮廓的变形,对比例1-2提供的配准图像中最外部轮廓局部出现失真,而实施例1方法可获得良好的配准图像。表3显示了每组真实ct图像与mr图像配准的tre的平均值,其中,组3图像即为图4中所用的真实ct图像与mr图像。
表2各方法在ct-mr图像配准时的tre(mm)对比
从表2中我们可以看出,与其它配准方法相比,本发明的方法对组1~组5的ct-mr图像进行配准皆可取得更低的tre,这说明本发明提出的方法在ct-mr图像配准方面较对比例的算法具有更高的配准精度。
本发明还提供了一种基于深度学习的多模医学图像非刚性配准系统,包括:
结构表征图构建模块,用于将参考图像输入到已训练的两层pcanet网络中,获得参考图像各级的图像特征,并将参考图像各级的图像特征进行合成得到参考图像的结构表征图,以及,将浮动图像输入到已训练的两层pcanet网络中,获得浮动图像各级的图像特征,并将浮动图像各级的图像特征进行合成得到浮动图像的结构表征图;
配准迭代模块,用于根据参考图像的结构表征图以及浮动图像的结构表征图建立目标函数,根据目标函数获得变换参数,基于变换参数变换浮动图像,并对变换后的浮动图像进行插值处理,获得配准图像。
在一个可选的实施方式中,该系统还包括:pcanet训练模块;
pcanet训练模块,用于对于n幅医学图像中的每幅图像的每个像素,无间隔的取k1×k2的块,将得到的所有块向量化,并将得到的所有向量进行组合得到目标矩阵,计算目标矩阵的特征向量,并将目标矩阵的特征值按从大到小进行排序,将前l1个特征值对应的特征向量进行矩阵化,得到第一层pcanet的l1个卷积模板;将各卷积模板分别与输入图像进行卷积得到nl1幅图像,将所述nl1幅图像输入到第二层pcanet中,得到第二层pcanet的l2个卷积模板,并得到nl1l2幅图像,其中,k1×k2表示块的大小,l1为第一层pcanet选取的特征个数,l2为第二层pcanet选取的特征个数。
在一个可选的实施方式中,结构表征图构建模块包括:第一层有效特征构建模块以及第二层有效特征构建模块;
第一层有效特征构建模块,用于基于已训练的两层pcanet网络中,得到参考图像的第一级图像特征,以及,基于已训练的两层pcanet网络中,得到浮动图像的第一级图像特征;
第二层有效特征构建模块,用于基于已训练的两层pcanet网络中,得到参考图像的第二级图像特征,以及,基于已训练的两层pcanet网络中,得到浮动图像的第二级图像特征。
在一个可选的实施方式中,配准模块包括求解模块和判定模块;
求解模块,用于得到参考图像的结构表征图和浮动图像的结构表征图之间的相似性度量,并加上正则化项构建目标函数,获得变换参数;
判定模块,用于判断目标函数是否符合迭代停止标准,如果不符合迭代停止标准,则根据变换参数变换浮动图像的结构表征图,对变换后的浮动图像的结构表征图进行插值处理,并更新浮动图像的原始结构表征图,直至满足迭代停止标准,最终将所得空间变换作用于浮动图像,获得配准图像。
本领域的技术人员容易理解,以上所述仅为本发明的较佳实施例而已,并不用以限制本发明,凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等,均应包含在本发明的保护范围之内。