一种扩散张量成像参数量化方法、系统、设备和介质

技术领域

本发明属于磁共振技术领域，尤其涉及一种扩散张量成像参数量化方法、系统、设备和介质。

背景技术

扩散磁共振成像(Diffusion Magnetic Resonance Imaging，dMRI)可以对活体组织中水分子的扩散运动进行无创映射，是临床重要检查手段之一。扩散张量成像(Diffusion tensor imaging,DTI)是一种常见的扩散磁共振成像方式，需要采集至少一副非扩散加权图像和六幅扩散加权(diffusion-weighted,DW)图像，目前已广泛应用于临床量化诊断以及脑白质纤维束可视化。b值是DW图像的重要参数，含义为扩散敏感系数，其大小影响着DW图像信号的衰减，b值越高，DW图像信号越弱。因此DW图像的噪声较常规的结构磁共振图像更为严重。噪声的存在会模糊图像细节，降低图像质量，导致DTI扩散张量场估计不准确，从而导致平均扩散率(mean diffusivity,MD)、各向异性分数(fractionalanisotropy,FA)等DTI量化参数计算不准确，影响临床诊断和科学研究。

因此，降低DW图像的噪声是十分必要的。在临床上，通常使用重复采集多次DW图像取平均值的方法来提高DW图像的信噪比，然而该方法会延长采集时间，进而增加采集成本和被试者在采集过程中发生运动的概率。为此，研究学者们相继提出了大量后处理技术以降低噪声对DTI参数量化的影响。这些后处理技术包括传统方法和深度学习方法。其中，MPPCA、NLM、BM3D等对DW图像进行去噪的方法属于传统方法。传统方法的缺点主要为计算时间长，算法复杂度高。现阶段的深度学习方法为：根据有噪声的DW图像，利用深度学习网络来预测无噪的DW图像或DTI扩散张量场或DTI量化参数。现有深度学习方法较传统方法效果好，速度快。然而现有的深度学习方法受限于指定的采集方案(扩散编码方向指向、方向个数、b值大小)，泛化能力差，难以在临床上推广使用。

发明内容

本发明的目的在于提供一种扩散张量成像参数量化方法、系统、设备和介质，能够提高DTI扩散张量场估计的准确性，进而提高DTI量化参数的准确性，为精准医疗提供可靠的量化依据。

本发明是通过以下技术方案实现的：

一种扩散张量成像参数量化方法，包括以下步骤：

获取dMRI数据{S

根据扩散编码方向g

对dMRI数据进行数据变换，得到扩散衰减图像，并根据扩散衰减图像和球谐函数的复共轭Y

将球谐系数图输入已构建的深度学习网络，得到DTI扩散张量场；

根据DTI扩散张量场，计算得到DTI量化参数图。

进一步地，对dMRI数据进行数据变换的步骤包括：

采用公式(1)对dMRI图像进行数据变换，公式(1)为：

式中，{S(g

进一步地，根据扩散衰减图像和球谐函数的复共轭Y

基于扩散衰减图像和球谐函数的复共轭Y

式中，

进一步地，深度学习网络包括：

编码器，用于对球谐系数图特征进行提取，得到特征图像；

解码器，用于对编码器提取到的特征图像进行解码，得到解码特征图；

跳跃连接层，将编码器输出的特征图像与解码器输出的解码特征图像连接在一起。

进一步地，DTI量化参数图至少包括FA图和MD图；

根据DTI扩散张量场，计算得到DTI量化参数图的步骤包括：

对DTI扩散张量场的扩散张量D进行特征值分解，得到特征值λ

根据公式(3)计算FA图，根据公式(4)计算MD图，其中，公式(3)如下：

公式(4)如下：

进一步地，深度学习网络训练过程包括以下步骤：

(1)训练数据准备，获取多例dMRI数据，并对每一例dMRI数据进行数据处理，其中数据处理的过程如下：

(1-1)挑选出b值为预设数值的扩散加权图像以及对应的非扩散加权图像，并进行预处理。预处理过程包括图像去噪，去吉布斯伪影以及偏场校正；

(1-2)对于预处理后的扩散加权图像以及对应的非扩散加权图像，采用加权线性最小二乘法进行逐像素拟合得到DTI扩散张量场，以得到的DTI扩散张量场作为dMRI数据的真实扩散张量场

(1-3)对DTI扩散张量场应用DTI模型，使用预设的扩散梯度表进行dMRI数据重建，得到扩散加权图像，DTI模型的公式为：

式中，s

(1-4)对于重建后的dMRI数据，加入莱斯噪声：

式中，I为重建后真实无噪的dMRI数据，N

(1-5)对于加入莱斯噪声的dMRI数据I

(2)深度学习网络训练；

(2-1)在多例dMRI数据中随机选择n例作为训练集，剩余部分作为验证集；

(2-2)采用训练集对深度学习网络进行训练，直到网络收敛，同时验证集损失达到稳定。

进一步地，深度学习网络的损失函数为：

式中，

本发明还提供了一种扩散张量成像参数量化系统，包括：

获取模块，用于获取dMRI数据{S

求取模块，用于根据扩散编码方向g

处理模块，用于对dMRI数据进行数据变换，得到扩散衰减图像，并根据扩散衰减图像和球谐函数的复共轭Y

输入模块，用于将球谐系数图输入已构建的深度学习网络，得到DTI扩散张量场；

计算模块，用于根据DTI扩散张量场，计算得到DTI量化参数图。

本发明还公开了一种电子设备，包括存储器和处理器，存储器存储有计算机程序，处理器执行计算机程序时实现上述任一项方法的步骤。

本发明还公开了一种计算机可读存储介质，其上存储有计算机程序，计算机程序被处理器执行时实现上述任一项的方法的步骤。

相比于现有技术，本发明的有益效果为：使用深度学习网络和球谐函数对扩散张量成像参数进行量化，较传统的方法计算速度快，可以得到更准确的DTI量化参数；并且较其他基于深度学习网络的扩散张量参数量化方法，不受限于指定的采集方案，能够准确的估计出不同采集方案的dMRI数据的扩散张量场，泛化能力强。

附图说明

图1为本发明扩散张量成像参数量化方法的步骤流程图；

图2为本发明扩散张量成像参数量化方法的整体框架图；

图3为本发明扩散张量成像参数量化方法中深度学习网络的结构示意图；

图4为实施例1仿真数据的实验方法对比结果，展示了DTI扩散张量场的六个元素图、FA图和MD图以及对应的误差图；

图5为实施例2高分辨率大脑真实数据的实验方法对比结果，展示了DTI扩散张量场的六个元素图、FA图和MD图以及局部放大图；

图6为本发明扩散张量成像参数量化系统的模块示意图；

图7为本发明电子设备一实施例的结构示意图；

图8为本发明计算机可读存储介质一实施例的结构示意框图。

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。通常在此处附图中描述和示出的本发明实施例的组件可以以各种不同的配置来布置和设计。

因此，以下对在附图中提供的本发明的实施例的详细描述并非旨在限制要求保护的本发明的范围，而是仅仅表示本发明的选定实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

应注意到：相似的标号和字母在下面的附图中表示类似项，因此，一旦某一项在一个附图中被定义，则在随后的附图中不需要对其进行进一步定义和解释。同时，在本发明的描述中，术语“第一”、“第二”等仅用于区分描述，而不能理解为指示或暗示相对重要性。

需要说明的是，在本文中，诸如第一和第二等之类的关系术语仅仅用来将一个实体或者操作与另一个实体或操作区分开来，而不一定要求或者暗示这些实体或操作之间存在任何这种实际的关系或者顺序。而且，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、方法、物品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、方法、物品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括所述要素的过程、方法、物品或者设备中还存在另外的相同要素。

在本发明的描述中，需要说明的是，术语“上”、“下”、“内”、“外”等指示的方位或位置关系为基于附图所示的方位或位置关系，或者是该发明产品使用时惯常摆放的方位或位置关系，仅是为了便于描述本发明和简化描述，而不是指示或暗示所指的装置或元件必须具有特定的方位、以特定的方位构造和操作，因此不能理解为对本发明的限制。

请参阅图1和图2，图1为本发明扩散张量成像参数量化方法的步骤流程图，图2为本发明扩散张量成像参数量化方法的整体框架图。一种扩散张量成像参数量化方法，包括以下步骤：

S1、获取dMRI数据{S

S2、根据扩散编码方向g

S3、对dMRI数据进行数据变换，得到扩散衰减图像，并根据扩散衰减图像和球谐函数的复共轭Y

S4、将球谐系数图输入已构建的深度学习网络，得到DTI扩散张量场；

S5、根据DTI扩散张量场，计算得到DTI量化参数图。

在上述步骤S1中，获取的dMRI数据{S

在上述步骤S2中，根据dMRI数据的扩散编码方向g

式中，P

在步骤S3中，对dMRI数据进行数据变换的步骤包括：

S31、采用公式(1)对dMRI数据进行数据变换，公式(1)为：

式中，{S(g

在步骤S31中，分别使用扩散加权图像S

在步骤S3中，根据扩散衰减图像和球谐函数的复共轭Y

S32、根据扩散衰减图像和球谐函数的复共轭Y

/>

式中，

在上述步骤S32中，球谐系数图是由球谐系数组成的图。每个体素对应一个球谐系数(标量)，那么这个球谐系数对应的就是一副图，称之为球谐系数图。具体地，在本实施例中，阶数l＝2，所以每个体素根据非扩散加权图像以及扩散加权图像，计算得到的球谐系数是1*6的向量，所以球谐系数图是六个图，如图2所示。另外，s(r)就是扩散衰减信号，r对应的就是相应的扩散编码方向g

在上述步骤S4中，预先构建深度学习网络，并采用预设的训练数据对深度学习网络进行训练，深度学习网络输入的是球谐系数图，深度学习网络输出的是DTI扩散张量场。因此，将计算得到的球谐系数图输入到预先构建的深度学习网络中，通过深度学习网络得到DTI扩散张量场。深度学习网络可以采用U-net网络，它是一个U型网络。

进一步地，深度学习网络包括：

编码器，用于对球谐系数图的特征进行提取，得到特征图像；

解码器，用于对编码器提取到的特征图像进行解码，得到解码特征图像；

跳跃连接层，，s(r)∈R。

深度学习网络包含一个编码器和一个解码器，即包含一个编码路径和一个解码路径，编码器由三维卷积层和三维最大值池化层组成。解码器由三维卷积层和三维上采样层组成。其中，三维卷积层的卷积核大小为3×3×3，激活函数为ReLU。跳跃连接层用于将编码器和解码器中维度相同的特征连接到一起，使得编码器和解码器中维度相同的特征在通道进行融合。应当注意，上述的深度学习网络的种类和布局用于说明而非限制，在此可以使用任何其他适当的配置。

进一步地，深度学习网络训练过程包括以下步骤：

(1)训练数据准备，获取多例dMRI数据，训练用的dMRI数据来自人类连接组计划(Human Connectome Project，HCP)中的健康成年人扩散磁共振数据，从中选取了20例dMRI数据用于深度学习网络的训练和验证；并对每一例dMRI数据进行数据处理，其中数据处理的过程如下：

(1-1)挑选出b值为预设数值的扩散加权图像以及对应的非扩散加权图像，并进行预处理；预处理的过程包括图像去噪、去吉布斯伪影以及偏场校正；预设数值可为1000s/mm

(1-2)对于预处理后的扩散加权图像以及对应的非扩散加权图像，采用加钱线性最小二乘法进行逐像素拟合得到DTI扩散张量场，以得到的DTI扩散张量场作为dMRI数据的真实扩散张量场

(1-3)对DTI扩散张量场应用DTI模型，使用预设的扩散梯度表进行dMRI数据重建，得到扩散加权图像，保证得到的扩散加权图像满足DTI模型，预设的扩散梯度表可以为给定15个扩散编码方向且b＝1000s/mm

/>

式中，s

(1-4)对于重建后的dMRI数据，加入莱斯噪声：

式中，I为重建后真实无噪的dMRI数据，N

(1-5)对于加入莱斯噪声的dMRI数据I

(2)深度学习网络训练；

(2-1)在多例dMRI数据中随机选择n例作为训练集，剩余部分作为验证集；如步骤(1-1)获取了20例dMRI数据中，则在20例dMRI数据中随机选择16例dMRI数据作为训练集，剩余4例dMRI数据作为验证集；

(2-2)采用训练集对深度学习网络进行训练，直到网络收敛，同时验证集损失达到稳定。

在上述步骤(1)和步骤(2)中，在进行深度网络训练时，以训练集中的dMRI数据在步骤(1)得到的球谐系数图以及该dMRI数据的真实扩散张量场

式中，

是由图像中每个体素对应的DTI扩散张量D∈R

进一步地，DTI量化参数图至少包括FA图和MD图；在步骤S5中，根据DTI扩散张量场，计算得到DTI量化参数图的步骤包括：

S51、对DTI扩散张量场的扩散张量D进行特征值分解，得到特征值λ

公式(4)如下：

在上述步骤S51中，DTI扩散张量场D

因此，对扩散张量D正定对称矩阵进行特征分解，得到特征值λ

为证明本发明扩散张量成像参数量化方法(以下简称本方法)的泛化能力，以下采用两个实施例进行说明。

实施例1

将本方法与现有的深度学习方法(dMRI数据作为深度学习网络的输入)以及传统去噪方法MPPCA(使用MPPCA进行去噪，然后使用WLLS对去噪后的dMRI数据进行逐像素张量拟合估计)进行了比较。本实施例对来自HCP的健康成年人扩散磁共振数据进行仿真。仿真数据包括一个b＝0s/mm

请参阅图4，图4为实施例1仿真数据的实验方法对比结果，展示了DTI扩散张量场的六个元素图、FA图和MD图以及对应的误差。如图4所示，本方法所得的DTI扩散张量场的六个元素图以及FA、MD图含有最少的噪声以及最接近真实值的细节，且RMSE值最低，证明了本方法在去噪和参数估计方面效果最佳，且泛化能力强。

实施例2，

将本方法与MPPCA算法(使用MPPCA进行去噪，然后使用WLLS对去噪后的dMRI数据进行逐像素张量拟合估计)进行了比较。现有的深度学习方法在15个扩散编码方向进行训练，不能直接应用于该真实数据(含有12个扩散编码方向)，因此未与现有深度学习方法作比较。本实施例采用的数据为真实数据，该数据为健康志愿者的大脑高分辨率dMRI数据(分辨率为0.9mm×0.9mm×0.9mm)，含有一个b＝0s/mm

请参阅图5，图5为实施例2高分辨率大脑真实数据的实验方法对比结果，展示了DTI扩散张量场的六个元素图、FA图和MD图以及局部放大图。如图5所示，使用本方法获得的DTI扩散张量场的六个元素图以及FA、MD图，其去噪效果和图像细节保留最好，证明了本方法在去噪和DTI参数量化方面效果好，且泛化能力强。

请参阅图6，图6为本发明扩散张量成像参数量化系统的模块示意图。本发明还提供了一种扩散张量成像参数量化系统，包括：

获取模块1，用于对扩散张量D正定对称矩阵进行特征分解，得到特征值λ

求取模块2，用于根据扩散编码方向g

处理模块3，用于对dMRI数据进行数据变换，得到扩散衰减图像，并根据扩散衰减图像和球谐函数的复共轭Y

输入模块4，用于将球谐系数图输入已构建的深度学习网络，得到DTI扩散张量场；

计算模块5，用于根据DTI扩散张量场，计算得到DTI量化参数图。

获取模块1获取的dMRI数据{S

求取模块2根据dMRI数据的扩散编码方向g

式中，P

处理模块3包括：

处理子模块，用于采用公式(1)对dMRI数据进行数据变换，公式(1)为：

式中，{(g

处理子模块分别使用扩散加权图像S

处理模块3还包括：

计算子模块，用于根据扩散衰减图像和球谐函数的复共轭Y

式中，

预先构建深度学习网络，并采用预设的训练数据对深度学习网络进行训练，深度学习网络输入的是球谐系数图，深度学习网络输出的是DTI扩散张量场。因此，输入模块4将计算得到的球谐系数图输入到预先构建的深度学习网络中，通过深度学习网络得到DTI扩散张量场。深度学习网络可以采用U-net网络，它是一个U型网络。

进一步地，DTI量化参数图至少包括FA图和MD图；计算模块5包括：

分解子模块，用于对DTI扩散张量场的扩散张量D进行特征值分解，得到特征值λ

公式(4)如下：

DTI扩散张量场D

因此，分解子模块对扩散张量D正定对称矩阵进行特征分解，得到特征值λ

请结合参阅图7，图7为本发明电子设备一实施例的结构示意框图。本发明一实施例还提出一种电子设备1001，包括存储器1003和处理器1002，存储器1003存储有计算机程序1004，处理器1002执行计算机程序1004时实现上述任一项扩散张量成像参数量化方法的步骤，包括：S1、获取dMRI数据{S

请结合参阅图8，图8为本发明计算机可读存储介质一实施例的结构示意框图。本发明一实施例还提供一种计算机可读存储介质2001，其上存储有计算机程序1004，计算机程序1004被处理器1002执行时实现上述任一项扩散张量成像参数量化方法的步骤，包括：S1、获取dMRI数据{S

相比于现有技术，本发明的有益效果为：使用深度学习网络和球谐函数对扩散张量成像参数进行量化，较传统的方法计算速度快，可以得到更准确的DTI量化参数；并且较其他基于深度学习网络的扩散张量参数量化方法，不受限于指定的采集方案，能够准确的估计出不同采集方案的dMRI数据的DTI扩散张量场，泛化能力强。

本领域普通技术人员可以理解实现上述实施例方法中的全部或部分流程，是可以通过计算机程序来指令相关的硬件来完成，所述的计算机程序可存储于一非易失性计算机可读取存储介质中，该计算机程序在执行时，可包括如上述各方法的实施例的流程。其中，本申请所提供的和实施例中所使用的对存储器、存储、数据库或其它介质的任何引用，均可包括非易失性和/或易失性存储器。非易失性存储器可以包括只读存储器(ROM)、可编程ROM(PROM)、电可编程ROM(EPROM)、电可擦除可编程ROM(EEPROM)或闪存。易失性存储器可包括随机存取存储器(RAM)或者外部高速缓冲存储器。作为说明而非局限，RAM以多种形式可得，诸如静态RAM(SRAM)、动态RAM(DRAM)、同步DRAM(SDRAM)、双速据率SDRAM(SSRSDRAM)、增强型SDRAM(ESDRAM)、同步链路(Synchlink)DRAM(SLDRAM)、存储器总线(Rambus)直接RAM(RDRAM)、直接存储器总线动态RAM(DRDRAM)、以及存储器总线动态RAM(RDRAM)等。

需要说明的是，在本文中，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的过程、装置、物品或者方法不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种过程、装置、物品或者方法所固有的要素。在没有更多限制的情况下，由语句“包括一个……”限定的要素，并不排除在包括该要素的过程、装置、物品或者方法中还存在另外的相同要素。

本发明并不局限于上述实施方式，如果对本发明的各种改动或变形不脱离本发明的精神和范围，倘若这些改动和变形属于本发明的权利要求和等同技术范围之内，则本发明也意图包含这些改动和变形。