CN104766304B - 一种基于多序列医学图像的血管配准方法 - Google Patents

一种基于多序列医学图像的血管配准方法 Download PDF

Info

Publication number
CN104766304B
CN104766304B CN201510088135.2A CN201510088135A CN104766304B CN 104766304 B CN104766304 B CN 104766304B CN 201510088135 A CN201510088135 A CN 201510088135A CN 104766304 B CN104766304 B CN 104766304B
Authority
CN
China
Prior art keywords
image
point
wavelet
vessel
multisequencing
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201510088135.2A
Other languages
English (en)
Other versions
CN104766304A (zh
Inventor
汪晓妍
李军伟
黄晓洁
张剑华
滕忠照
陈胜勇
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhejiang University of Technology ZJUT
Original Assignee
Zhejiang University of Technology ZJUT
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Zhejiang University of Technology ZJUT filed Critical Zhejiang University of Technology ZJUT
Priority to CN201510088135.2A priority Critical patent/CN104766304B/zh
Publication of CN104766304A publication Critical patent/CN104766304A/zh
Application granted granted Critical
Publication of CN104766304B publication Critical patent/CN104766304B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Landscapes

  • Image Processing (AREA)
  • Image Analysis (AREA)

Abstract

一种基于多序列医学图像的血管配准方法,包括以下步骤:1)获取多序列医学血管图像,选定一个序列图像为配准的参考图像,其他序列图像作为浮动图像;2)对所述多序列医学血管图像进行去噪处理;3)血管边缘的shape context描述;4)边缘点匹配;5)滤除误匹配;6)边缘校正与插值:根据匹配点对血管壁做迭代演化,使形变后的参考图像与浮动图像血管轮廓能量函数E(pi,qj)达到最小,根据边缘校正后的结果,采用样条插值的方法对浮动图像的局部血管区域做插值,最终得到配准后的结果。本发明提供一种精度较高的基于多序列医学图像的血管配准方法。

Description

