CN113298855B - 基于自动勾画的图像配准方法 - Google Patents

基于自动勾画的图像配准方法 Download PDF

Info

Publication number
CN113298855B
CN113298855B CN202110585065.7A CN202110585065A CN113298855B CN 113298855 B CN113298855 B CN 113298855B CN 202110585065 A CN202110585065 A CN 202110585065A CN 113298855 B CN113298855 B CN 113298855B
Authority
CN
China
Prior art keywords
oar
point
image
matrix
points
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
CN202110585065.7A
Other languages
English (en)
Other versions
CN113298855A (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.)
Perception Vision Medical Technology Co ltd
Original Assignee
Perception Vision Medical Technology Co ltd
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 Perception Vision Medical Technology Co ltd filed Critical Perception Vision Medical Technology Co ltd
Priority to CN202110585065.7A priority Critical patent/CN113298855B/zh
Publication of CN113298855A publication Critical patent/CN113298855A/zh
Priority to PCT/CN2021/136311 priority patent/WO2022247218A1/zh
Application granted granted Critical
Publication of CN113298855B publication Critical patent/CN113298855B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/30—Determination of transform parameters for the alignment of images, i.e. image registration
    • G06T7/33—Determination of transform parameters for the alignment of images, i.e. image registration using feature-based methods
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N3/00—Computing arrangements based on biological models
    • G06N3/02—Neural networks
    • G06N3/04—Architecture, e.g. interconnection topology
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N3/00—Computing arrangements based on biological models
    • G06N3/02—Neural networks
    • G06N3/08—Learning methods
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/10—Segmentation; Edge detection
    • G06T7/12—Edge-based segmentation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V10/00—Arrangements for image or video recognition or understanding
    • G06V10/70—Arrangements for image or video recognition or understanding using pattern recognition or machine learning
    • G06V10/74—Image or video pattern matching; Proximity measures in feature spaces
    • G06V10/75—Organisation of the matching processes, e.g. simultaneous or sequential comparisons of image or video features; Coarse-fine approaches, e.g. multi-scale approaches; using context analysis; Selection of dictionaries
    • G06V10/752—Contour matching
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10072—Tomographic images
    • G06T2207/10084—Hybrid tomography; Concurrent acquisition with multiple different tomographic modalities
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20076—Probabilistic image processing
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20081—Training; Learning
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/20—Special algorithmic details
    • G06T2207/20084—Artificial neural networks [ANN]

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Evolutionary Computation (AREA)
  • Artificial Intelligence (AREA)
  • General Health & Medical Sciences (AREA)
  • Health & Medical Sciences (AREA)
  • Software Systems (AREA)
  • Computing Systems (AREA)
  • Biophysics (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Physics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Biomedical Technology (AREA)
  • General Engineering & Computer Science (AREA)
  • Computational Linguistics (AREA)
  • Molecular Biology (AREA)
  • Medical Informatics (AREA)
  • Databases & Information Systems (AREA)
  • Multimedia (AREA)
  • Image Analysis (AREA)
  • Apparatus For Radiation Diagnosis (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)
  • Image Processing (AREA)

Abstract

本发明公开了一种基于自动勾画的图像配准方法，包括：输入两个任意模态的医学图像；采用域自适应的迁移学习策略训练神经网络，在训练数据仅有CT的OAR勾画金标准时，使其能对除CT外的其它模态也具备识别OAR的能力；对两个输入图像分别输入上述神经网络，得到各自的全身OAR分割结果；基于OAR轮廓点采用迭代式优化方法得到刚性配准的变换矩阵和位移向量，从而得到刚性配准后的医学图像warped image；在刚性配准的基础上，针对每个OAR，对fixed image中该OAR的轮廓点与对应moving image中该OAR的轮廓点进行匹配，进而得到轮廓点的位移值；通过基于控制点的插值法得到完整的位移场三维矩阵，从而得到非刚性配准后的医学图像warped image。借此，能有效解决OAR金标准缺乏的问题。

Description

基于自动勾画的图像配准方法
技术领域
本发明涉及图像处理领域、深度学习领域、医疗领域，尤其是一种基于自动勾画的图像配准方法。
背景技术
图像配准在医学图像处理与分析中有众多具有实用价值的应用。随着医学成像设备的进步，对于同一患者，可以采集含有准确解剖信息的多种不同模态的图像，如CT、CBCT、MRI、PET等。然而，通过观察不同图像进行诊断需要凭着空间想象和医生的主观经验。采用正确的图像配准方法则可以将多种多样的信息准确地融合到同一图像中，使医生更方便更精确地从各个角度观察病灶和结构。同时，通过对不同时刻采集的动态图像的配准，可以定量分析病灶和器官的变化情况，使得医疗诊断、制定手术计划、放射治疗计划更准确可靠。
传统的图像配准方法基于相似度目标函数的优化求解问题，容易收敛至局部极小值，尤其对不同模态图像的配准效果较差，且迭代求解的过程耗时较长。而基于危及器官(OAR)勾画的图像配准方法能解决上述问题，但OAR金标准的获得需要耗费医生、专家的大量时间，成本较高。近年来，人们对探索利用人工智进行诊断产生了浓厚的兴趣，并在某些领域利用AI算法建立了表现优于人类医学专家的数学模型。因此有理由相信，利用AI算法对传统图像配准方法进行改进能有效提高图像配准的效果。
公开于该背景技术部分的信息仅仅旨在增加对本发明的总体背景的理解，而不应当被视为承认或以任何形式暗示该信息构成已为本领域一般技术人员所公知的现有技术。
发明内容
本发明的目的在于提供一种基于自动勾画的图像配准方法，其能够利用AI算法对传统图像配准方法进行改进能有效提高图像配准的效果。
为实现上述目的，本发明提供了一种基于自动勾画的图像配准方法，包括以下步骤：输入两个任意模态(CT、CBCT、MRI、PET等)的医学图像，一个作为fixed image(参考图像)，另一个作为moving image(待配准图像)；采用域自适应的迁移学习策略训练神经网络，在训练数据仅有CT的OAR(危及器官)勾画金标准时，使其能对除CT外的其它模态(CBCT、MRI、PET等)也具备识别OAR的能力；对两个输入图像分别输入上述神经网络，得到各自的全身OAR分割结果；基于OAR轮廓点采用迭代式优化方法得到刚性配准的变换矩阵和位移向量，从而得到刚性配准后的医学图像warped image；在刚性配准的基础上，针对每个OAR，对fixed image中该OAR的轮廓点与对应moving image中该OAR的轮廓点进行匹配，进而得到轮廓点的位移值；以及通过基于控制点的插值法得到完整的位移场三维矩阵，从而得到非刚性配准后的医学图像warped image。
在一优选的实施方式中，基于自动勾画的图像配准方法还包括：采用域自适应的迁移学习策略训练神经网络，在训练数据仅有CT的全身危及器官(OAR)的勾画标记时，使其能对除CT外的其它模态(CBCT、MRI、PET等)也具备识别OAR的能力，包括：神经网络采用GAN架构，包含生成器网络与判别器网络。生成器用于输出OAR分割结果与混淆判别器，判别器用于判断生成器输出的特征与结果是否属于CT。生成器网络的结构包含若干个域适配网络(每种模态一个)及公共的主网络。以Unet为基础，域适配网络由2个Residual Modules(残差模块)和下采样交替组成，主网络的前半部分由1个Residual Modules和下采样交替组成，与域适配网络组成U型结构的左侧encoder(编码器)，而主网络的后半部分由3个Residual Modules和上采样交替组成U型结构的右侧decoder(解码器)，并通过skip-connection将主网络的编码器中部分较浅层的高分辨率特征与解码器中对应层的特征进行融合，以补充由下采样带来的细节损失，主网络最后一个Residual Modules后接一个channel(通道)数为OAR个数+1(加上背景)的卷积层，输出多channel概率图，表示每个像素属于某个OAR或背景的概率，进而得到所有OAR的分割结果。判别器网络的结构由4个卷积模块和下采样交替组成，后接一个全局池化层与全连接层，输出结果表示输入特征属于CT的概率。训练时，选择CT与另一种模态的图像一起进行训练。两个图像分别送入各自的域适配网络及公共主网络，提取主网络中所有Residual Modules输出的feature map(特征图)并与最终输出的概率图合并为一个总的feature map，输入判别器判断该特征是否属于CT。即生成器与判别器交替训练，当训练生成器时任务为降低判别器的准确率，当训练判别器时任务为提高判断准确率。同时，CT经过生成器后输出的概率图与专家勾画的OAR金标准进行有监督学习。测试及应用时，输入图像并根据其模态选择对应的域适配网络，与公共主网络组成生成器，输出OAR分割结果。神经网络能同时识别全身的OAR。
在一优选的实施方式中，基于自动勾画的图像配准方法还包括：基于OAR轮廓点采用迭代式优化方法得到刚性配准的变换矩阵和位移向量，从而得到刚性配准后的医学图像warped image，包括：对每个OAR，分割结果为二值图的三维矩阵，把矩阵与对该矩阵进行腐蚀操作后做差得到轮廓二值图，从而得到该OAR所有轮廓点的物理坐标。取不易发生形变的一个OAR(如，骨性结构)的轮廓点用于下述的刚性配准步骤：对fixed image中该OAR的每个轮廓点，找到在moving image中该OAR的最近轮廓点，形成匹配点对。得到所有匹配点对后，通过最小化下式，得到变换矩阵和位移向量的最优解：
Figure GDA0003343327740000041
最优解为：
R＝(PTP)-1PTQ
A＝R[0:3,0:3]
b＝R[0:3,3]
其中，N为匹配点对个数，pn为fixed image的第n个匹配点，qn为对应的movingimage中的像素点。P为fixed image所有匹配点组成的矩阵，大小为[N,4]，即N个四维行向量组成的矩阵，四维的前三维是像素点的物理坐标，第四维是固定值1。Q为moving image所有匹配点组成的矩阵，大小为[N,4]。矩阵R的大小为[4,4]，R[0:3,0:3]指取矩阵R的前3行、前3列组成的大小为[3,3]的矩阵，R[0:3,3]指取矩阵R的前3行第3列的三维列向量。A和b分别为变换矩阵和位移向量的最优解。
用A和b作用于fixed image中的每个轮廓点，得到变换后的轮廓点用于下一轮迭代。
重复上述步骤，即找到在moving image中该OAR的最近轮廓点，形成匹配点对后求解得到第t轮迭代的解At和bt。当符合下式时，迭代结束：
‖At-At-1‖+‖bt-bt-1‖<10-6
最后通过A和b得到刚性配准后的医学图像warped image。
在一优选的实施方式中，基于自动勾画的图像配准方法还包括：在刚性配准的基础上，针对每个OAR，对fixed image中该OAR的轮廓点与对应moving image中该OAR的轮廓点进行匹配，进而得到轮廓点的位移值，包括：针对每个OAR，用刚性配准后的变换矩阵和位移向量作用于fixed image中该OAR的每个轮廓点，得到刚性配准后新的轮廓点位置，并在此基础上与对应moving image中该OAR的轮廓点进行匹配，步骤如下：用Kuhn-Munkres算法得到所有匹配点对。循环fixed image上的每一个匹配点pn，执行以下优化匹配点的步骤：若pn的相邻两点pn-1与pn+1在moving image中的对应匹配点qn-1与qn+1不相邻，则qn-1与qn+1间所有轮廓点，连同pn，从上述识别OAR的神经网络中提取主网络所有Residual Modules输出的feature map(特征图)合并为一个总的feature map，从而得到相应点的特征向量。计算pn的特征向量与qn-1与qn+1间所有轮廓点的特征向量的相似度，选取相似度最高的点作为新的匹配点。相似度的计算公式为：
Figure GDA0003343327740000051
其中，F(pn)为pn的特征向量，F(qk)为qk的特征向量。
最后得到最终的fixed image上每个OAR的轮廓点的moving image中的匹配点。
在一优选的实施方式中，基于自动勾画的图像配准方法还包括：通过基于控制点的插值法得到完整的位移场三维矩阵，从而得到非刚性配准后的医学图像warped image，包括：初始化位移场的三维矩阵为零矩阵，其大小与fixed image的矩阵相同，对fixedimage中轮廓点的位移赋值为与该点对应的在moving image中的匹配点的物理坐标与该点物理坐标的差。然后采用如下基于控制点的样条插值法得到位移场矩阵剩余像素点的值。
Figure GDA0003343327740000061
B0(t)＝(1-t)3/6
B1(t)＝(3t3-6t2+4)/6
B2(t)＝(-3t3+3t2+3t+1)/6
B3(t)＝t3/6
假设由一组m×n×l个控制点组成的网格作为新的像素坐标系，其中(x,y,z)为位移场矩阵中某个已赋值的像素点在此新坐标系下的坐标位置，
Figure GDA0003343327740000063
表示x的往下取整，
Figure GDA0003343327740000064
为坐标位置在
Figure GDA0003343327740000065
Figure GDA0003343327740000066
的控制点，f(x,y,z)为该像素点的位移值，f′(x,y,z)为通过该像素点附近的16个控制点拟合的近似值。
Figure GDA0003343327740000062
对于某个像素点，满足f(x，y，z)＝f′(x，y，z)的控制点的值有多种解，因此加上最小化上式的约束条件后求解，得到控制点的最优解为：
Figure GDA0003343327740000071
Figure GDA0003343327740000074
其中，φi，j，k为
Figure GDA0003343327740000075
的简写。
当能求解控制点φi，j，k的像素点不只一个时，各像素点通过上式求解得到该控制点的最优解可能不同，通过最小化下式，得到φi，j，k的最终解：
min e(φi，j，k)＝∑s(wsφi，j，k-wsφs)2
最终解为：
Figure GDA0003343327740000072
其中：
Figure GDA0003343327740000073
Figure GDA0003343327740000076
Si，j，k＝{(xs，ys，zs)|i-2≤xs＜i+2，j-2≤ys＜j+2，k-2≤zs＜k+2}
Si，j，k为参与求解控制点φi，j，k的所有已赋值的像素点的集合。φs为通过坐标为(xs，ys，zs)的像素点得到的控制点最优解。
当能求解控制点φi，j，k的像素点个数为0时，对其赋值为0。
在得到所有控制点的值后，位移场中剩余像素点的位移值均通过其附近的16个控制点拟合得到。由于位移值是三维向量，即x、y、z方向，故上述插值过程需重复3次，即每个方向分别进行一次。最后得到位移场三维矩阵，从而得到非刚性配准后的医学图像warpedimage。
与现有技术相比，本发明的基于自动勾画的图像配准方法具有以下有益效果：本发明能对任意两个模态的图像进行配准，且采用域自适应的迁移学习策略训练神经网络，在训练数据仅有CT的OAR勾画金标准时，使其能对除CT外的其它模态也具备识别OAR的能力，能有效解决OAR金标准缺乏的问题。
附图说明
图1是根据本发明一实施方式的基于自动勾画的图像配准方法的流程示意图；
图2是根据本发明一实施方式的基于自动勾画的图像配准方法的模型结构示意图。
具体实施方式
下面结合附图，对本发明的具体实施方式进行详细描述，但应当理解本发明的保护范围并不受具体实施方式的限制。
除非另有其它明确表示，否则在整个说明书和权利要求书中，术语“包括”或其变换如“包含”或“包括有”等等将被理解为包括所陈述的元件或组成部分，而并未排除其它元件或其它组成部分。
如图1所示，根据本发明优选实施方式的一种基于自动勾画的图像配准方法，主要包括以下步骤：
输入两个任意模态(CT、CBCT、MRI、PET等)的医学图像，一个作为fixed image，另一个作为moving image。
采用域自适应的迁移学习策略训练神经网络，在训练数据仅有CT的OAR(危及器官)勾画金标准时，使其能对除CT外的其它模态(CBCT、MRI、PET等)也具备识别OAR的能力。
对两个输入图像分别输入上述神经网络，得到各自的全身OAR勾画结果。
基于OAR轮廓点采用迭代式优化方法得到刚性配准的变换矩阵和位移向量，从而得到刚性配准后的医学图像warped image。
在刚性配准的基础上，针对每个OAR，对fixed image中该OAR的轮廓点与对应moving image中该OAR的轮廓点进行匹配，进而得到轮廓点的位移值。
通过基于控制点的插值法得到完整的位移场三维矩阵，从而得到非刚性配准后的医学图像warped image。
在一些实施方式中，基于自动勾画的图像配准方法还包括以下步骤：
S1、参照图1，采用域自适应的迁移学习策略训练神经网络，在训练数据仅有CT的全身危及器官(OAR)的勾画标记时，使其能对除CT外的其它模态(CBCT、MRI、PET等)也具备识别OAR的能力；
步骤S1具体包括：
S11、神经网络采用GAN架构，包含生成器网络与判别器网络。生成器用于输出OAR分割结果与混淆判别器，判别器用于判断生成器输出的特征与结果是否属于CT。
S12、生成器网络的结构包含若干个域适配网络(每种模态一个)及公共的主网络。以Unet为基础，域适配网络由2个Residual Modules(残差模块)和下采样交替组成，主网络的前半部分由1个Residual Modules和下采样交替组成，与域适配网络组成U型结构的左侧encoder(编码器)，而主网络的后半部分由3个Residual Modules和上采样交替组成U型结构的右侧decoder(解码器)，并通过skip-connection将主网络的编码器中部分较浅层的高分辨率特征与解码器中对应层的特征进行融合，以补充由下采样带来的细节损失，主网络最后一个Residual Modules后接一个channel(通道)数为OAR个数+1(加上背景)的卷积层，输出多channel概率图，表示每个像素属于某个OAR或背景的概率，进而得到所有OAR的分割结果。
S13、判别器网络的结构由4个卷积模块和下采样交替组成，后接一个全局池化层与全连接层，输出结果表示输入特征属于CT的概率。
S14、训练时，选择CT与另一种模态的图像一起进行训练。两个图像分别送入各自的域适配网络及公共主网络，提取主网络中所有Residual Modules输出的feature map(特征图)并与最终输出的概率图合并为一个总的feature map，输入判别器判断该特征是否属于CT。即生成器与判别器交替训练，当训练生成器时任务为降低判别器的准确率，当训练判别器时任务为提高判断准确率。同时，CT经过生成器后输出的概率图与专家勾画的OAR金标准进行有监督学习。
S15、测试及应用时，输入图像并根据其模态选择对应的域适配网络，与公共主网络组成生成器，输出OAR分割结果。
S16、神经网络能同时识别全身的OAR。
在一些实施方式中，基于自动勾画的图像配准方法还包括下列步骤：
S2、对两个输入图像分别输入上述神经网络，得到各自的全身OAR分割结果。
在一些实施方式中，基于自动勾画的图像配准方法还包括下列步骤：
S3、基于OAR轮廓点采用迭代式优化方法得到刚性配准的变换矩阵和位移向量，从而得到刚性配准后的医学图像warped image；
步骤S3具体包括以下步骤：
S31、对每个OAR，分割结果为二值图的三维矩阵，把矩阵与对该矩阵进行腐蚀操作后做差得到轮廓二值图，从而得到该OAR所有轮廓点的物理坐标。
S32、取不易发生形变的骨性结构OAR的轮廓点用于下述的刚性配准步骤：
S33、对fixed image中该OAR的每个轮廓点，找到在moving image中该OAR的最近轮廓点，形成匹配点对。得到所有匹配点对后，通过最小化下式，得到变换矩阵和位移向量的最优解：
Figure GDA0003343327740000111
最优解为：
R＝(PTP)-1PTQ
A＝R[0：3，0：3]
b＝R[0：3，3]
其中，N为匹配点对个数，pn为fixed image的第n个匹配点，qn为对应的movingimage中的像素点。P为fixed image所有匹配点组成的矩阵，大小为[N，4]，即N个四维行向量组成的矩阵，四维的前三维是像素点的物理坐标，第四维是固定值1。Q为moving image所有匹配点组成的矩阵，大小为[N，4]。矩阵R的大小为[4，4]，R[0：3，0：3]指取矩阵R的前3行、前3列组成的大小为[3，3]的矩阵，R[0：3，3]指取矩阵R的前3行第3列的三维列向量。A和b分别为变换矩阵和位移向量的最优解。
S34、用A和b作用于fixed image中的每个轮廓点，得到变换后的轮廓点用于下一轮迭代。
S35、重复步骤S23～S24得到第t轮迭代的解At和bt。
S36、当符合下式时，迭代结束：
||At-At-1||+||bt-bt-1||＜10-6
S37、最后通过A和b得到刚性配准后的医学图像warped image。
在一些实施方式中，基于自动勾画的图像配准方法还包括下列步骤：
S4、在刚性配准的基础上，针对每个OAR，对fixed image中该OAR的轮廓点与对应moving image中该OAR的轮廓点进行匹配，进而得到轮廓点的位移值；
步骤S4具体包括以下步骤：
S41、针对每个OAR，用刚性配准后的变换矩阵和位移向量作用于fixedimage中该OAR的每个轮廓点，得到刚性配准后新的轮廓点位置，并在此基础上与对应moving image中该OAR的轮廓点进行匹配，步骤如下：
S42、用Kuhn-Munkres算法得到所有匹配点对。
S43、循环fixed image上的每一个匹配点pn，执行以下优化匹配点的步骤：
(1)若pn的相邻两点pn-1与pn+1在moving image中的对应匹配点qn-1与qn+1不相邻，则qn-1与qn+1间所有轮廓点，连同pn，从上述识别OAR的神经网络中提取主网络所有ResidualModules输出的feature map(特征图)合并为一个总的feature map，从而得到相应点的特征向量。
(2)计算pn的特征向量与qn-1与qn+1间所有轮廓点的特征向量的相似度，选取相似度最高的点作为新的匹配点。相似度的计算公式为：
Figure GDA0003343327740000131
其中，F(pn)为pn的特征向量，F(qk)为qk的特征向量。
S44、对每个OAR执行步骤S41～S43，最后得到最终的fixed image上每个OAR的每个轮廓点的在moving image中的匹配点。
在一些实施方式中，基于自动勾画的图像配准方法还包括下列步骤：
S5、通过基于控制点的插值法得到完整的位移场三维矩阵，从而得到非刚性配准后的医学图像warped image；
步骤S5具体包括以下步骤：
S51、初始化位移场的三维矩阵为零矩阵，其大小与fixed image的矩阵相同，对fixed image中轮廓点的位移赋值为与该点对应的在moving image中的匹配点的物理坐标与该点物理坐标的差。然后采用如下基于控制点的样条插值法得到位移场矩阵剩余像素点的值。
(1)求解控制点：
Figure GDA0003343327740000141
B0(t)＝(1-t)3/6
B1(t)＝(3t3-6t2+4)/6
B2(t)＝(-3t3+3t2+3t+1)/6
B3(t)＝t3/6
假设由一组m×n×l个控制点组成的网格作为新的像素坐标系，其中(x，y，z)为位移场矩阵中某个已赋值的像素点在此新坐标系下的坐标位置，
Figure GDA0003343327740000144
表示x的往下取整，
Figure GDA0003343327740000145
为坐标位置在
Figure GDA0003343327740000146
Figure GDA0003343327740000147
的控制点，f(x，y，z)为该像素点的位移值，f′(x，y，z)为通过该像素点附近的16个控制点拟合的近似值。
Figure GDA0003343327740000142
对于某个像素点，满足f(x，y，z)＝f′(x，y，z)的控制点的值有多种解，因此加上最小化上式的约束条件后求解，得到控制点的最优解为：
Figure GDA0003343327740000143
Figure GDA0003343327740000148
其中，φi，j，k为
Figure GDA0003343327740000149
的简写。
当能求解控制点φi，j，k的像素点不只一个时，各像素点通过上式求解得到该控制点的最优解可能不同，通过最小化下式，得到φi，j，k的最终解：
min e(φi，j，k)＝∑s(wsφi，j，k-wsφs)2
最终解为：
Figure GDA0003343327740000151
其中：
Figure GDA0003343327740000152
Figure GDA0003343327740000153
Si，j，k＝{(xs，ys，zs)|i-2≤xs＜i+2，j-2≤ys＜j+2，k-2≤zs＜k+2}
Si，j，k为参与求解控制点φi，j，k的所有已赋值的像素点的集合。φs为通过坐标为(xs，ys，zs)的像素点得到的控制点最优解。
当能求解控制点φi，j，k的像素点个数为0时，对其赋值为0。
(2)在得到所有控制点的值后，位移场中剩余像素点的位移值均通过其附近的16个控制点拟合得到。
S52、由于位移值是三维向量，即x、y、z方向，故上述插值过程需重复3次，即每个方向分别进行一次。
S53、最后得到位移场三维矩阵，从而得到非刚性配准后的医学图像warpedimage。
综上所述，本发明的基于自动勾画的图像配准方法具有以下优点：本发明能对任意两个模态的图像进行配准，且采用域自适应的迁移学习策略训练神经网络，在训练数据仅有CT的OAR勾画金标准时，使其能对除CT外的其它模态也具备识别OAR的能力，能有效解决OAR金标准缺乏的问题。
前述对本发明的具体示例性实施方案的描述是为了说明和例证的目的。这些描述并非想将本发明限定为所公开的精确形式，并且很显然，根据上述教导，可以进行很多改变和变化。对示例性实施例进行选择和描述的目的在于解释本发明的特定原理及其实际应用，从而使得本领域的技术人员能够实现并利用本发明的各种不同的示例性实施方案以及各种不同的选择和改变。本发明的范围意在由权利要求书及其等同形式所限定。

Claims (5)

1.一种基于自动勾画的图像配准方法，其特征在于，包括以下步骤：
输入两个任意模态的医学图像，一个作为fixed image，另一个作为moving image，其中所述两个任意模态的医学图像为CT、CBCT、MRI及PET中的任意两个；
采用域自适应的迁移学习策略训练神经网络，在训练数据仅有CT的OAR勾画金标准时，使其能对除CT外的CBCT、MRI及PET的所述医学图像也具备识别OAR的能力；
对两个输入的所述医学图像分别输入所述迁移学习策略训练神经网络，得到各自的全身OAR分割结果；
基于OAR轮廓点采用迭代式优化方法得到刚性配准的变换矩阵和位移向量，从而得到刚性配准后的医学图像warped image；
在刚性配准的基础上，针对每个OAR，对fixed image中该OAR的轮廓点与对应movingimage中该OAR的轮廓点进行匹配，进而得到轮廓点的位移值；以及
通过基于控制点的插值法得到完整的位移场三维矩阵，从而得到非刚性配准后的医学图像warped image。
2.根据权利要求1所述的一种基于自动勾画的图像配准方法，其特征在于，还包括：
采用域自适应的迁移学习策略训练神经网络，在训练数据仅有CT的OAR的勾画标记时，使其能对除CT外的其它模态的所述医学图像也具备识别OAR的能力，包括：
神经网络采用GAN架构，包含生成器网络与判别器网络；其中所述生成器网络用于输出OAR分割结果与混淆判别器，所述判别器网络用于判断生成器输出的特征与结果是否属于CT；
所述生成器网络的结构包含若干个域适配网络及公共的主网络；以Unet为基础，所述域适配网络由2个Residual Modules和下采样交替组成，所述主网络的前半部分由1个Residual Modules和下采样交替组成，与所述域适配网络组成U型结构的左侧encoder，而所述主网络的后半部分由3个Residual Modules和上采样交替组成U型结构的右侧decoder，并通过skip-connection将所述主网络的编码器中部分较浅层的高分辨率特征与解码器中对应层的特征进行融合，以补充由下采样带来的细节损失，主网络最后一个Residual Modules后接一个channel数为OAR个数+1的卷积层，输出多channel概率图，表示每个像素属于某个OAR或背景的概率，进而得到所有OAR的分割结果；
所述判别器网络的结构由4个卷积模块和下采样交替组成，后接一个全局池化层与全连接层，输出结果表示输入特征属于CT的概率；
训练时，选择CT与另一种模态的图像一起进行训练；两个图像分别送入各自的域适配网络及公共主网络，提取主网络中所有Residual Modules输出的feature map并与最终输出的概率图合并为一个总的feature map，输入判别器判断该特征是否属于CT；即生成器与判别器交替训练，当训练生成器时任务为降低判别器的准确率，当训练判别器时任务为提高判断准确率；同时，CT经过生成器后输出的概率图与专家勾画的OAR金标准进行有监督学习；以及
测试及应用时，输入图像并根据其模态选择对应的所述域适配网络，与公共的所述主网络组成生成器，输出OAR分割结果；以及
神经网络能同时识别全身的OAR。
3.根据权利要求1所述的一种基于自动勾画的图像配准方法，其特征在于，还包括：
基于OAR轮廓点采用迭代式优化方法得到刚性配准的变换矩阵和位移向量，从而得到刚性配准后的医学图像warped image，包括：
对每个OAR，分割结果为二值图的三维矩阵，把矩阵与对该矩阵进行腐蚀操作后做差得到轮廓二值图，从而得到该OAR所有轮廓点的物理坐标；
取不易发生形变的一个OAR的轮廓点用于下述的刚性配准步骤：
对fixed image中该OAR的每个轮廓点，找到在moving image中该OAR的最近轮廓点，形成匹配点对，得到所有匹配点对后，通过最小化下式，得到变换矩阵和位移向量的最优解：
Figure FDA0003343327730000031
最优解为：
R＝(PTP)-1PTQ
A＝R[0:3,0:3]
b＝R[0:3,3]
其中，N为匹配点对个数，pn为fixed image的第n个匹配点，qn为对应的moving image中的像素点；P为fixed image所有匹配点组成的矩阵，大小为[N,4]，即N个四维行向量组成的矩阵，四维的前三维是像素点的物理坐标，第四维是固定值1；Q为moving image所有匹配点组成的矩阵，大小为[N,4]；矩阵R的大小为[4,4]，R[0:3,0:3]指取矩阵R的前3行、前3列组成的大小为[3,3]的矩阵，R[0:3,3]指取矩阵R的前3行第3列的三维列向量；A和b分别为变换矩阵和位移向量的最优解；
用A和b作用于fixed image中的每个轮廓点，得到变换后的轮廓点用于下一轮迭代；
重复上述步骤，即找到在moving image中该OAR的最近轮廓点，形成匹配点对后求解得到第t轮迭代的解At和bt；以及
当符合下式时，迭代结束：
||At-At-1||+||bt-bt-1||＜10-6
最后通过A和b得到刚性配准后的医学图像warped image。
4.根据权利要求1所述的一种基于自动勾画的图像配准方法，其特征在于，还包括：
在刚性配准的基础上，针对每个OAR，对fixed image中该OAR的轮廓点与对应movingimage中该OAR的轮廓点进行匹配，进而得到轮廓点的位移值，包括：
针对每个OAR，用刚性配准后的变换矩阵和位移向量作用于fixed image中该OAR的每个轮廓点，得到刚性配准后新的轮廓点位置，并在此基础上与对应moving image中该OAR的轮廓点进行匹配，步骤如下：
用Kuhn-Munkres算法得到所有匹配点对；
循环fixed image上的每一个匹配点pn，执行以下优化匹配点的步骤：
若pn的相邻两点pn-1与pn+1在moving image中的对应匹配点qn-1与qn+1不相邻，则qn-1与qn+1间所有轮廓点，连同pn，从上述识别OAR的神经网络中提取主网络所有ResidualModules输出的feature map(特征图)合并为一个总的feature map，从而得到相应点的特征向量；
计算pn的特征向量与qn-1与qn+1间所有轮廓点的特征向量的相似度，选取相似度最高的点作为新的匹配点；相似度的计算公式为：
Figure FDA0003343327730000041
其中，F(pn)为pn的特征向量，F(qk)为qk的特征向量；以及
最后得到最终的fixed image上每个OAR的轮廓点的moving image中的匹配点。
5.根据权利要求1所述的一种基于自动勾画的图像配准方法，其特征在于，还包括：
通过基于控制点的插值法得到完整的位移场三维矩阵，从而得到非刚性配准后的医学图像warped image，包括：
初始化位移场的三维矩阵为零矩阵，其大小与fixed image的矩阵相同，对fixedimage中轮廓点的位移赋值为与该点对应的在moving image中的匹配点的物理坐标与该点物理坐标的差；然后采用如下基于控制点的样条插值法得到位移场矩阵剩余像素点的值；
Figure FDA0003343327730000051
B0(t)＝(1-t)3/6
B1(t)＝(3t3-6t2+4)/6
B2(t)＝(-3t3+3t2+3t+1)/6
B3(t)＝t3/6
假设由一组m×n×l个控制点组成的网格作为新的像素坐标系，其中(x，y，z)为位移场矩阵中某个已赋值的像素点在此新坐标系下的坐标位置，
Figure FDA0003343327730000053
表示x的往下取整，
Figure FDA0003343327730000054
为坐标位置在
Figure FDA0003343327730000055
Figure FDA0003343327730000056
的控制点，f(x，y，z)为该像素点的位移值，f′(x，y，z)为通过该像素点附近的16个控制点拟合的近似值；
Figure FDA0003343327730000052
对于某个像素点，满足f(x，y，z)＝f′(x，y，z)的控制点的值有多种解，因此加上最小化上式的约束条件后求解，得到控制点的最优解为：
Figure FDA0003343327730000061
Figure FDA0003343327730000065
其中，φi，j，k为
Figure FDA0003343327730000064
的简写；
当能求解控制点φi，j，k的像素点不只一个时，各像素点通过上式求解得到该控制点的最优解可能不同，通过最小化下式，得到φi，j，k的最终解：
min e(φi，j，k)＝∑s(wsφi，j，k-wsφs)2
最终解为：
Figure FDA0003343327730000062
其中：
Figure FDA0003343327730000063
Figure FDA0003343327730000066
Si，j，k＝{(xs，ys，zs)|i-2≤xs＜i+2，j-2≤ys＜j+2，k-2≤zs＜k+2}
Si，j，k为参与求解控制点φi，j，k的所有已赋值的像素点的集合；φs为通过坐标为(xs，ys，zs)的像素点得到的控制点最优解；
当能求解控制点φi，j，k的像素点个数为0时，对其赋值为0；
在得到所有控制点的值后，位移场中剩余像素点的位移值均通过其附近的16个控制点拟合得到；
由于位移值是三维向量，即x、y、z方向，故上述插值过程需重复3次，即每个方向分别进行一次；以及
最后得到位移场三维矩阵，从而得到非刚性配准后的医学图像warped image。
CN202110585065.7A 2021-05-27 2021-05-27 基于自动勾画的图像配准方法 Active CN113298855B (zh)

Priority Applications (2)

Application Number Priority Date Filing Date Title
CN202110585065.7A CN113298855B (zh) 2021-05-27 2021-05-27 基于自动勾画的图像配准方法
PCT/CN2021/136311 WO2022247218A1 (zh) 2021-05-27 2021-12-08 基于自动勾画的图像配准方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110585065.7A CN113298855B (zh) 2021-05-27 2021-05-27 基于自动勾画的图像配准方法

Publications (2)

Publication Number Publication Date
CN113298855A CN113298855A (zh) 2021-08-24
CN113298855B true CN113298855B (zh) 2021-12-28

Family

ID=77325620

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110585065.7A Active CN113298855B (zh) 2021-05-27 2021-05-27 基于自动勾画的图像配准方法

Country Status (2)

Country Link
CN (1) CN113298855B (zh)
WO (1) WO2022247218A1 (zh)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113298855B (zh) * 2021-05-27 2021-12-28 广州柏视医疗科技有限公司 基于自动勾画的图像配准方法
CN113920178B (zh) * 2021-11-09 2022-04-12 广州柏视医疗科技有限公司 一种基于标记点的多视觉2d-3d图像配准方法及系统
CN114187338B (zh) * 2021-12-08 2023-04-28 卡本(深圳)医疗器械有限公司 一种基于估算2d位移场的器官形变配准方法
CN114796901B (zh) * 2022-05-30 2024-08-27 北京大学第一医院 一种腰骶神经根的自动勾画方法、设备及存储介质

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111862022A (zh) * 2020-07-13 2020-10-30 中山大学 全身多部位放射治疗危及器官自动勾画方法

Family Cites Families (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20170337682A1 (en) * 2016-05-18 2017-11-23 Siemens Healthcare Gmbh Method and System for Image Registration Using an Intelligent Artificial Agent
CN107403201A (zh) * 2017-08-11 2017-11-28 强深智能医疗科技（昆山）有限公司 肿瘤放射治疗靶区和危及器官智能化、自动化勾画方法
CN109934861B (zh) * 2019-01-22 2022-10-18 广东工业大学 一种头颈部多模态医学图像自动配准方法
US20220180964A1 (en) * 2019-03-28 2022-06-09 Phase Genomics, Inc. Systems and methods for karyotyping by sequencing
CN110210486B (zh) * 2019-05-15 2021-01-01 西安电子科技大学 一种基于素描标注信息的生成对抗迁移学习方法
US11651487B2 (en) * 2019-07-12 2023-05-16 The Regents Of The University Of California Fully automated four-chamber segmentation of echocardiograms
CN111127444B (zh) * 2019-12-26 2021-06-04 广州柏视医疗科技有限公司 基于深度语义网络的ct影像中自动识别放疗危及器官的方法
CN111199550B (zh) * 2020-04-09 2020-08-11 腾讯科技（深圳）有限公司 图像分割网络的训练方法、分割方法、装置和存储介质
CN111784706B (zh) * 2020-06-28 2021-06-04 广州柏视医疗科技有限公司 鼻咽癌原发肿瘤图像自动识别方法及系统
CN111862174B (zh) * 2020-07-08 2023-10-03 清华大学深圳国际研究生院 一种跨模态医学图像配准方法及装置
CN112419320B (zh) * 2021-01-22 2021-04-27 湖南师范大学 基于sam与多层uda的跨模态心脏分割方法
CN112733859B (zh) * 2021-01-25 2023-12-19 重庆大学 一种组织病理学图像的深度迁移半监督域自适应分类方法
CN113298855B (zh) * 2021-05-27 2021-12-28 广州柏视医疗科技有限公司 基于自动勾画的图像配准方法

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111862022A (zh) * 2020-07-13 2020-10-30 中山大学 全身多部位放射治疗危及器官自动勾画方法

Also Published As

Publication number Publication date
WO2022247218A1 (zh) 2022-12-01
CN113298855A (zh) 2021-08-24

Similar Documents

Publication Publication Date Title
CN113298855B (zh) 基于自动勾画的图像配准方法
CN113298854B (zh) 基于标记点的图像配准方法
US11514573B2 (en) Estimating object thickness with neural networks
CN108416802B (zh) 一种基于深度学习的多模医学图像非刚性配准方法及系统
CN111932550B (zh) 一种基于深度学习的3d心室核磁共振视频分割系统
CN110363802B (zh) 基于自动分割和骨盆对齐的前列腺图像配准系统及方法
CN107403201A (zh) 肿瘤放射治疗靶区和危及器官智能化、自动化勾画方法
CN113674330B (zh) 一种基于生成对抗网络的伪ct影像生成系统
CN115457020B (zh) 一种融合残差图像信息的2d医学图像配准方法
CN112598649B (zh) 基于生成对抗网络的2d/3d脊椎ct非刚性配准方法
Jin et al. Object recognition in medical images via anatomy-guided deep learning
CN116258732A (zh) 一种基于pet/ct图像跨模态特征融合的食管癌肿瘤靶区分割方法
CN110853048A (zh) 基于粗、精训练的mri图像分割方法、装置和存储介质
He et al. Cephalometric landmark detection by considering translational invariance in the two-stage framework
CN113870327A (zh) 基于预测多层次变形场的医学图像配准方法
CN114581453A (zh) 基于多轴面特征融合二维卷积神经网络的医学图像分割方法
Dai et al. CAN3D: Fast 3D medical image segmentation via compact context aggregation
CN115830016A (zh) 医学图像配准模型训练方法及设备
CN118279361A (zh) 基于无监督深度学习和模态转换的多模态医学图像配准方法
CN117853547A (zh) 一种多模态的医学图像配准方法
CN117333751A (zh) 一种医学图像融合方法
Gan et al. Probabilistic modeling for image registration using radial basis functions: Application to cardiac motion estimation
CN113850710B (zh) 一种跨模态医学影像精准转换方法
CN108596900B (zh) 甲状腺相关性眼病医学影像数据处理装置、方法、计算机可读存储介质及终端设备
CN116309754A (zh) 一种基于局部-全局信息协作的大脑医学图像配准方法及系统

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant