CN115511932A - 基于医学影像图像的配准方法、可读存储介质及电子设备 - Google Patents
基于医学影像图像的配准方法、可读存储介质及电子设备 Download PDFInfo
- Publication number
- CN115511932A CN115511932A CN202211198336.4A CN202211198336A CN115511932A CN 115511932 A CN115511932 A CN 115511932A CN 202211198336 A CN202211198336 A CN 202211198336A CN 115511932 A CN115511932 A CN 115511932A
- Authority
- CN
- China
- Prior art keywords
- polygonal mesh
- medical image
- vertex
- mesh
- registration
- 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.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 98
- 238000003860 storage Methods 0.000 title claims abstract description 23
- 238000013528 artificial neural network Methods 0.000 claims abstract description 50
- 230000009466 transformation Effects 0.000 claims abstract description 33
- 238000013507 mapping Methods 0.000 claims abstract description 8
- 210000004556 brain Anatomy 0.000 claims description 31
- 210000000056 organ Anatomy 0.000 claims description 29
- 238000011176 pooling Methods 0.000 claims description 28
- 238000012952 Resampling Methods 0.000 claims description 24
- 230000006870 function Effects 0.000 claims description 21
- 230000008569 process Effects 0.000 claims description 18
- 238000010606 normalization Methods 0.000 claims description 12
- 230000004913 activation Effects 0.000 claims description 11
- 210000005013 brain tissue Anatomy 0.000 claims description 8
- 238000005070 sampling Methods 0.000 claims description 8
- 238000010586 diagram Methods 0.000 claims description 7
- 210000001638 cerebellum Anatomy 0.000 claims description 5
- 210000000133 brain stem Anatomy 0.000 claims description 3
- 230000009191 jumping Effects 0.000 claims description 2
- 238000012545 processing Methods 0.000 abstract description 14
- 230000001054 cortical effect Effects 0.000 description 63
- 210000003710 cerebral cortex Anatomy 0.000 description 35
- 238000012549 training Methods 0.000 description 21
- 238000004364 calculation method Methods 0.000 description 12
- 238000003384 imaging method Methods 0.000 description 11
- 239000011159 matrix material Substances 0.000 description 11
- 238000013461 design Methods 0.000 description 9
- 238000005192 partition Methods 0.000 description 9
- 230000008901 benefit Effects 0.000 description 7
- 238000012360 testing method Methods 0.000 description 6
- 238000009826 distribution Methods 0.000 description 5
- 230000000694 effects Effects 0.000 description 5
- 238000002595 magnetic resonance imaging Methods 0.000 description 5
- 238000012916 structural analysis Methods 0.000 description 4
- 239000013598 vector Substances 0.000 description 4
- 230000005540 biological transmission Effects 0.000 description 3
- 238000002591 computed tomography Methods 0.000 description 3
- 238000012937 correction Methods 0.000 description 3
- 238000002599 functional magnetic resonance imaging Methods 0.000 description 3
- 238000009499 grossing Methods 0.000 description 3
- 230000033001 locomotion Effects 0.000 description 3
- 238000003062 neural network model Methods 0.000 description 3
- PTLRDCMBXHILCL-UHFFFAOYSA-M sodium arsenite Chemical compound [Na+].[O-][As]=O PTLRDCMBXHILCL-UHFFFAOYSA-M 0.000 description 3
- 210000001519 tissue Anatomy 0.000 description 3
- 238000012795 verification Methods 0.000 description 3
- 238000012935 Averaging Methods 0.000 description 2
- 241001465754 Metazoa Species 0.000 description 2
- 230000002490 cerebral effect Effects 0.000 description 2
- 210000004720 cerebrum Anatomy 0.000 description 2
- 238000006243 chemical reaction Methods 0.000 description 2
- 238000004891 communication Methods 0.000 description 2
- 238000013527 convolutional neural network Methods 0.000 description 2
- 210000004884 grey matter Anatomy 0.000 description 2
- 230000003287 optical effect Effects 0.000 description 2
- 230000011218 segmentation Effects 0.000 description 2
- 238000004904 shortening Methods 0.000 description 2
- 238000002604 ultrasonography Methods 0.000 description 2
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 description 1
- RYGMFSIKBFXOCR-UHFFFAOYSA-N Copper Chemical compound [Cu] RYGMFSIKBFXOCR-UHFFFAOYSA-N 0.000 description 1
- LFQSCWFLJHTTHZ-UHFFFAOYSA-N Ethanol Chemical compound CCO LFQSCWFLJHTTHZ-UHFFFAOYSA-N 0.000 description 1
- 241000251539 Vertebrata <Metazoa> Species 0.000 description 1
- 238000002583 angiography Methods 0.000 description 1
- 210000001557 animal structure Anatomy 0.000 description 1
- 238000003491 array Methods 0.000 description 1
- 210000000988 bone and bone Anatomy 0.000 description 1
- 210000004027 cell Anatomy 0.000 description 1
- 210000002421 cell wall Anatomy 0.000 description 1
- 238000012512 characterization method Methods 0.000 description 1
- 238000013135 deep learning Methods 0.000 description 1
- 238000002598 diffusion tensor imaging Methods 0.000 description 1
- 238000002597 diffusion-weighted imaging Methods 0.000 description 1
- 210000002249 digestive system Anatomy 0.000 description 1
- 230000010339 dilation Effects 0.000 description 1
- 238000006073 displacement reaction Methods 0.000 description 1
- 210000000750 endocrine system Anatomy 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 238000004880 explosion Methods 0.000 description 1
- 239000000835 fiber Substances 0.000 description 1
- 230000005484 gravity Effects 0.000 description 1
- 210000002216 heart Anatomy 0.000 description 1
- 229910052739 hydrogen Inorganic materials 0.000 description 1
- 239000001257 hydrogen Substances 0.000 description 1
- 125000004435 hydrogen atom Chemical class [H]* 0.000 description 1
- 238000003709 image segmentation Methods 0.000 description 1
- 210000003734 kidney Anatomy 0.000 description 1
- 210000004072 lung Anatomy 0.000 description 1
- 210000002751 lymph Anatomy 0.000 description 1
- 230000000877 morphologic effect Effects 0.000 description 1
- 210000003205 muscle Anatomy 0.000 description 1
- 210000000653 nervous system Anatomy 0.000 description 1
- 239000005416 organic matter Substances 0.000 description 1
- 230000010412 perfusion Effects 0.000 description 1
- 230000035790 physiological processes and functions Effects 0.000 description 1
- 238000002600 positron emission tomography Methods 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 210000004994 reproductive system Anatomy 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 210000002345 respiratory system Anatomy 0.000 description 1
- 230000002441 reversible effect Effects 0.000 description 1
- 210000001202 rhombencephalon Anatomy 0.000 description 1
- 238000002603 single-photon emission computed tomography Methods 0.000 description 1
- 210000003625 skull Anatomy 0.000 description 1
- 238000005563 spheronization Methods 0.000 description 1
- 210000002784 stomach Anatomy 0.000 description 1
- 238000011426 transformation method Methods 0.000 description 1
- 230000001131 transforming effect Effects 0.000 description 1
- 230000002485 urinary effect Effects 0.000 description 1
- 210000004885 white matter Anatomy 0.000 description 1
- 238000002759 z-score normalization Methods 0.000 description 1
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
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
- G06T17/20—Finite element generation, e.g. wire-frame surface description, tesselation
-
- 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/10081—Computed x-ray tomography [CT]
-
- 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/10088—Magnetic resonance imaging [MRI]
-
- 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]
-
- 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/30—Subject of image; Context of image processing
- G06T2207/30004—Biomedical image processing
- G06T2207/30016—Brain
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Graphics (AREA)
- Geometry (AREA)
- Software Systems (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Apparatus For Radiation Diagnosis (AREA)
Abstract
本发明公开了一种基于医学影像图像的配准方法、可读存储介质和电子设备，属于医学图像处理领域。该配准方法包括以下步骤：根据该医学影像图像进行三维模型重建获得多边形网格和该多边形网格中的每个顶点的属性特征；将该多边形网格膨胀并映射至球面以获得具有该属性特征的球化后的多边形网格；将该球化后的多边形网格和所对应的球化模版输入刚性配准图神经网络中，以获得该球化后的多边形网格配准至该所对应的球化模版的刚性配准欧拉角；将该球化后的多边形网格基于该刚性配准欧拉角旋转变换以获得刚性配准后的多边形网格。
Description
技术领域
本发明涉及医学图像处理领域，特别涉及一种基于医学影像图像的配准方法、可读存储介质及电子设备。
背景技术
在对大脑皮层进行结构分析时，通常需要经过大脑皮层表面重建、大脑皮层表面配准和大脑皮层表面的结构分区等多个环节。大脑皮层表面配准是指同一区域内以不同成像手段所获得的不同图像图形的地理坐标的匹配，配准通常包括几何纠正、投影变换与统一比例尺三方面的处理。大脑皮层表面结构分区是根据大脑的结构或连接信息，将大脑皮层表面划分为多个不同的具有神经生物学意义的区域。
大脑皮层的结构分析主要依据大脑皮层表面的结构分区来进行，因此大脑皮层表面的结构分区是大脑皮层结构分析中最主要的环节。而由于大脑皮层表面的结构分区是通过非线性配准将模版上标记的结构分区投射到个体大脑皮层上完成的，因此大脑皮层的配准的准确性将会影响到结构分区的准确性。
同时，由于在结构分析过程中，大脑皮层表面重建和大脑皮层表面配准由于技术较为复杂，因此在处理过程中往往将耗费十几个小时才能完成图像重建和配准。
因此，确有必要提供一种准确性高且能够快速实现配准的新的基于医学影像图像的配准方法、可读存储介质及电子设备。
发明内容
为了解决现有技术中存在的上述问题和缺陷的至少一个方面，本发明提供了一种基于医学影像图像的配准方法、可读存储介质及电子设备，能够至少部分地实现快速进行三维模型配准、三维模型配准的准确性高的技术效果。所述技术方法如下：
根据本发明的一个方面，提供了一种基于医学影像图像的配准方法，包括以下步骤：
根据所述医学影像图像进行三维模型重建获得多边形网格和所述多边形网格中的每个顶点的属性特征；
将所述多边形网格膨胀并映射至球面以获得具有所述属性特征的球化后的多边形网格；
将所述球化后的多边形网格和所对应的球化模版输入刚性配准图神经网络中，以获得所述球化后的多边形网格配准至所述所对应的球化模版的刚性配准欧拉角；
将所述球化后的多边形网格基于所述刚性配准欧拉角旋转变换以获得刚性配准后的多边形网格。
根据本发明的另一方面，提供了一种可读存储介质，其中，所述可读存储介质上存储程序或指令，所述程序或指令被处理器执行时以执行上述的医学影像图像配准方法。
根据本发明的还一方面，提供了一种电子设备，所述电子设备包括处理器和存储器，其中，
所述存储器上存储程序，所述程序被所述处理器执行时以执行上述任一项所述的基于医学影像图像的配准方法。
根据本发明实施例的基于医学影像图像的配准方法、可读存储介质及电子设备具有以下优点中的至少一个：
(1)本发明提供的基于医学影像图像的配准方法、可读存储介质及电子设备能够通过刚性配准图神经网络的设计，能够快速实现基于医学影像图像的三维模型配准，由此显著缩短了基于医学影像的三维模型的配准时间，大大提升了基于医学影像图像的三维模型处理效率；
(2)本发明提供的基于医学影像图像的配准方法、可读存储介质及电子设备通过刚性配准图神经网络和非刚性配准图神经网络的设计，在配准过程中依次通过粗配准对每个顶点位置的统一调整和精配准对每个顶点的位置分别进行的精细调整，从而能够获得更加精准的大脑皮层表面配准三维模型。
附图说明
本发明的这些和/或其他方面和优点从下面结合附图对优选实施例的描述中将变得明显和容易理解，其中：
图1是根据本发明一个实施例的基于医学影像图像的配准方法的流程图；
图2A-图2C是图1所示的多边形网格球化过程的示意图，其中图2A为医学影像图像经三维模型重建后得到的多边形网格，图2B为多边形网格膨胀后得到的多边形网格，图2C为膨胀后的三维模型映射至球面得到的球化后的多边形网格；
图3是图1所示的刚性配准图神经网络的结构示意图；
图4是图1所示的球化后的多边形网格基于刚性配准欧拉角旋转变换的流程图；
图5是图1所示的欧拉角的原理图；
图6是图1所示的刚性配准后的多边形网格基于非刚性配准欧拉角旋转变换的流程图；
图7是图1所示的非刚性配准图神经网络的结构示意图；
图8A是图1所示的刚性配准后的多边形网格的局部示意图；
图8B是图8A所示的刚性配准后的多边形网格发生折叠时的效果示意图；
图8C是图8A所示的刚性配准后的多边形网格未发生折叠时的效果示意图。
具体实施方式
下面通过实施例，并结合附图，对本发明的技术方案作进一步具体的说明。在说明书中，相同或相似的附图标号指示相同或相似的部件。下述参照附图对本发明实施方式的说明旨在对本发明的总体发明构思进行解释，而不应当理解为对本发明的一种限制。
本文中所使用的术语“医学影像图像”应当被广义地理解为为了医疗或医学研究，对人体或人体某部分，以非侵入方式取得的内部组织影像图像。
本文中所使用的术语“多边形网格表示”，即mesh表示，应当被广义地理解为在三维空间中，通过一系列顶点和一组轮廓线连接构成的曲面并彼此组合以表示三维曲面的一种表示方式。通常包括三角形网格、四边形网格或者其他凸多边形网格。
本文中所使用的术语“模版”，应当被广义地理解为是指通过将正常人大脑皮层表面进行平均得到的具有代表性的大脑结构或者特性分布。例如通过将大脑的结构、沟回深度和曲率进行平均得到的模版。
本文中所使用的术语“欧拉角”，应当被广义地理解为用来确定定点转动刚体位置的3个一组独立角参量(α,β,γ)。欧拉角(α,β,γ)的三个角度分别表示依次围绕Z、Y和X轴旋转的角度表示方式。欧拉角表示与旋转矩阵具有等价的表示能力，且两者之间能够互相转换。
本文中所使用的术语“三维模型”，应当被广义地理解为物体的多边形表示，通常用计算机或者其它视频设备进行显示。
参见图1，示出了根据本发明的一个实施例的基于医学影像图像的配准方法。该配准方法包括以下步骤：
步骤S1根据医学影像图像进行三维模型重建获得多边形网格和多边形网格中的每个顶点的属性特征；
步骤S2将多边形网格膨胀并映射至球面以获得具有属性特征的球化后的多边形网格；
步骤S3将球化后的多边形网格和所对应的球化模版输入刚性配准图神经网络中，以获得球化后的多边形网格配准至所对应的球化模版的刚性配准欧拉角；
步骤S4将球化后的多边形网格基于刚性配准欧拉角旋转变换以获得刚性配准后的多边形网格。
在一个示例中，采集的医学影像图像包括解剖成像图像和功能成像图像。解剖成像图像包括X射线摄影(Digital Radiography，DR)成像图像、电子计算机断层扫描(Computed Tomography，CT)成像图像、超声成像(Ultrasound，US)图像、磁共振成像(Magnetic Resonance Imaging，MRI)图像、内窥镜成像图像以及数字减影血管造影(Digital Subtraction Angiography，DSA)图像。功能成像图像包括正电子发射型计算机断层显像(Positron Emission Computed Tomography，PET)图像、单光子发射计算机断层成像(Single-Photon Emission Computed Tomography，SPECT)图像、功能性磁共振成像(functional magnetic resonance imaging，FMRI)图像。
优选地，医学影像图像包括CT图像、MRI图像和FMRI图像。更加优选地，医学影像图像为MRI图像。MRI图像可以为T1加权成像(T1 weighted imaging，T1WI)图像、T2加权成像(T2 weighted imaging，T2WI)图像、弥散加权成像(Diffusion Weighted Imaging，DWI)图像、弥散张量成像(Diffusion Tensor Imaging，DTI)图像和/或灌注成像图像。优选地，MRI图像为T1WI图像。
在一个示例中，医学影像图像可以为人体器官的医学影像图像，也可以为动物器官的医学影像图像。人体器官包括人体八大系统中的各个器官，人体八大系统包括运动系统、神经系统、内分泌系统、循环系统、呼吸系统、消化系统、泌尿系统以及生殖系统。器官可以是几种不同类型的组织经发育分化并相互结合构成具有一定形态和功能的结构；也可以是几种组织相互结合，组成具有一定形态和功能的结构；还可以是生物体内能担任某种独立的生理机能的部分。例如，人体器官可以是胃、肾、心脏、肺部、脑组织、骨骼、淋巴或者肌肉等。优选地，人体器官为脑组织，脑组织包括大脑、小脑和脑干中的至少一种。例如，脑组织可以为大脑，也可以为小脑，还可以为大脑和小脑的组合，当然还可以为大脑、小脑和脑干的组合。在一个示例中，动物应当被理解为不能将无机物合成有机物，由细胞构成的、且细胞具有细胞核但没有细胞壁的一类生命体。动物包括无脊椎动物、脊椎动物以及尾索动物与头索动物。
在一个示例中，多边形网格可以为三角形网格、四边形网格等。优选地，多边形网格为三角形网格。在使用除三角形网格以外的其他类型多边形网格时，可以直接使用该其他类型多边形网格，也可以将该其他类型多边形网格进行曲面细化(例如三角剖分)转换成三角形网格之后再使用。
在一个示例中，属性特征包括多边形网格中的每个顶点的形态学特征，例如厚度、曲率、密度、深度、面积和各个顶点的位置坐标中的一个或者它们的任意组合。
下面将以大脑作为人体器官、以T1WI作为医学影像图像、以三角形网格作为多边形网格对基于医学影像图像的配准方法进行示例性说明，大脑的其他类型的基于医学影像图像的配准方法以及其他人体器官的基于医学影像图像的配准方法与大脑皮层的T1WI图像配准方法的原理相同或者相似，在此不再一一赘述。
在一个示例中，大脑皮层的属性特征包括脑沟回深度、曲率、皮层厚度、皮层面积、灰质密度和/或各个顶点的位置坐标等。优选地，属性特征包括各个顶点的位置坐标与脑沟回深度和曲率的任意组合。也就是说，属性特征可以为各个顶点的第一位置信息(例如坐标信息)和各个顶点处的第一脑沟回深度特征，也可以为各个顶点的第一位置信息和各个顶点处的皮层的第一曲率特征，还可以为各个顶点的第一位置信息、各个顶点处的第一脑沟回深度特征和第一曲率特征。
在一个示例中，采集的大脑T1WI图像为个体脑图像。在一个示例中，在进行三维模型重建之前可以对大脑的T1WI图像进行图像预处理，例如可以进行头动校正、图像分割(例如颅骨分割、灰质与白质分割等)、灰度不均匀矫正等。
在步骤S1中，可以采用FreeSurfer软件、CAT12软件、FSL软件、SPM软件、Mimics软件、Simpleware软件、3D-DOCTOR软件、Amira软件或者现有的基于深度学习的三维模型重建方法(例如FastSurfer软件、FastCSR软件)等进行大脑皮质表面的三维模型重建。对大脑皮质表面进行三维模型重建时，需要分别对右脑和左脑的皮质表面进行三维模型重建，由此得到的多边形网格为右脑皮质表面重建多边形网格和左脑皮质表面重建多边形网格，且右脑、左脑皮质表面重建多边形网格均为三角形网格。在获得右脑皮质表面重建多边形网格的同时，得到了右脑皮质表面的各个顶点的第一脑沟回深度特征和/或第一曲率特征的计算结果。例如在三角形网格表示的右脑皮层表面多边形网格上的每个顶点处均显示出相应顶点的第一脑沟回深度特征和/或第一曲率特征。在获得左脑皮质表面重建多边形网格的同时，得到了左脑皮质表面的各个顶点的第一脑沟回深度特征和/或第一曲率特征的计算结果。
结合图2A-2C所示，将多边形网格(即三维模型)(如图2A所示)球化得到具有第一脑沟回深度特征和/或第一曲率特征的球化后的大脑皮层表面多边形网格(如图2C所示)。球化是将个体大脑皮层表面的三维模型(例如三角形网格)进行膨胀，之后将膨胀得到的大脑皮层表面多边形网格(如图2B所示)映射至具有固定半径的球面的过程，通过球化的多边形网格仍然保留了三维模型(例如三角形网格)的网格表示特征。换句话说，球化后的大脑皮层表面多边形网格依然保留了三维模型(例如三角形网格)表示的特性，即：在映射过程中，大脑皮层表面三角网格中的每个顶点所包含的位置信息、脑沟回深度、曲率等属性在球化后的大脑皮层表面多边形网格上同样被体现，其中构成大脑皮层表面的三维模型(例如三角形网格)的顶点(vertex)的数量不变，三角面片(即三角形图元)的连接关系也不变，同时每个顶点上对应的第一脑沟回深度特征、第一曲率特征也不会改变，球化的多边形网格上改变的是各个顶点的空间位置(例如，各顶点的坐标)。通过这样的设计，使得后续输入刚性配准图神经网络的多边形网格能够具有第一位置信息、第一脑沟回深度特征和/或第一曲率特征的属性。
在一个示例中，通过膨胀处理可以将脑沟展平，这样便于在大脑皮层表面上展示大脑结构、功能分区等属性。
结合图3和图4所示，刚性配准图神经网络包括7个图卷积块和1个输出层，每个图卷积块包括2个图卷积子块和1个图池化层，每个图卷积子块均包括依次连接的图卷积子层、批标准化层和激活函数(例如，LeakyRelu激活函数)。本领域技术人员可以明白，图卷积块的个数根据输入多边形网格上的顶点数来确定，图卷积子块数量可以根据实验测试调整，例如可以设计为1个、2个、3个或者更多个。
在一个示例中，由于刚性配准图神经网络(Graph Neural Network，GNN)的输入为个体和模版的图结构和图属性，其为单帧数据，因此刚性配准图神经网络可以基于图卷积网络(Graph Convolutional Network，GCN)为架构设计，也可以基于图自编码器(GraphAuto-encoder，GAE)为架构设计当然还可以基于图注意力网络(Graph AttentionNetwork，GAT)为架构设计。优选地，刚性配准图神经网络可以是基于图卷积网络(GCN)的架构来设计的。由于本发明的刚性配准图神经网络是基于图卷积网络来设计的网络结构，因此其具有节点数量固定且两层之间具有固定的节点计算对应关系，而无需通过随机采样来确定计算关系，使得在配准任务中仅需要改变节点的坐标而无需预测节点之间的连接关系等优点。
在一个示例中，当刚性配准图神经网络为图卷积网络时，网络所输入的多边形网格可以是具有40962个顶点(即图节点)的球化后的大脑皮层表面多边形网格。当然本领域技术人员可以根据需要选择顶点数(即图节点数)为163842个的球化后的大脑皮层表面多边形网格作为输入，当然还可以根据需要选择顶点数为10242个或2562个的球化后的大脑皮层表面多边形网格作为输入。优选地，选择顶点数为40962个的球化后的大脑皮层表面多边形网格作为输入。
在一个示例中，当刚性配准图神经网络为图卷积网络时，输入的多边形网格为具有第一脑沟回深度特征的多边形网格或第一曲率特征的多边形网格时，每个顶点的输入的通道可以设计为2个，即输入图卷积网络的多边形网格的特征数量为40962×2。在每个顶点的2个通道中，其中一个通道用于个体的球化后的大脑皮层表面多边形网格的第一脑沟回深度特征或者第一曲率特征，另外一个通道用于大脑球化模版的第二脑沟回深度特征或者第二曲率特征。
当输入的多边形网格为具有第一脑沟回深度特征的多边形网格和具有第一曲率特征的多边形网格时，每个顶点的输入通道可以设计为4个，即输入图卷积网络的多边形网格的特征数量为40962×4。在每个顶点的4个通道中，其中2个通道分别用于个体的球化后的大脑皮层表面多边形网格的第一脑沟回深度特征和第一曲率特征，另外2个通道用于大脑球化模版的第二脑沟回深度特征和第二曲率特征。
在一个示例中，图卷积子层对输入的多边形网格上的每个顶点和与该顶点连接其它顶点进行卷积运算，以学习输入多边形网格中的属性特征分布。例如当输入的多边形网格的特征数量为40962×4时，图卷积子层可以对当前顶点以及与该顶点相邻连接的6个顶点一起进行卷积运算。
在一个示例中，7个图卷积块中的经过两个图卷积子层、批标准化(BatchNorm)层、LeakyRelu激活函数单元处理，变换为40962*8，其中8代表每个顶点的通道数，该刚性配准图神经网络是为球面刚性任务专门设计的，具体表现为在该任务中每层的图节点数量(即顶点数量)由该任务来确定，例如，在图3中经每层网络输入的多边形网格的顶点数量(即图节点数量)40962，10242，2562...的排列遵从了不同特征数量的脑模版的三角形网格节点数量(即三角形网格中顶点数量)的变化规律，且相邻层之间的池化计算关系与脑模版的网格简化关系相一致，即此任务决定了这部分图卷积网络的结构及计算关系。同时由于刚性配准的欧拉角表达了三维空间中个体的球化后的大脑皮层表面多边形网格到大脑球化模版的变换关系，因此图神经网络的输出的多边形网格特征数量应当设置为1×3，也就是说，图神经网络的输出层所输出的节点数目为1，该节点具有3个维度的角度(即欧拉角)。之后将该1个节点的属性特征应用至重采样前的球化后的大脑皮层表面多边形网格以进行统一调整。
在一个示例中，大脑球化模版Fixed(F)为通过正二十面体离散化得到的球面，即是正二十面体的分形体所表示的球面。例如，FreeSurfer软件包中提供的fsaverage模版，该模版顶点数量为163842个，其提供了T1WI图像的球化后脑模版。当然本领域技术人员应当理解，大脑球化模版还可以是FreeSurfer软件包中提供的fsaverage6模版，其顶点数量为40962个；fsaverage5模版，其顶点数量为10242个。优选地，采用fsaverage6模版作为大脑球化模版。
在一个示例中，由于每个个体球化后的大脑皮层表面多边形网格会由于个体之间大脑皮层结构、多边形网格质量、重建方法等因素影响而生成顶点数量和三角面片连接关系不一致的大脑皮层表面三角形网格，这种顶点数量和连接关系不恒定的特性不便于后续处理，因此需要先使用球面采样方法将个体大脑皮层三角形网格采样成具有统一顶点数量及连接关系的表示。在一个示例中，将球化后的大脑皮层表面多边形网格Moving(M)采用FreeSurfer的fsaverage6模版作为基准对个体大脑皮层三角形网格及在其顶点上对应的第一脑沟回深度特征、第一曲率特征通过重心差值法进行采样，以生成与大脑球化模版相同顶点数量和连接关系的个体大脑皮层表面重采样多边形网格。在该大脑皮层表面重采样多边形网格中各个顶点的第一脑沟回深度特征、第一曲率特征与重采样之前相应位置的第一脑沟回深度特征、第一曲率特征一致。
在一个示例中，当球化后的大脑皮层表面多边形网格为除三角形网格以外的多边形网格时，需要将该大脑皮层表面多边形网格先行进行曲面细分为三角形网格。在一个示例中，曲面细分的方法可以为Loop细分法、蝶形细分法、√3细分法、Catmull-Clark细分法、分段平滑细分法以及位移细分法中的至少一种，或者基于上述细分方法改进的细分方法。在将除三角形网格以外的多边形网格进行曲面细分以后，再对其进行重采样。
当球化后的大脑皮层表面多边形网格为球化后的大脑皮层表面三角形网格或者经曲面细分为球化后的大脑皮层表面三角形网格时，对该球化后的大脑皮层表面三角形网格进行重采样以保持输入刚性配准图神经网络中的三角形网格与输入通道中的大脑球化模版具有相同的顶点数量和三角面片连接关系，同时降低了输入刚性配准图神经网络的处理量。
在一个示例中，如图4所示，设将球化后的大脑皮层表面三角形网格Moving(M)的属性特征采样通过重心差值法至大脑球化模版Fixed(F)(例如fsaverage6模版)上，。所述重心差值法包括以下步骤：
i：对大脑球化模版Fixed(F)的一个顶点在球化后的大脑皮层表面三角形网格Moving(M)的顶点中找到三个与之最接近的顶点。
ii：由第i步所使用的大脑球化模版Fixed(F)上的顶点与所找到的三个顶点中的任意两个顶点构成三角形，由此构成了三个三角形且每个三角形对应三角形网格Moving(M)中未参与构成该三角形的顶点分别计算三个三角形的面积，之后将三个面积进行归一化，从而将三个三角形的面积转换为三个权重，该权重与未参与组成对应三角形的三角形网格Moving(M)上的顶点相对应；例如大脑球化模版Fixed(F)上的顶点A(未示出)与所找到的三个顶点B、C和D中的任意两个顶点构成三角形，即可以构成三角形ABC、ABD、ACD，之后分别求得上述三个三角形的面积并归一化，得到权重W1、W2和W3，其中权重W1是三角形ABC的面积归一化的结果，其与之对应的顶点为D，权重W2是三角形ABD的面积归一化的结果，其与之对应的顶点为C，权重W3是三角形ACD的面积归一化的结果，其与之对应的顶点为B。
iii：通过使用第ii步求得的权重和三角形网格Moving(M)中对应的三个顶点的属性值加权求得模版Fixed(F)上的顶点对应的属性值，即完成模版Fixed(F)上的一个顶点的差值。例如，设顶点B、C和D的属性值分别为P1、P2和P3，模版Fixed(F)上的顶点A对应的属性值P0，则P0＝W3P1+W2P2+W1P3。
重复步骤i～iii，遍历三角形网格Moving(M)的顶点，即可完成采样。
在一个示例中，刚性配准图神经网络的输入为大脑皮层表面重采样多边形网格。在一个示例中，依次连接的7个图卷积块中，最后一个图卷积块的池化层为平均池化层。剩余的图卷积块中，图池化层可以设计为最大池化层(Max Pooling)、平均池化层、重叠池化层(Overlapping Pooling)或空金字塔池化层(Spatial Pyramid Pooling)。优选地，图池化层为平均池化层。通过平均池化层的设计，可以将经过该图卷积块输出的多边形网格进行网格简化，经上一层图卷积块的图池化层简化网格后的多边形网格作为下一层图卷积块的输入。也就是说，通过平均池化能够实现边坍塌操作，即能够通过平均池化将多边形网格中的例如2个顶点变为1个顶点，由此来删除1条边，从而实现了网格的简化。所述输出层包括线性层，所述线性层的输入为经上一层神经网络中的图池化层平均池化后的多边形网格。
在一个示例中，刚性配准图神经网络模型通过非监督方法训练获得。采用fsaverage6模版作为大脑球化模版，其左脑、右脑均包含40962个顶点。在一个示例中，刚性配准图神经网络模型采用了同时输入第一脑沟回深度特征和第一曲率特征的方案，因此将网络模型的输入设置为40962×4，在进行模型训练时，训练集采用了723例数据，测试集采用了181例数据，且测试集同时用作验证集。为了与FreeSurfer软件保持兼容，采用了FreeSurfer软件的球化处理多边形网格。网络模型训练的损失函数loss(pindi，patlas)的表达式为：
loss(pindi，patlas)＝wcorrLcorr(pindi，patlas)+wL1LL1(pindi，patlas)+wL2LL2(pindi，patlas) (1)
在式(1)中，pindi表示基于刚性配准欧拉角进行刚性配准(即旋转变换)后重采样得到的个体的大脑皮层表面的第一脑沟回深度特征或第一曲率特征，该第一脑沟回深度与第一曲率特征为40962维向量；patlas表示大脑球化模版的大脑皮层表面的第二脑沟回深度特征或第二曲率特征，该第二脑沟回深度与第二曲率特征为40962维向量；Lcorr表示变换后的个体的大脑皮层表面的第一脑沟回深度特征或第一曲率特征与模版计算的相关性损失，LL1、LL2分别表示变换后的个体的大脑皮层表面的第一脑沟回深度特征或第一曲率特征与模版计算的L1损失和L2损失；wcorr、wL1、wL2分别表示损失函数loss(pindi，patlas)的权重。
在一个示例中，wcorr、wL1、wL2分别取值为1、10、1。在计算损失时，对第一脑沟回深度特征、第一曲率特征分别计算了损失并进行了加权求和，两者的权重分别是0.95和0.05。
Lcorr＝1-corr(pindi，patlas) (2)
在式(2)-(5)中，Lcorr表示变换后的个体的大脑皮层表面的第一脑沟回深度特征或第一曲率特征与模版计算的相关性损失；N表示顶点数量，N取值为40962；
分别表示个体大脑皮层表面重采样多边形网格第一脑沟回深度特征或第一曲率特征与大脑球化模版第二脑沟回深度特征与第二曲率特征的均值；pindi表示基于刚性配准欧拉角进行刚性配准(即旋转变换)后重采样得到的个体的大脑皮层表面的第一脑沟回深度特征或第一曲率特征；patlas表示大脑球化模版的大脑皮层表面的第二脑沟回深度特征或第二曲率特征；LL1、LL2分别表示变换后的个体的大脑皮层表面的第一脑沟回深度特征或第一曲率特征与模版计算的L1损失和L2损失。
在一个示例中，网络模型训练的优化器为AdamW，初始学习率设置为0.001，在训练集上最多训练200轮，模型训练结束后，选择在验证集上损失函数值最低的模型作为最终模型。
在步骤S4中，如4所示，重采样前的球化后的大脑皮层表面多边形网格基于刚性配准图神经网络所输出的欧拉角r(α，β，γ)进行旋转变换以获得刚性配准后的大脑皮层表面多边形网格Rigid Moved。如图5所示，刚性配准图神经网络输出的欧拉角r(α，β，γ)的三个角度α，β，γ分别表示依次围绕Z轴、Y轴和X轴进行旋转的角度。当然本领域技术人员可以明白，在旋转时旋转顺序也可以设置为Y→X→Z；还可以设置为X→Y→Z等，本领域技术人员可以根据需要来设置。在将重采样前的球化后的大脑皮层表面多边形网格进行基于欧拉角的旋转变换过程中，多边形网格上的每一个顶点能够基于同一个欧拉角在球面上移动，通过这样一次性地曲面移动就实现了大脑皮层表面的配准。且由于欧拉角值的大小与球面上的移动距离成正比，因此在表示球面上任意两点之间的旋转变换时，不仅能够表示0～90°的旋转变换，而且还能表示90°以上的旋转变换，也就是说当顶点的移动距离较大(例如对应的欧拉角大于或等于90°)时，欧拉角同样能够表示出来，进而能够实现每一个顶点一次性移动至配准位置。同时，通过欧拉角的这种设计旋转得到的多边形网格还非常便于施加正则约束。
在一个示例中，重采样前的球化后的大脑皮层表面多边形网格基于旋转矩阵(如图4所示的3×3矩阵)进行旋转变换，该旋转矩阵为基于欧拉角r(α，β，γ)转换得到的。欧拉角与转换矩阵的转换关系为：
r(α，β，γ)＝rZ(α)rY′(β)rX″(γ) (9)
式中，rZ(α)表示绕Z轴旋转α角度的得到的弧度；rY′(β)表示绕Y′轴旋转β角度得到的弧度，Y′轴为将Y轴按rZ(α)进行旋转后的坐标轴；rX″(γ)表示绕X″轴旋转γ角度得到的弧度，X″轴为将X轴按rZ(α)旋转，再按rY′(β)进行旋转后的坐标轴；r(α，β，γ)为由三次旋转合成的旋转矩阵，其旋转结果与依次按rZ(α)，rY′(β)，rX″(γ)三个旋转矩阵旋转等价。
对重采样前的球化后的大脑皮层表面多边形网格上的任意一点P3×1＝[x，y，z]′应用旋转矩阵r(α，β，γ)得到配准后的点P′3×1＝[x′，y′，z′]′的变换关系如下所示：
P′3×1＝r(α，β，γ)P3×1
例如，如图5所示，球面上的点vn依次围绕Z轴、Y′轴(经过第一次旋转后的Y轴)、X″轴(经过两次旋转后的X轴)旋转到v′n点的变换。
在一个示例中，将所述重采样多边形网格和球化模版进行数据标准化处理，之后将数据标准化处理后的重采样多边形网格和球化模版输入所述刚性配准图神经网络中。在一个示例中，重采样多边形网格中的各个顶点的第一脑沟回深度特征、第一曲率特征使用Z-Score进行标准正态化。同样地，大脑皮层表面球化模版使用相同方法来处理。在一个示例中，重采样多边形网格和球化模版的各个顶点坐标(即第一位置信息、第二位置信息)通过除以100进行标准化。通过数据标准化处理可以使网络模型训练过程更加稳定。
结合图1和图6所示，基于医学影像图像的配准方法还包括以下步骤：
步骤S5将刚性配准后的人体器官多边形网格输入非刚性配准图神经网络以获得在所述刚性配准后的人体器官多边形网格配准至所述球化模版时的每个顶点的非刚性配准欧拉角；
步骤S6将所述刚性配准后的人体器官多边形网格中的每个顶点基于其所对应的非刚性配准欧拉角旋转变换以获得非刚性配准后的人体器官多边形网格。
在将经过球化后的人体器官多边形网格基于同一个欧拉角旋转变换之后，再基于每个顶点的自身的非刚性配准欧拉角进行旋转变换，得到非刚性配准后的人体器官多边形网格。换句话说，球化后的人体器官多边形网格依次经历对各个顶点的统一的粗调整之后，再对各个顶点分别进行进一步的精调整，使得球化后的人体器官多边形网格的属性特征与对应人体器官的球化模版更加接近。
下面将以脑组织作为人体器官对医学影像配准方法中的步骤S5-S6(以下简称“精调整”)进行示例性说明，其他人体器官的精调整与脑组织的精调整的方法和原理相同或者相似，在此不再一一赘述。
在一个示例中，在输入非刚性配准图神经网络之前，需要依次对刚性配准后的大脑皮层表面多边形网格进行重采样和数据标准化。在一个示例中，重采样为重心差值法。在一个示例中，数据标准化为将重采样后的刚性配准后的大脑皮层表面多边形网格(即，重采样刚性配准多边形网格)中的各个顶点的第三脑沟回深度特征、第三曲率特征使用Z-Score进行标准正态化。同样地，大脑皮层表面球化模版使用相同方法来处理。在一个示例中，重采样刚性配准多边形网格和大脑皮层表面的球化模版的各个顶点坐标(即第三位置信息、第二位置信息)通过除以100进行标准化。通过数据标准化处理可以使网络模型训练过程更加稳定。
在一个示例中，如图7所示，非刚性配准图神经网络为基于U-net网络架构设计的网络模型。当然也可以采用U-net++、Runet、attention Unet以及根据U-net网络改进的其它网络作为模型框架，当然还可以采用基于图自编码器(GAE)作为模型架构。
在一个示例中，非刚性配准图神经网络包括编码路径、解码路径和跳跃连接。所述编码路径包括多个第一图卷积块，所述多个第一图卷积块中的每一个第一图卷积块包括至少一个第一图卷积子块和用于下采样的第一图池化层，所述至少一个第一图卷积子块中的每一个第一图卷积子块包括第一图卷积层、第一批标准化层和第一激活函数(例如，LeakyRelu函数)。
所述解码路径包括多个第二图卷积块，所述多个第二图卷积块中的每一个第二图卷积块包括至少一个第二图卷积子块和用于上采样的第二反图卷积层，所述至少一个第二图卷积子块中的每一个第二图卷积子块包括第二图卷积层、第二批标准化层和第二激活函数(例如，LeakyRelu函数)。
同一层的第一图卷积块和第二图卷积块之间通过跳线进行跳跃连接，且同一层次的第一图卷积块的经过下采样获得的特征图与第二图卷积块的经过上采样获得的特征图进行拼接作为下一层的第二图卷积块的输入。即通过跳线将同层的编码结构和解码结构连接，这样可以补充不同层次模型需要的细节。
在一个示例中，多个第一图卷积块、多个第二图卷积块均可以设计为5个，也就是说，输入的多边形网格需要经历5次下采样和5次上采样。本领域技术人员可以明白，第一图卷积块和第二图卷积块的数量需要根据输入的多边形网格的顶点数目来决定，例如可以根据需要设计为3个第一、第二图卷积块，也可以设计为6个第一、第二图卷积块。
在一个示例中，第一图卷积块中设置有2个第一图卷积子块，第二图卷积块中设置有2个第二图卷积子块，当然本领域技术人员可以根据需要将第一、第二图卷积子块的个数均设置为1个或者3个。
在一个示例中，非刚性配准图卷积神经网络输入的多边形网格具有第三脑沟回深度特征或第三曲率特征时，每个顶点的网络输入通道数量设置为2个，其中一个通道用于重采样刚性配准多边形网格的第三位置信息(即顶点坐标)和第三脑沟回深度特征或者第三位置信息和第三曲率特征的输入，另一个通道用于大脑皮层表面球化脑模版的第二位置信息和第二脑沟回深度特征或者第二位置信息和第二曲率特征的输入。
当非刚性配准图卷积神经网络输入的多边形网格具有第三脑沟回深度特征和第三曲率特征时，每个顶点的网络输入通道数量设置为4个，其中2个通道分别用于重采样刚性配准多边形网格的第三位置信息与第三脑沟回深度特征和第三位置信息与第三曲率特征的输入，另外2个通道分别用于大脑皮层表面球化脑模版的第二位置信息与第二脑沟回深度特征和第二位置信息与第二曲率特征的输入。
在一个示例中，第三脑沟回深度特征与第一脑沟回深度特征相同，第三曲率特征与第一曲率特征相同。也就是说，不论是刚性配准多边形网格还是非刚性配准多边形网格，在配准过程中改变的仅是球化后的人体器官(例如脑组织)多边形网格上的各个顶点的位置，而每个顶点的属性特征并没有发生改变。
在一个示例中，大脑皮层表面球化模版为刚性配准时使用的球化模版。因此，在重采样时，将刚性配准后的大脑皮层表面多边形网格的顶点数量重采样为40962个。当输入通道设置为4个时，输入非刚性配准图神经网络的多边形网格的特征数量为40962×4。也就是说，网络模型的输入通道数量基于需要输入的顶点的属性特征的类别数量来决定。
在一个示例中，编码路径中通过逐层卷积、逐层池化显著增加了图卷积操作的感受野，且通过这样的设计还不会引起计算量和显存的爆炸。
在一个示例中，解码路径中的逐层的第一图卷积层基于个体的刚性配准重采样多边形网格与球化模版中的每个顶点的属性特征进行图卷积操作，实现了逐层对感兴趣信息的添加，进而逐层完成对各个顶点的欧拉角的预测。
在一个示例中，第一图池化层将同一网络层中的至少一个第一图卷积子块中的最后一个第一图卷积子块输出的多边形网格进行边坍塌(Edge Collapse)操作以简化激活后的多边形网格中的网格。例如，当同一层网络中的至少一个第一图卷积子块为2个第一图卷积子块时，第一图池化层对第二个第一图卷积子块中的第一激活函数激活后的多边形网格进行边坍塌操作。
在一个示例中，由于边坍塌过程具有可逆性，因此每一层网络的第一图池化层在进行边坍塌操作过程中，顺序地将边坍塌信息进行存储，由此使得多边形网格输入解码路径后，能够顺序地恢复、重构多边形网格。
在一个示例中，在解码路径中，第二反图卷积层同一网络层中的至少一个第二图卷积子块中的最后一个第二图卷积子块输出的多边形网格进行曲面细分，以加密激活后的多边形网格。
在一个示例中，第二反卷积层将输入其的多边形网格例如通过增加新的顶点的方式来对多边形网格进行曲面细分，从而增加图结构特征数量。
在一个示例中，第二反卷积层将输入的多边形网格的顶点数量按照1:5～9的比例关系扩展，例如1:5的关系进行扩展。在一个示例中，第二反卷积层依据解码路径中第一图池化层所存储的边坍塌信息(例如各个顶点被删除的顺序)来依次恢复各个顶点，重构多边形网格。
在一个示例中，结合图6-7所示，刚性配准重采样多边形网格Moving(M’)经过非刚性配准图神经网络的预测获得了刚性配准重采样多边形网格中的每一个顶点待调整的非刚性配准欧拉角。例如，刚性配准重采样多边形网格上的顶点数量为40962个，则需要非刚性配准图神经网络输出的欧拉角为40962个。由于欧拉角是从3个维度的角度来描述的，因此需要非刚性配准图神经网络的每个顶点的输出通道设置为3个，这样通过非刚性配准图神经网络输出的图结构的特征数量为40962×3。
在一个示例中，非刚性配准图神经网络模型通过非监督训练获得。在训练中，模型的训练集采用了723例数据，模型测试集采用了181例数据，并采用球面刚性配准并进行球面重采样后的个体的球化大脑皮层多边形网格、FreeSurfer的fsaverage6模版作为训练输入，且输入的属性特征包括脑沟回深度特征和曲率特征。模型训练的损失函数loss(pindi，patlas，Smoved)的表达式为：
loss(pindi，patlas，Smoved)＝wcorrLcorr(pindi，patlas)+wL1LL1(pindi，patlas)+wL2LL2(pindi，patlas)+wsmoothLsmooth(Smoved) (10)
在式(10)中，pindi表示基于非刚性配准欧拉角进行非刚性配准(即旋转变换)后重采样多边形网格的脑沟回深度特征和/或曲率特征，在训练过程中采用的40962维的向量；patlas表示球化模版的脑沟回深度特征和/或曲率特征，在训练过程中采用的40962维的向量；Smoved表示基于非刚性配准欧拉角进行非刚性配准后的大脑皮层表面多边形网格；Smoving每个顶点表示基于非刚性配准欧拉角进行非刚性配准后的三角形网格；Lcorr表示非刚性配准后的个体大脑皮层多边形网格的脑沟回深度特征和/或曲率特征与模版计算的相关损失；LL1、LL2分别表示非刚性配准后的个体大脑皮层多边形网格的脑沟回深度特征和/或曲率特征与模版计算的L1损失和L2损失，Lsmooth表示旋转场平滑损失；wcorr、wL1、wL2、wsmooth分别是损失函数的权重。
在一个示例中，wcorr、wL1、wL2、wsmooth的取值选择为1、10、1、1000。在一个示例中，在计算损失时，对第三脑沟回深度特征、第三曲率特征分别进行了计算并加权求和，两者的权重分别选择为0.95，0.05。
在一个示例中，旋转场平滑损失采用三角形网格版本的拉普拉斯矩阵计算方法进行计算，该项损失的引入是为了对非线性变换进行约束，避免球形三维曲面在旋转变换中发生局部折叠的问题。
如图8A-图8B所示，图8A是刚性配准重采样多边形网格的曲面状态。在图8B中，当未使用拉普拉斯算子进行约束时，刚性配准重采样多边形网格在基于非刚性配准欧拉角旋转变换后，将发生折叠现象。在图8C中，通过拉普拉斯算子的约束，刚性配准重采样多边形网格在基于非刚性配准欧拉角旋转变换后，未发生折叠现象，由此可以看出，通过拉普拉斯算子的约束能够有效避免折叠现象的发生，由此保障了个体大脑皮层和模版球面位置映射的一一对应关系。
在一个示例中，Lcorr(pindi，patlas)、LL1(pindi，patlas)和LL2(pindi，patlas)的计算方法的原理与上文中式(2)-(5)相同，在此不再一一赘述。Lsmooth(Smoved)采用了python第三方库pytroch3d中的计算来实现。
在一个示例中，网络模型训练的优化器为AdamW，初始学习率设置为0.003，在训练集上最多训练200轮，待模型训练结束后，选择在验证集上损失函数值最低的模型作为最终模型。
在一个示例中，我们将使用FreeSurfer软件进行大脑皮层表面的非刚性配准与本发明一个实施例提供的医学影像图像配准方法(依次进行了刚性配准和非刚性配准)进行了比较。测试指标是在测试集181例数据上计算并统计的，并采用了皮尔逊相关系数(PCC，Pearson correlation coefficient)、平均绝对误差(MAE，Mean Absolute Error)和Dice系数作为配准评价指标来衡量配准后个体大脑皮层表面脑沟回深度特征(sulc)和曲率特征(curv)分布以及脑结构分区一致性与模版的接近程度。指标值均以“均值±标准差”的形式给出，且PCC数值越高表示与脑模板的特征分布的相似性越高，MAE数值越高表示与脑模板的特征分布的差异越大；Dice系数越接近于1则代表分区一致性高。配准效果指标对比如表一所示。
表一 配准效果指标对比表
由于FreeSurfer软件无法在GPU上运行，因此只做了其与本发明非刚性配准方法在单纯的CPU环境下的对比。
从表中可以得出：在单纯的CPU环境下和单纯的GPU环境下，本发明的一个示例的刚性配准方法在处理速度(即运行时间)上相比现有的FreeSurfer软件的非刚性配准方法具有很大优势，在单纯CPU环境下约有14637倍的速度提升。且本发明另一示例的非刚性配准方法在处理速度(即运行时间)上相比现有的FreeSurfer软件的非刚性配准方法同样具有很大优势，在单纯CPU环境下约有54倍的速度提升。由此可见，当依次使用本发明一个实施例提供的刚性配准方法和非刚性配准方法进行配准时，相比于现有的FreeSurfer软件的非刚性配准方法在保证了配准精度的同时，还能够大大缩短处理时间。
在一个示例中，刚性配准的大脑皮层表面多边形网格Moving(M’)的每一个顶点基于其对应的非配准欧拉角转换的旋转矩阵N×3×3(如图6所示)旋转变换以实现非刚性配准得到非刚性配准的大脑皮层表面多边形网格Non-rigid Moved。其中，N为刚性配准的大脑皮层表面多边形网格的顶点数量。在刚性配准的大脑皮层表面多边形网格进行基于欧拉角的旋转变换过程中，多边形网格上的每一个顶点能够基于自身对应的欧拉角在球面上移动，通过每个顶点自身进行一次性地曲面移动就实现了大脑皮层表面的精确配准。且由于欧拉角值的大小与球面上的移动距离成正比，因此在表示球面上任意两点之间的旋转变换时，不仅能够表示0～90°的旋转变换，而且还能表示90°以上的旋转变换，也就是说当顶点的移动距离较大(例如对应的欧拉角大于或等于90°)时，欧拉角同样能够表示出来，进而能够实现每一个顶点基于自身的欧拉角即可一次性移动至精确配准位置。
在一个示例中，欧拉角与旋转矩阵的转换方法和原理如上文所示的式(6)-(9)所表达的方法和原理。
对刚性配准的大脑皮层表面多边形网格Moving(M’)上的任意一点K3×1＝[m，n，l]′应用式(9)的旋转矩阵r(α，β，γ)得到配准后的点K′3×1＝[m’，n′，l′]′的变换关系如下所示：
K′3×1＝r(α，β，γ)K3×1
在一个示例中，根据本发明的另一实施例提供了一种可读存储介质。本发明的实施例的“可读存储介质”是指参与向处理器提供程序或指令以供执行的任何介质。所述介质可以采用多种形式，包括但不限于非易失性介质、易失性介质和传输介质。非易失性介质包括例如光盘或磁盘，诸如存储设备。易失性介质包括动态存储器，诸如主存储器。传输介质包括同轴电缆、铜线和光纤，包括包含总线的导线。传输介质还可以采用声波或光波的形式，诸如在射频(RF)和红外(IR)数据通信期间生成的声波或光波。可读存储介质的常见形式包括例如软盘、柔性盘、硬盘、磁带、任何其他磁介质、CD-ROM、DVD、任何其他光学介质、穿孔卡、纸带、任何具有孔图案的其他物理介质、RAM、PROM和EPROM、FLASH-EPROM、任何其他存储器芯片或盒、如下所述的载波、或计算机可从其中进行读取的任何其他介质。
该可读存储介质上存储程序或指令，所述程序或指令被处理器执行时以执行上述的基于医学影像图像的配准方法。
在一个示例中，根据本发明的还一实施例提供了一种电子设备。该电子设备(未示出)包括处理器(未示出)和存储器(未示出)。存储器上存储程序，该程序被处理器执行时，可以执行上述任一示例的基于医学影像图像的配准方法以实现三维模型的配准。
在一个示例中，处理器可以为微处理器，例如图形处理器(GPU)、中央处理器(CPU)、数字信号处理器(Digital Signal Processor，DSP)等通用处理器。在一个示例中，处理器还可以为通过硬件电路实现的微处理器核，例如通过可重新配置逻辑在硬件逻辑组件中实现的微处理器核，所述硬件逻辑组件包括现场可编程门阵列(FPGA)、复杂可编程逻辑器件(CPLD)、程序专用集成电路(ASIC)、程序专用标准产品(ASSP)以及片上系统(SOC)以及等。
在一个示例中，处理器还可以是虚拟处理器，该虚拟处理器可以为带有英特尔x86处理器特性的虚拟处理器，还可以可为具有PowerPC处理器的特性的虚拟处理器。优选地，处理器为图形处理器。在一个示例中，所述处理器可以为单核处理器，也可以为多核处理器。
在一个示例中，存储器包括易失性存储器(即随机存取存储器)和非易失性存储器。易失性存储器包括主存储器、高速缓存(Cache)等，非易失性存储器包括辅助存储器等。在一个示例中，存储器可以设置为远程存储器，远程存储器可以通过网络(有线网络或无线网络)连接至处理器。所述网络包括但不限于广域网、局域网、城域网、个域网、互联网、卫星通信网及其任意组合。
在一个示例中，处理器执行基于存储器中获取的程序创建对应的任务线程，并执行线程。在一个示例中，处理器基于内存中的读取指令从外存中获取程序以创建对应的任务线程并执行线程。上述程序用于实现基于医学影像图像的配准方法。
尽管此处所述的主题是在结合操作系统和应用程序在计算机系统上的执行而执行的一般上下文中提供的，但本领域技术人员可以认识到，还可结合其他类型的程序模块来实现。一般而言，程序模块包括执行特定任务或实现特定抽象数据类型的例程、程序、组件、数据结构和其他类型的结构。本领域技术人员可以明白，结合本文中任何一个示例所描述的方法步骤，均能够以电子硬件、或者计算机软件和电子硬件的结合来实现。以硬件或者软件方式来执行主要取决于技术方案的特定应用和设计约束条件。本领域技术人员可以对每个特定的应用来使用不同方法以实现所描述的功能，但是这种实现不应当认为超出了本申请的范围。
当方法步骤以软件功能单元的形式实现并作为独立的产品销售或使用时，可以存储在一个计算机可读存储介质中。因此，本发明的技术方案本质上或者说对原有技术做出贡献的部分或者该技术方案的部分可以以软件产品的形式体现出来，该计算机软件产品存储在一个存储介质中，包括若干指令用以使得一台计算机设备(可以是个人计算机，服务器，或者网络设备等)执行本发明各个示例所述方法的全部或部分步骤。
根据本发明各个示例提供的基于医学影像图像的配准方法、可读存储介质及电子设备具有以下优点中的至少一个：
(1)本发明提供的基于医学影像图像的配准方法、可读存储介质及电子设备能够通过刚性配准图神经网络的设计，能够快速实现基于医学影像图像的三维模型配准，由此显著缩短了基于医学影像图像的三维模型的配准时间，大大提升了基于医学影像图像的三维模型的处理效率；
(2)本发明提供的基于医学影像图像配准方法、可读存储介质及电子设备通过刚性配准图神经网络和非刚性配准图神经网络的设计，在配准过程中依次通过粗配准对每个顶点位置的统一调整和精配准对每个顶点的位置分别进行的精细调整，从而能够获得更加精准的大脑皮层表面配准三维模型。
虽然本总体发明构思的一些实施例已被显示和说明，本领域普通技术人员将理解，在不背离本总体发明构思的原则和精神的情况下，可对这些实施例做出改变，本发明的范围以权利要求和它们的等同物限定。
Claims (18)
1.一种基于医学影像图像的配准方法，包括以下步骤：
根据所述医学影像图像进行三维模型重建获得多边形网格和所述多边形网格中的每个顶点的属性特征；
将所述多边形网格膨胀并映射至球面以获得具有所述属性特征的球化后的多边形网格；
将所述球化后的多边形网格和所对应的球化模版输入刚性配准图神经网络中，以获得所述球化后的多边形网格配准至所述所对应的球化模版的刚性配准欧拉角；
将所述球化后的多边形网格基于所述刚性配准欧拉角旋转变换以获得刚性配准后的多边形网格。
2.根据权利要求1所述的基于医学影像图像的配准方法，其中，
所述球化后的多边形网格上的所有顶点基于同一个刚性配准欧拉角旋转变换。
3.根据权利要求2所述的基于医学影像图像的配准方法，其中，
所述医学影像图像为人体器官的医学影像图像，
所述人体器官为脑组织，所述脑组织包括大脑、小脑和脑干中的至少一种。
4.根据权利要求3所述的基于医学影像图像的配准方法，其中，
所述属性特征包括每个顶点的第一位置信息和每个顶点处的第一脑沟回深度特征，
所述刚性配准图神经网络的输入通道分别用于输入具有所述第一脑沟回深度特征与所述第一位置信息的球化后的多边形网格和具有每个模版顶点的第二脑沟回深度特征与第二位置信息的球化模版。
5.根据权利要求1-4中任一项所述的基于医学影像图像的配准方法，其中，
所述属性特征包括每个顶点的第一曲率特征，
所述刚性配准图神经网络的输入通道还分别用于输入具有第一曲率特征的球化后的多边形网格和具有每个模版顶点的第二曲率特征的球化模版。
6.根据权利要求3所述的基于医学影像图像的配准方法，其中，
所述刚性配准图神经网络包括多个图卷积块和输出层，所述多个图卷积块中的每一个图卷积块均包括至少一个图卷积子块和图池化层，所述至少一个图卷积子块中的每一个图卷积子块包括图卷积子层、批标准化层和激活函数，所述图卷积子层对每个顶点和与每个顶点中的同一顶点相邻连接的其它顶点进行卷积运算。
7.根据权利要求6所述的基于医学影像图像的配准方法，其中，
所述多边形网格为三角形网格，
在获得具有所述属性特征的球化后的多边形网格之后，将所述球化后的多边形网格进行球面重采样以获得与所述所对应的球化模版的顶点数量相同的重采样多边形网格，输入所述刚性配准图神经网络中的球化后的多边形网格为所述重采样多边形网格，
在球化后的多边形网格基于所述刚性配准欧拉角旋转变换过程中，基于刚性配准欧拉角旋转变换的多边形网格为重采样前的所述球化后的多边形网格。
8.根据权利要求7所述的基于医学影像图像的配准方法，其中，
将所述重采样多边形网格和球化模版进行数据标准化处理，之后将数据标准化处理后的重采样多边形网格和球化模版输入所述刚性配准图神经网络中。
9.根据权利要求3-8中任一项所述的基于医学影像图像的配准方法，其中，
在获得所述刚性配准后的多边形网格之后，所述基于医学影像图像的配准方法还包括以下步骤：
将刚性配准后的多边形网格输入非刚性配准图神经网络以获得在所述刚性配准后的人体器官多边形网格配准至所述球化模版时的每个顶点的非刚性配准欧拉角；
将所述刚性配准后的多边形网格中的每个顶点基于其所对应的非刚性配准欧拉角旋转变换以获得非刚性配准后的人体器官多边形网格。
10.根据权利要求9所述的基于医学影像图像的配准方法，其中，
输出的非刚性配准欧拉角的数量与刚性配准后的人体器官多边形网格中的顶点数量相同。
11.根据权利要求9所述的基于医学影像图像的配准方法，其中，
所述刚性配准后的人体器官多边形网格的属性特征包括旋转后的每个顶点的第三位置信息和旋转后的每个顶点处的第三脑沟回深度，
所述非刚性配准图神经网络的输入通道分别用于输入具有第三脑沟回深度特征与第三位置信息的刚性配准后的人体器官多边形网格和具有每个模版顶点的第二脑沟回深度特征与每个模版顶点的第二位置信息的球化模版。
12.根据权利要求11所述的基于医学影像图像的配准方法，其中，
所述第三脑沟回深度与第一脑沟回深度相同。
13.根据权利要求12所述的基于医学影像图像的配准方法，其中，
所述刚性配准后的人体器官多边形网格的属性特征还包括旋转后的每个顶点处的第三曲率特征，
所述非刚性配准图神经网络的输入通道还分别用于输入具有第三曲率特征的人体器官多边形网格和具有第二曲率特征的球化模版。
14.根据权利要求13所述的基于医学影像图像的配准方法，其中，
所述第三曲率特征与所述第一曲率特征相同。
15.根据权利要求9所述的基于医学影像图像的配准方法，其中，
所述非刚性配准图神经网络包括编码路径、解码路径和跳跃连接，所述编码路径包括多个第一图卷积块，所述多个第一图卷积块中的每一个第一图卷积块包括至少一个第一图卷积子块和用于下采样的第一图池化层，所述至少一个第一图卷积子块中的每一个第一图卷积子块包括第一图卷积层、第一批标准化层和第一激活函数，
所述解码路径包括多个第二图卷积块，所述多个第二图卷积块中的每一个第二图卷积块包括至少一个第二图卷积子块和用于上采样的第二反图卷积层，所述至少一个第二图卷积子块中的每一个第二图卷积子块包括第二图卷积层、第二批标准化层和第二激活函数，
同一层的第一图卷积块和第二图卷积块之间通过跳线进行跳跃连接，且所述同一层次的第一图卷积块的经过下采样获得的特征图与第二图卷积块的经过上采样获得的特征图进行拼接作为下一层的第二图卷积块的输入。
16.根据权利要求15所述的基于医学影像图像的配准方法，其中，
在编码路径中，第一图池化层将经过第一激活函数激活后的多边形网格进行边坍塌操作以简化激活后的多边形网格，
在解码路径中，第二反图卷积层将经过第二激活函数激活后的多边形网格进行曲面细分以加密激活后的多边形网格。
17.一种可读存储介质，其中，
所述可读存储介质上存储程序或指令，所述程序或指令被处理器执行时以执行权利要求1-16中任一项所述的基于医学影像图像的配准方法。
18.一种电子设备，所述电子设备包括处理器和存储器，其中，
所述存储器上存储程序，所述程序被所述处理器执行时以执行权利要求1-16中任一项所述的基于医学影像图像的配准方法。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202211198336.4A CN115511932B (zh) | 2022-09-29 | 2022-09-29 | 基于医学影像图像的配准方法、可读存储介质及电子设备 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202211198336.4A CN115511932B (zh) | 2022-09-29 | 2022-09-29 | 基于医学影像图像的配准方法、可读存储介质及电子设备 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN115511932A true CN115511932A (zh) | 2022-12-23 |
| CN115511932B CN115511932B (zh) | 2024-02-13 |
Family
ID=84508388
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202211198336.4A Active CN115511932B (zh) | 2022-09-29 | 2022-09-29 | 基于医学影像图像的配准方法、可读存储介质及电子设备 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN115511932B (zh) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116778020A (zh) * | 2023-08-17 | 2023-09-19 | 山东大学 | 一种基于深度学习的柔性超声聚束成像方法及系统 |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111457930A (zh) * | 2020-04-02 | 2020-07-28 | 武汉中海庭数据技术有限公司 | 一种利用车载Lidar与无人机联合的高精度测图定位方法 |
| CN113643336A (zh) * | 2021-07-26 | 2021-11-12 | 之江实验室 | 基于球极坐标系深度神经网络的三维图像刚性匹配方法 |
| US20220122250A1 (en) * | 2020-10-19 | 2022-04-21 | Northwestern University | Brain feature prediction using geometric deep learning on graph representations of medical image data |
| CN114972366A (zh) * | 2022-07-27 | 2022-08-30 | 山东大学 | 基于图网络的大脑皮层表面全自动分割方法及系统 |
-
2022
- 2022-09-29 CN CN202211198336.4A patent/CN115511932B/zh active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111457930A (zh) * | 2020-04-02 | 2020-07-28 | 武汉中海庭数据技术有限公司 | 一种利用车载Lidar与无人机联合的高精度测图定位方法 |
| US20220122250A1 (en) * | 2020-10-19 | 2022-04-21 | Northwestern University | Brain feature prediction using geometric deep learning on graph representations of medical image data |
| CN113643336A (zh) * | 2021-07-26 | 2021-11-12 | 之江实验室 | 基于球极坐标系深度神经网络的三维图像刚性匹配方法 |
| CN114972366A (zh) * | 2022-07-27 | 2022-08-30 | 山东大学 | 基于图网络的大脑皮层表面全自动分割方法及系统 |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116778020A (zh) * | 2023-08-17 | 2023-09-19 | 山东大学 | 一种基于深度学习的柔性超声聚束成像方法及系统 |
| CN116778020B (zh) * | 2023-08-17 | 2023-11-14 | 山东大学 | 一种基于深度学习的柔性超声聚束成像方法及系统 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN115511932B (zh) | 2024-02-13 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109754394B (zh) | 三维医学图像处理装置及方法 | |
| Eppenhof et al. | Pulmonary CT registration through supervised learning with convolutional neural networks | |
| CN113344799A (zh) | 用于使用深度神经网络减少医学图像中的有色噪声的系统和方法 | |
| Cheng et al. | A morphing-Based 3D point cloud reconstruction framework for medical image processing | |
| Akkus et al. | Robust brain extraction tool for CT head images | |
| Habijan et al. | Whole heart segmentation from CT images using 3D U-net architecture | |
| CN110251231A (zh) | 超声三维重建的方法及装置 | |
| JP7349158B2 (ja) | 機械学習装置、推定装置、プログラム及び学習済みモデル | |
| CN115511932B (zh) | 基于医学影像图像的配准方法、可读存储介质及电子设备 | |
| CN112884819A (zh) | 一种影像配准及神经网络的训练方法、装置和设备 | |
| Chai et al. | A voxel‐based finite element model for the prediction of bladder deformation | |
| Chatterjee et al. | Ddos-unet: Incorporating temporal information using dynamic dual-channel unet for enhancing super-resolution of dynamic mri | |
| US20240105328A1 (en) | Neural network simulator for ultrasound images and clips | |
| CN109741439A (zh) | 一种二维mri胎儿图像的三维重建方法 | |
| Menchón-Lara et al. | Fast 4D elastic group-wise image registration. Convolutional interpolation revisited | |
| CN115861399B (zh) | 基于医学影像图像的配准方法、可读存储介质及电子设备 | |
| Amirrajab et al. | A framework for simulating cardiac MR images with varying anatomy and contrast | |
| Harders | Surgical scene generation for virtual reality-based training in medicine | |
| Zhong et al. | 3D‐2D Deformable Image Registration Using Feature‐Based Nonuniform Meshes | |
| CN108064148A (zh) | 用于体外受精的图像引导的胚胎转移 | |
| Sang et al. | Inter-phase 4D cardiac MRI registration with a motion prior derived from CTA | |
| Qi et al. | Multi-task MR imaging with iterative teacher forcing and re-weighted deep learning | |
| CN111476888A (zh) | 一种基于三维空间体拟合的医学图像层间插值方法，设备及可读存储介质 | |
| Cai et al. | [Retracted] Detection of 3D Arterial Centerline Extraction in Spiral CT Coronary Angiography | |
| Chen et al. | Multi-resolution statistical shape models for multi-organ shape modelling |
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 |