一种基于多序列医学图像的血管配准方法
技术领域
本发明属于医学图像处理技术领域,尤其是涉及到多序列血管图像的配准方法。
背景技术
当前在MR血管图像配准领域主要是采用传统的自然图像配准方法,即建立一种图像相似度量函数或者误差函数,通过一定的优化策略使得参考图像与浮动图像的相似度达到最大或者误差最小,最终得到参考图像与浮动图像的变换矩阵,然而这些方法没有考虑到MR血管图像的结构复杂,对比度、信噪比(SNR)较低等情况,使得相似度度量函数或者误差函数容易陷入局部极大或者极小值,而达不到全局精确配准的目的。另外人体血管是一种非刚性组织,其形态会受到肌肉收缩,人体颤抖,脉动等因素的影响,这些因素使得一些面向刚形体配准的方法和有限自由度的配准模型得不到理想结果。当前基于像素的配准方法利用参考图像和浮动图像所包含的信息,直接利用图像的灰度值做配准,通过优化图像的互信息(MI)、归一化互信息(NMI)、联合熵等手段获取图像之间的变换矩阵,但是在对比度较低的医学图像对这些度量信息的并不敏感,极容易使得图像陷入局部极大值而达不到配准精度。基于特征点集(不连续点,转折点,交叉点,线交叉点,角点等)的配准方法首先要精确的提取图像的特征点,然后做精确匹配,而医学图像的收到噪声的干扰,极容易提取到伪特征点而错误匹配,从而影响变换参数精度性。
发明内容
为了克服已有医学图像处理技术的血管配准精度较低的不足,本发明提供一种精度较高的基于多序列医学图像的血管配准方法。
本发明解决其技术问题所采用的技术方案是:
一种基于多序列医学图像的血管配准方法,所述配准方法包括以下步骤:
1)获取多序列医学血管图像,选定一个序列图像为配准的参考图像,其他序列图像作为浮动图像;
2)对所述多序列医学血管图像进行去噪处理;
3)血管边缘的shape context描述,血管内外壁分别记为Clumen(X),Cwall(X),X代表边缘点坐标,边缘上每一点分别用Shape Contex描述子统计其直方图;
4)边缘点匹配
根据参考图像血管边缘点与浮动图像血管边缘点的代价函数Cij
其中,pi和qj分别代表颈动脉血管参考图像序列和浮动图像序列的边缘点集,hi(k)、hj(k)分别代表参考图像边缘点直方图和浮动图像边缘点直方图,K代表直方图总共的等级数,k是其中的第k级;
匹配准则是最小化匹配代价函数H(π):
式中π表示匹配点对的排列;
匹配流程如下:
(4.1)现有pi和qj两个边缘点集,对于pi中的点i,分别寻找qj中Cost值最小的点j;
(4.2)将匹配的信息保存;
(4.3)重复(4.1),对剩余的点进行匹配,直到匹配完所有点。
(4.4)判断H(π)是否为最小,否则重复步骤(4.1)~(4.3);
5)滤除误匹配
步骤4)的匹配过程存在误匹配点对,消除误匹配过程如下:
(5.1)计算匹配点对(pi,qj)之间的距离;
D(i,j)=||pi|-|qj||,(i,j)∈π (3)
其中,π表示匹配点对的排列;
(5.2)统计D(i,j)分布情况,求其概率密度函数,式(4),(5)分别表示其距离分布的均值mean(D(i,j))和方差sigma(D(i,j)):
式中L是匹配点对个数,μ代表式(4)的均值mean(D(ij));
根据式(4),(5)求其概率密度函数。
(5.3)滤除误匹配,依据小概率事件理论把匹配点间的距离大于Γ的情况视为误匹配点对,予以滤除;
Γ=mean(D(i,j))+2.58*sigma(D(i,j)) (7)
6)边缘校正与插值
为了实现参考图像与浮动图像血管组织的像素级配准,根据匹配点对血管壁做迭代演化,使形变后的参考图像与浮动图像血管轮廓能量函数E(pi,qj)达到最小。
式中pi,qj代表参考图像血管壁边缘和浮动图像血管壁边缘,π表示匹配点对的排列;
根据边缘校正后的结果,采用样条插值的方法对浮动图像的局部血管区域做插值,最终得到配准后的结果。
进一步,所述步骤2)中,采用小波变换去噪方法,在小波分解的每一层都选自适应的选取阈值,在小波分解的低频域采用中值滤波器滤除分布在低频信号中的噪声。
再进一步,所述步骤2)中,阈值的函数如式(9):
其中αn是噪声的标准差,σg,j是无噪声图像g在小波域内第j层的标准差,M即为小波域内小波系数的总体个数,可调参数k1,k2满足关系k1+k2=1,αj为小波域第j层的小波系数取为αj=1/2j-1,j为小波分解层数;
σn取值是小波域内的最高频小波系数,其值由噪声图像小波分解后最高频子带的中值计算得到:
为噪声图像小波分解的最高频子带;
小波系数在j层的标准差由噪声小波系数的标准差计算得到:
本发明的有益效果主要表现在:基于shape context描述子的多序列医学血管图像的精确配准,以在图像结构复杂,质量(对比度,信噪比)不高的情况下达到像素级别的精度,以满足计算机辅助诊断对血管图像配准精度的要求。
附图说明
图1是基于多序列医学图像的血管配准方法的流程图。
图2是shape context描述子的示意图,其中,(a)表示我们要描述的目标图形,(b)是提取其边缘轮廓的结果,并对边缘轮廓实施降采样,(c)把边缘轮廓映射到以各边缘点为中心的极坐标内,(d)图反映了边缘上每一点与其他点之间的关系,(e)图为边缘点的shapecontext描述子统计直方图。
图3是血管边缘匹配与矫正过程示意图,其中,(a)表示参考图像与浮动图像血管之间的差异,(b)是用shape context描述子做匹配后的结果(匹配点对中直线连接),(c)图是经过迭代后参考图像与浮动图像血管边缘之间的差异,(d)是最后边缘校正后结果。
图4是配准结果示意图,其中,(a)、(d)、(g)分别表示三种浮动图像序列T1WI、T2WI、T1GD的血管边缘,(b)、(e)、(h)表示浮动图像的血管边缘与参考图像血管之间的差异,(c)、(f)、(i)是用本文配准方法配准的最终结果。
具体实施方式
下面结合附图对本发明作进一步描述。
参照图1~图4,一种基于多序列医学图像的血管配准方法,包括以下步骤:
1)获取多序列医学血管图像(如多序列磁共振图像(MRI)、多序列CT图像等)。本说明以多序列磁共振图像为例说明其技术流程,多序列采用的是1.5T磁场强度下采集的T1(T1加权成像)、T1GD(T1造影剂成像)、T2(T2加权成像)、PD(质子密度成像)、STIR(短翻转时间翻转恢复成像)序列图像。在配准过程中把T1图像作为配准的参考图像,其他序列图像作为浮动图像;
2)由于医学图像会受到各种干扰而产生噪声影响图像质量,如MR图像会受到外界电磁设备等干扰获取的图像信号收到随机干扰,进而影响图像质量,为此在做血管配准之前有必要做相应的去噪处理。本文在传统的小波变换去噪的基础上改进了阈值函数,在小波分解的每一层都选自适应的选取阈值,在小波分解的低频域采用中值滤波器滤除分布在低频信号中的噪声;
(1)阈值函数的选取
小波去噪过程中阈值函数的好坏直接影响去噪效果,当阈值选择较小的时候,一部分大于该阈值的噪声系数就会被保留下来,去噪后的图像中仍然会存在大量的噪声,如果阈值函数选择的过大,就会使很多有用的图像细节和边缘信息被过滤掉,导致图像退化。Chang等人提出了一种最优阈值选择法,这种阈值函数式(12)是由贝叶斯最大后验概率推导得到。
其中是噪声的方差,σg,j是无噪声图像g在小波域内第j层的标准差。αj是小波域内第j层的可调节系数,通常取为1。本发明在结合经典小波阈值函数的基础上提出了一种更加适用于医学图像噪声模型的阈值函数式(10)。
其中αn是噪声的标准差,M即为小波域内小波系数的总体个数,σg,j是无噪声图像g在小波域内第j层的标准差,可调参数k1,k2满足关系k1+k2=1。αj为小波域第j层的小波系数取为αj=1/2j-1,j为小波分解层数。
σn取值是小波域内的最高频小波系数,其值由噪声图像小波分解后最高频子带的中值计算得到:
为噪声图像小波分解的最高频子带。
小波系数在j层的标准差由噪声小波系数的标准差计算得到:
根据小波变换的特点,图像的噪声更多的分布在小波变换的高频子带,因此高频子带的阈值应该相较于低频子带大。通过在式(2)中引入项,可以实现对每一层的小波系数选择合适的值,另一方面也弱化了统一阈值函数对小波系数总体个数M的影响,将使阈值函数的选取更具有针对性。
(2)结合中值滤波
经过上述软阈值去噪处理后并不能完全抑制噪声对图像的影响,为将原始图像真实完整的还原出来本文采用对图像细节信息和边缘信息都能很好保护的中值滤波器对阈值后的图像做进一步处理。
3)血管边缘的shape context描述,血管内外壁分别记为Clumen(X),Cwall(X),X代表边缘点坐标。为了实现浮动图像动脉血管边缘点和参考图像边缘点做匹配,边缘上每一点分别用Belongie等人提出的Shape Context[1]描述子统计其直方图。
附图(2)简要说明shape context描述子。如图(2-b)有n个点,然后进行极坐标变换。某点Pi与其余的(n-1)个点均存在关系如图(2-d),即产生(n-1)个向量,这(n-1)个向量描述了丰富轮廓信息,决定了目标的形状特征,然后分别统计与其关联的向量得到第i个点的直方图如图(3-e)所示。如果n越大,信息量越大,描述就越准确。极坐标直方图hi(k)内包含了点与点的全局关系,即shape context所需的描述。
我们采用同样的方法描述血管的边缘点。
4)边缘点匹配
边缘点匹配是指T1序列动脉血管边缘点与T1GD,PD,STIR序列血管边缘点做匹配,匹配准则是根据参考图像血管边缘点与浮动图像血管边缘点的代价函数Cij
其中pi和qj分别代表颈动脉血管参考图像序列(T1WI)和浮动图像序列(T1GD、PD、STIR)的边缘点集,hi(k)、hj(k)分别代表参考图像边缘点直方图和浮动图像边缘点直方图,K代表直方图总共的等级数,k是其中的第k级;
匹配准则是最小化匹配代价函数H(π):
式中π表示匹配点对的排列。
匹配流程如下:
(1)现有pi和qj两个边缘点集,对于pi中的点i,分别寻找qj中Cost值最小的点j;
(2)将匹配的信息保存;
(3)重复第(1)步,对剩余的点进行匹配,直到匹配完所有点。
(4)判断H(π)是否为最小,否则重复步骤1~3。
5)滤除误匹配
上述的匹配过程存在一定的误匹配点对,误匹配点对于浮动图像颈动脉血管壁的迭代演化变形有巨大影响,为此本文提出了以下消除误匹配方法:
1)计算匹配点对(pi,qj)之间的距离;
D(i,j)=||pi|-|qj||,(i,j)∈π (3)
式中π表示匹配点对的排列。
2)统计D(i,j)分布情况,求其概率密度函数,式(4),(5)分别表示其距离分布的均值和方差。
式中L是匹配点对个数,μ代表式(8)的均值mean(D(ij))。
根据式(4),(5)求其概率密度函数。
3)滤除误匹配,本文依据小概率事件理论把匹配点间的距离大于Γ的情况视为误匹配点对,予以滤除。
Γ=mean(D(i,j))+2.58*sigma(D(i,j)) (7)
边缘校正和迭代的结果如图(3)所示。
6)边缘校正与插值
为了实现参考图像与浮动图像血管组织的像素级配准,下一步根据匹配点对血管壁做迭代演化,使形变后的参考图像与浮动图像血管轮廓能量函数E(pi,qj)达到最小。
式中pi,qj代表参考图像血管壁边缘和浮动图像血管壁边缘,π表示匹配点对的排列;
根据边缘校正后的结果,采用样条插值的方法对浮动图像的局部血管区域做插值,最终得到配准后的结果。
插值结果如附图(4)所示。

Claims (3)

1.一种基于多序列医学图像的血管配准方法,其特征在于:所述配准方法包括以下步骤:
1)获取多序列医学血管图像,选定一个序列图像为配准的参考图像,其他序列图像作为浮动图像;
2)对所述多序列医学血管图像进行去噪处理;
3)血管边缘的shape context描述,血管内外壁分别记为Clumen(X),Cwall(X),X代表边缘点坐标,边缘上每一点分别用Shape Context描述子统计其直方图;
4)边缘点匹配
根据参考图像血管边缘点与浮动图像血管边缘点的代价函数Ci,j
其中,pi和qj分别代表颈动脉血管参考图像序列和浮动图像序列的边缘点集,hi(k)、hj(k)分别代表参考图像边缘点直方图和浮动图像边缘点直方图,K代表直方图总共的等级数,k是其中的第k级;
匹配准则是最小化匹配代价函数H(π):
式中π表示匹配点对的排列;
匹配流程如下:
(4.1)现有pi和qj两个边缘点集,对于pi中的点i,分别寻找qj中Cost值最小的点j;
(4.2)将匹配的信息保存;
(4.3)重复(4.1),对剩余的点进行匹配,直到所有点全部匹配完成;
(4.4)判断H(π)是否为最小,否则重复步骤(4.1)~(4.3);
5)滤除误匹配
步骤4)的匹配过程存在误匹配点对,消除误匹配过程如下:
(5.1)计算匹配点对(pi,qj)之间的距离D(i,j)
D(i,j)=||pi|-|qj||,(i,j)∈π (3)
其中,π表示匹配点对的排列;
(5.2)统计D(i,j)分布情况,求其概率密度函数,式(4),(5)分别表示其距离分布的均值mean(D(i,j))和方差sigma(D(i,j)):
式中L是匹配点对个数,μ代表式(4)的均值mean(D(ij));
根据式(4),(5)求其概率密度函数P(D(i,j)):
(5.3)滤除误匹配,依据小概率事件理论把匹配点间的距离大于Γ的情况视为误匹配点对,予以滤除;
Γ=mean(D(i,j))+2.58*sigma(D(i,j)) (7)
6)边缘校正与插值
为了实现参考图像与浮动图像血管组织的像素级配准,根据匹配点对血管壁做迭代演化,使形变后的参考图像与浮动图像血管轮廓能量函数E(pi,qj)达到最小;
式中,pi和qj分别代表颈动脉血管参考图像序列和浮动图像序列的边缘点集,π表示匹配点对的排列;
根据边缘校正后的结果,采用样条插值的方法对浮动图像的局部血管区域做插值,最终得到配准后的结果。
2.如权利要求1所述的一种基于多序列医学图像的血管配准方法,其特征在于:所述步骤2)中,采用小波变换去噪方法,在小波分解的每一层都自适应的选取阈值,在小波分解的低频域采用中值滤波器滤除分布在低频信号中的噪声。
3.如权利要求2所述的一种基于多序列医学图像的血管配准方法,其特征在于:所述步骤2)中,阈值T的函数式(9):
其中αn是噪声的标准差,σg,j是无噪声图像g在小波域内第j层的标准差,M即为小波域内小波系数的总体个数,可调参数k1,k2满足关系 k1+k2=1,αj为小波域第j层的小波系数取为αj=1/2j-1,j为小波分解层数;
σn取值是小波域内的最高频小波系数,其值由噪声图像小波分解后最高频子带的中值计算得到:
为噪声图像小波分解的最高频子带;
小波系数在j层的标准差由噪声小波系数的标准差计算得到:
CN201510088135.2A 2015-02-26 2015-02-26 一种基于多序列医学图像的血管配准方法 Active CN104766304B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201510088135.2A CN104766304B (zh) 2015-02-26 2015-02-26 一种基于多序列医学图像的血管配准方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201510088135.2A CN104766304B (zh) 2015-02-26 2015-02-26 一种基于多序列医学图像的血管配准方法

Publications (2)

Publication Number Publication Date
CN104766304A CN104766304A (zh) 2015-07-08
CN104766304B true CN104766304B (zh) 2017-12-05

Family

ID=53648113

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201510088135.2A Active CN104766304B (zh) 2015-02-26 2015-02-26 一种基于多序列医学图像的血管配准方法

Country Status (1)

Country Link
CN (1) CN104766304B (zh)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105205812A (zh) * 2015-09-01 2015-12-30 哈尔滨工业大学 一种基于微小卫星星座的多帧图像重构方法
CN107067420A (zh) * 2017-04-28 2017-08-18 上海联影医疗科技有限公司 图像处理方法、装置及设备
CN107133959B (zh) * 2017-06-12 2020-04-28 上海交通大学 一种快速的血管边界三维分割方法及系统
CN109697702A (zh) * 2018-11-02 2019-04-30 浙江工业大学 基于弯曲波变换的医学超声图像去噪方法
CN114782358A (zh) * 2022-04-18 2022-07-22 上海博动医疗科技股份有限公司 一种血管形变自动计算的方法、装置及存储介质

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102446366A (zh) * 2011-09-14 2012-05-09 天津大学 时空联合多视角视频插值及三维建模方法
CN102679871A (zh) * 2012-05-07 2012-09-19 上海交通大学 亚像素精度工业物体快速检测方法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102446366A (zh) * 2011-09-14 2012-05-09 天津大学 时空联合多视角视频插值及三维建模方法
CN102679871A (zh) * 2012-05-07 2012-09-19 上海交通大学 亚像素精度工业物体快速检测方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
Adaptive wavelet thresholding for image denoising and compression;Chang S G et al;《IEEE Transactions on Image Processing》;20020806;第9卷(第9期);第1532-1546页 *
Shape matching and object recognition using shape contexts;Belongie S et al;《IEEE Transactions on Pattern Analysis and Machine Intelligence》;20020807;第24卷(第4期);第509-522页 *
基于薄板样条和形状内容的医学图像非刚性配准方法研究;吴月娥 等;《航天医学与医学工程》;20070228;第20卷(第1期);第43-46页 *

Also Published As

Publication number Publication date
CN104766304A (zh) 2015-07-08

Similar Documents

Publication Publication Date Title
CN104766304B (zh) 一种基于多序列医学图像的血管配准方法
US11308587B2 (en) Learning method of generative adversarial network with multiple generators for image denoising
CN110097512B (zh) 基于Wasserstein生成对抗网络的MRI图像去噪模型的构建方法及应用
CN109191399B (zh) 基于改进的多路径匹配追踪算法的磁共振图像去噪方法
CN107194912B (zh) 基于稀疏表示的改进耦合字典学习的脑部ct/mr图像融合方法
CN115496771A (zh) 一种基于脑部三维mri图像设计的脑肿瘤分割方法
WO2022121100A1 (zh) 一种基于darts网络的多模态医学图像融合方法
CN110599530B (zh) 基于双正则约束的mvct图像纹理增强方法
CN115018728A (zh) 基于多尺度变换和卷积稀疏表示的图像融合方法及系统
CN110163825A (zh) 一种人类胚胎心脏超声图像去噪和增强方法
CN116645283A (zh) 基于自监督感知损失多尺度卷积神经网络的低剂量ct图像去噪方法
Xu et al. A novel multi-modal fundus image fusion method for guiding the laser surgery of central serous chorioretinopathy
CN107845081B (zh) 一种磁共振图像去噪方法
Mentl et al. Noise reduction in low-dose ct using a 3D multiscale sparse denoising autoencoder
Qiao et al. Automatic liver segmentation method based on improved region growing algorithm
CN102968779A (zh) 基于FCM聚类的非下采样Contourlet域MRI图像增强方法
Li et al. Skin lesion segmentation via dense connected deconvolutional network
CN116894783A (zh) 基于时变约束的对抗生成网络模型的金属伪影去除方法
CN106570837B (zh) 基于高斯多尺度空间的脑组织mri图像偏移场校正方法
CN111477304A (zh) 一种pet和mri图像相融合的肿瘤照射成像组合方法
CN111127344A (zh) 一种基于bp神经网络的自适应双边滤波的超声图像降噪的方法
CN110570369B (zh) 一种甲状腺结节超声图像去噪方法
CN112560808B (zh) 一种基于灰度信息的体内静脉识别方法及装置
CN109377461A (zh) 一种基于nsct的乳腺x射线图像自适应增强方法
CN114387173A (zh) 一种oct图像降噪方法、电子设备及存储介质

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
EXSB Decision made by sipo to initiate substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant