活体心脏复杂组织学特征精准解析的超融合体素内不相干运动张量磁共振成像方法

1.本发明属于医学磁共振成像领域，具体涉及一种超融合体素内不相干运动张量磁共振成像方法。


背景技术：

2.心血管疾病作为人类健康的头号杀手，一直被医学界高度关注。研究表明，多数心血管疾病都同时伴有心肌纤维和心肌微循环系统结构及功能异常的现象，例如，冠脉粥样硬化导致的缺血性心肌病，由于心肌供氧和需氧失衡，致使心肌细胞变性、坏死，心肌纤维化瘢痕形成；同时，粥样硬化导致的冠脉狭窄也使得下级小动脉、小静脉和毛细血管的血液流速及血管分布发生改变。因此，能够获悉心肌纤维和心肌微循环系统的微观结构和生理状况，明确心肌纤维和微循环系统之间的关系，对进一步认识各种心血管疾病成因，以及心血管疾病的早期诊断和治疗具有重要的临床价值。
3.目前，获取心肌纤维和微循环系统组织学特征的手段主要包括离体解剖和医学影像。其中，前者通过在显微镜下观察心肌组织的切片样本，获取心肌纤维和微循环系统的结构信息。但由于被解剖的心肌组织已丧失活性，所以该方法不能反映心肌纤维和微循环系统的生理状况。相比而言，医学影像能够实现对活体心脏结构及功能的探测。例如，扩散张量磁共振成像能够获得活体心肌组织的纤维走向；单光子发射计算机断层成像，正电子发射断层成像以及动脉自旋标记磁共振成像能够检测心肌微循环的血液灌注率和血流分布情况。然而，遗憾的是，这些医学成像技术均无法同时获取心肌纤维和心肌微循环的组织学特征。这意味着为了全面、准确评估心脏的生理或病理状况，心血管疾病患者不得不接受多种成像设备的检查。这不仅增加了检查的复杂度和诊疗的相关费用，而且引入了不同类型图像间的配准问题以及数据的同一性问题。
4.体素内不相干运动(intra-voxel incoherent motion，ivim)成像是目前唯一不依赖造影剂便能在一次成像中同时获取组织细胞和微血管生理或病理信息的磁共振检查方法。ivim成像是扩散加权成像(diffusion weighted imaging，dwi)向微循环领域的拓展，通过建模细胞性水分子和血液性水分子的扩散和伪扩散行为间接描述细胞和微血管的组织学特征。细胞性水分子的扩散运动与细胞密度，细胞膜通透性，细胞外间隙以及扩散介质的黏性度有关，能够在一定程度上反映组织细胞活性和纤维化程度；血液性水分子的伪扩散运动与血管壁透性，血液流速以及微血管拓扑结构等有关，能够在一定程度上反映血管内皮细胞损伤及毛细血管扩张等血管和血液相关病变。
5.然而，常规的ivim成像方法将血液性水分子狭义地认定为在毛细血管内流动的水分子，并且假设组织细胞和毛细血管在单位体素内随机分布，即水分子运动对扩散加权信号衰减的贡献与扩散敏感梯度的方向无关，表现为各向同性。但就心肌组织而言，心肌细胞(或纤维)与毛细血管均服从各向异性分布。也就是说，水分子沿不同方向表现出不同的扩散或伪扩散运动行为。尽管最近提出的张量化ivim模型能够在一定程度上描述心肌细胞和
毛细血管各向异性的分布特征，但忽视了各向同性分布的小动脉和小静脉中水分子运动对扩散加权信号衰减的贡献。小动脉和小静脉作为心肌微循环系统的重要组成部分，其组织学异常与多种心血管疾病的发生、发展密切相关。此外，如果ivim建模时没有考虑小动脉和小静脉中水分子的伪扩散行为，那么这些水分子的运动将被错误地编码到反映心肌纤维和毛细血管组织学特征的参数中，进而导致严重的估计偏差。


技术实现要素：

6.本发明要解决的技术问题是：
7.本发明的目的是为了解决现有医学成像技术无法在一次检查中全面且准确地获得各向异性分布的心肌纤维、毛细血管以及各向同性分布的小动脉、小静脉组织学特征的问题，而提出了活体心脏复杂组织学特征精准解析的超融合体素内不相干运动张量磁共振成像方法。
8.本发明为解决上述技术问题所采用的技术方案为：
9.本发明活体心脏复杂组织学特征精准解析的超融合体素内不相干运动张量磁共振成像方法按下列步骤实现：
10.一、采用张量模型表达活体心脏单位体素内心肌纤维间隙和纤维中水分子的各向异性扩散运动以及毛细血管中水分子的各向异性伪扩散运动对扩散加权信号衰减的贡献，具体表达式如下：
[0011][0012]
式中：s
dc
(b)为给定b值下心肌纤维间隙和纤维中水分子以及毛细血管中水分子运动所产生的扩散加权信号；
[0013]
s0为未施加扩散敏感梯度(即b＝0)的磁共振信号；
[0014]
fd为单位体素内心肌纤维间隙及纤维中水分子容积张量(二阶矩阵)，fc为单位体素内毛细血管中水分子容积张量(二阶矩阵)，d为单位体素内纤维性水分子扩散张量(二阶矩阵)，d
c*
单位体素内毛细血管中水分子伪扩散张量(二阶矩阵)，张量fd，fc，d，d
c*
均为待求解的未知量；
[0015]
为扩散敏感梯度单位矢量(人为设定的已知量)；t为转置运算符；b为衰减因子(人为设定的已知量)，可通过如下表达式设定：
[0016][0017]
式中γ——氢核的旋磁比(常数)；
[0018]
δ——扩散敏感梯度持续时间(人为设定的已知量)；
[0019]
δ——两个扩散敏感梯度的时间间隔(人为设定的已知量)；
[0020]
二、使用一个连续的伪扩散变量d
′
和一个连续的容积分数函数f
′
(d
′
)描述单位体素内小动脉流和小静脉流在扩散敏感梯度作用下产生的扩散加权信号：
[0021][0022]
式中，s
av
(b)表示为给定b值下小动静脉内水分子运动所产生的扩散加权信号；
[0023]fav
为小动脉和小静脉血管内水分子容积分数(标量)，其中，f
′
(d
′
)/f
av
为小动静脉容积分数密度函数，记为ρ
av
(d
′
)；
[0024]
ρ
av
(d
′
)在某区间上的积分表示拥有某类伪扩散属性的小动静脉流占全部小动静脉容积和的百分比；ρ
av
(d
′
)在区间(0,+∞)上的积分等于1；代入ρ
av
(d
′
)，等式(3)变为：
[0025][0026]
令ρ
av
(d
′
)服从高斯函数，则
[0027][0028]
式中：为小动静脉流的平均伪扩散系数(标量)，σ
av
为小动静脉流的伪扩散系数方差(标量)，标量和σ
av
均为待求解的未知量；
[0029]
n为归一化常数，表示为：
[0030][0031]
三、将表达式(1)和(4)加权求和，便可得到最终形态的超融合ivim张量模型(即离散指数加权方程)，其表达式为：
[0032][0033]
式中：s(b)为s
dc
(b)与s
av
(b)的和；s(b)和s0(即s(b＝0))均通过测量获得；
[0034]
四、建立超融合ivim张量模型后，采用围绕非负最小二乘法(non-negative least square，nnls)的步进递归策略求解超融合ivim张量模型(7)中各未知参数：fd张量(二阶矩阵)、fc张量(二阶矩阵)、d张量(二阶矩阵)、dc*张量(二阶矩阵)、f
av
标量、标量、σ
av
标量，并利用这七个参数表征活体心脏单位体素内心肌的复杂组织学特征；
[0035]
五、通过对活体心脏的每个体素建立超融合ivim模型并求得其相应的七个参数，从而完成对整个活体心脏的体素内不相干运动张量磁共振成像。
[0036]
进一步地，在步骤四中，求解超融合ivim张量模型参数的具体过程为：
[0037]
首先，针对每个体素，测得人为设定的多个b值所对应的s(b)(包含s0)信号，并将这些信号按b值大小排序(从小到大)，起始b值记为0；
[0038]
然后，基于所有b值对应的s(b)信号执行nnls，并记录nnls输出的谱峰数量(记为p)；
[0039]
接下来，从小到大依次剔除b值(以下称为步进)，并基于剩余b值所对应s(b)信号执行nnls；每步进一次，nnls将输出一个p值，待p＝1时停止步进，并将此时最小的b值记为bd；基于bd～b
max
(b
max
代表最大b值)范围内的s(b)信号，利用非线性最小二乘曲线拟合方法，求解出纤维性水分子的扩散相关参数，即张量fd和d；
[0040]
再以bd为起点，向b＝0的方向回归(以下称为递归)；每次递归后，再次执行nnls，并记录其输出的p值；待p＝2时停止递归，将此时的最小b值记为bc；基于bc～b
max
范围内的s
(b)信号以及已求解产生的fd和d，利用非线性最小二乘曲线拟合方法计算出毛细血管内水分子的伪扩散相关参数，即张量fc和d
c*
；
[0041]
最后，基于全b值范围的s(b)信号以及已求解产生的fd、d、fc、d
c*
，利用非线性最小二乘曲线拟合，求解出小动脉和小静脉的伪扩散相关参数，即标量f
av
、和σ
av
。
[0042]
进一步地，步骤一中描述各向异性毛细血管的指数表达式为：
[0043][0044]
式中j——虚数因子；
[0045]
u0——毛细血管内平均血流速度；
[0046]
μ——毛细血管的俯仰角；
[0047]
v——毛细血管的方位角；
[0048]
δ——两个扩散敏感梯度的时间间隔；
[0049]
——倒数位移矢量，如式(9)所示；
[0050][0051]
式中δ——扩散敏感梯度宽度；
[0052]
——扩散敏感梯度矢量。
[0053]
进一步地，步骤二中令ρ
av
(d
′
)为狄拉克函数。
[0054]
进一步地，利用莱温伯格-马夸尔特(levenberg-marquardt)方法拟合非线性衰减的s(b)信号。
[0055]
活体心脏复杂组织学特征精准解析的超融合体素内不相干运动张量磁共振成像系统，该系统具有与上述步骤对应的程序模块，运行时执行上述的活体心脏复杂组织学特征精准解析的超融合体素内不相干运动张量磁共振成像方法中的步骤。
[0056]
本发明具有以下有益技术效果：
[0057]
本发明采用张量模型表达心肌纤维间隙和纤维内水分子的各向异性扩散运动以及毛细血管内水分子的各向异性伪扩散运动对扩散减权信号衰减的贡献；使用一个连续的伪扩散变量和一个连续的容积分数密度函数描述单位体素内小动静脉流在扩散敏感梯度作用下产生的扩散加权信号；将各表达式按容积分数加权求和，构建超融合ivim张量模型；采用围绕nnls的步进递归策略求解超融合ivim张量模型参数。
[0058]
本发明首先利用二阶张量和概率密度函数的表达潜能，将采自心肌组织的超维扩散加权信号中所蕴含的各类水分子运动特征分别建模为扩散张量指数式、伪扩散张量指数式和伪扩散系数概率密度函数积分指数式；然后，将三个表达式按各自的容积分数加权求和，形成能够同时反映各向异性心肌纤维和毛细血管以及各向同性小动脉和小静脉组织学特征的超融合ivim张量模型；最后，采用围绕nnls的步进递归策略求解超融合ivim张量模型参数。
[0059]
本发明提出用参数fd、fc、d、d
c*
、f
av
、σ
av
表征活体心脏单位体素的复杂心肌组织学特征，通过对活体心脏的每个体素建立超融合ivim模型并求得其相应的七个参数，从而完成对整个活体心脏的体素内不相干运动张量磁共振成像。采用围绕nnls的步进递归策略求解超融合ivim张量模型参数能实现对复杂数学模型的准确求解，大大提高求解精度。
使用本发明能够获取心肌细胞完整性，细胞膜透性，心肌纤维走向，以及微循环(毛细血管、小动脉和小静脉)网络拓扑结构、血液流速和容积分数等重要的心肌组织学特征。
[0060]
本发明所述活体心脏复杂组织学特征精准解析的超融合体素内不相干运动张量磁共振成像方法的有益技术效果表现在以下几个方面：
[0061]
(1)通过一次成像，可同时获得，心肌细胞完整性，细胞膜透性，心肌纤维走向，以及微循环(毛细血管、小动脉和小静脉)网络拓扑结构、血液流速和容积分数等重要的心肌组织学特征，避免了多种检查手段相结合产生的累加费用以及由此带来的图像配准难和数据同一性差等问题；
[0062]
(2)上述特征能够准确反映心肌细胞坏死、水肿、凋亡，微血管内皮细胞损伤、血管畸形以及局灶性缺血等心肌组织的病理学改变，为心血管疾病的临床诊断和治疗提供了重要的参考依据；
[0063]
(3)采用围绕nnls的步进递归策略解决了常规算法不能解析复杂ivim张量方程的问题，实现了对心肌组织数学模型的准确求解，杜绝了经验主义导致的误判；
[0064]
(4)在获取心肌微循环灌注特征方面，相比传统增强磁共振成像手段，该方法不依赖造影剂，能够有效避免造影剂注射可能引发的过敏风险，如恶心、呕吐、喉头水肿、肾功能损伤等，可显著减轻患者的身体痛苦和心理负担。
附图说明
[0065]
图1为活体心脏复杂组织学特征精准解析的超融合体素内不相干运动张量磁共振成像方法的流程图；
[0066]
图2为实施例一中计算生成的活体心脏纤维性水分子扩散张量d的三维可视图；
[0067]
图3为实施例一中计算生成的活体心脏心肌纤维间隙及纤维内水分子容积张量fd的三维可视图；
[0068]
图4为实施例一中计算生成的活体心脏毛细血管内水分子伪扩散张量d
c*
的三维可视图；
[0069]
图5为实施例一中计算生成的活体心脏毛细血管内水分子容积张量fc的三维可视图；
[0070]
图6为实施例一中计算生成的活体心脏小动脉和小静脉血管内水分子容积分数f
av
的伪彩图；
[0071]
图7为实施例一中计算生成的活体心脏小动脉和小静脉流平均伪扩散系数的伪彩图；
[0072]
图8为实施例一中计算生成的活体心脏小动脉和小静脉流伪扩散系数方差σ
av
的伪彩图；
具体实施方式
[0073]
具体实施方式一：本实施方式活体心脏复杂组织学特征精准解析的超融合体素内不相干运动张量磁共振成像方法按下列步骤实现：
[0074]
一、采用张量模型表达活体心脏单位体素内心肌纤维间隙和纤维中水分子的各向异性扩散运动以及毛细血管中水分子的各向异性伪扩散运动对扩散加权信号衰减的贡献，
具体表达式如下：
[0075][0076]
式中：s
dc
(b)为给定b值下心肌纤维间隙和纤维中水分子以及毛细血管中水分子运动所产生的扩散加权信号；
[0077]
s0为未施加扩散敏感梯度(即b＝0)的磁共振信号；
[0078]
fd为单位体素内心肌纤维间隙及纤维中水分子容积张量(二阶矩阵)，fc为单位体素内毛细血管中水分子容积张量(二阶矩阵)，d为单位体素内纤维性水分子扩散张量(二阶矩阵)，d
c*
单位体素内毛细血管中水分子伪扩散张量(二阶矩阵)，张量fd，fc，d，d
c*
均为待求解的未知量。
[0079]
为扩散敏感梯度单位矢量(人为设定的已知量)；t为转置运算符；b为衰减因子(人为设定的已知量)，可通过如下表达式设定：
[0080][0081]
式中γ——氢核的旋磁比(常数)；
[0082]
δ——扩散敏感梯度持续时间(人为设定的已知量)；
[0083]
δ——两个扩散敏感梯度的时间间隔(人为设定的已知量)；
[0084]
二、使用一个连续的伪扩散变量d
′
和一个连续的容积分数函数f
′
(d
′
)描述单位体素内小动脉流和小静脉流在扩散敏感梯度作用下产生的扩散加权信号：
[0085][0086]
式中，s
av
(b)表示为给定b值下小动静脉内水分子运动所产生的扩散加权信号；
[0087]fav
为小动脉和小静脉血管内水分子容积分数(标量)，其中，f
′
(d
′
)/f
av
为小动静脉容积分数密度函数，记为ρ
av
(d
′
)；
[0088]
ρ
av
(d
′
)在某区间上的积分表示拥有某类伪扩散属性的小动静脉流占全部小动静脉容积和的百分比；ρ
av
(d
′
)在区间(0,+∞)上的积分等于1；代入ρ
av
(d
′
)，等式(10)变为：
[0089][0090]
令ρ
av
(d
′
)服从高斯函数，则
[0091][0092]
式中：为小动静脉流的平均伪扩散系数(标量)，σ
av
为小动静脉流的伪扩散系数方差(标量)，标量和σ
av
均为待求解的未知量；
[0093]
n为归一化常数，表示为：
[0094][0095]
三、将表达式(1)和(4)加权求和，便可得到最终形态的超融合ivim张量模型(即离
散指数加权方程)，其表达式为：
[0096][0097]
式中：s(b)为s
dc
(b)与s
av
(b)的和。s(b)和s0(即s(b＝0))均通过测量获得。
[0098]
四、建立超融合ivim张量模型后，采用围绕非负最小二乘法(non-negative least square，nnls)的步进递归策略求解超融合ivim张量模型(7)中各未知参数：fd张量(二阶矩阵)、fc张量(二阶矩阵)、d张量(二阶矩阵)、dc*张量(二阶矩阵)、f
av
标量、标量、σ
av
标量。
[0099]
本发明采用围绕非负最小二乘法(non-negative least square，nnls)的步进递归策略求解超融合ivim张量模型参数是本发明重要发明点。首先，针对每个体素，测得人为设定的多个b值所对应的s(b)(包含s0)信号，并将这些信号按b值大小排序(从小到大)，起始b值记为0；然后，基于所有b值对应的s(b)信号执行nnls，并记录nnls输出的谱峰数量(记为p)；接下来，从小到大依次剔除b值(以下称为步进)，并基于剩余b值所对应s(b)信号执行nnls；每步进一次，nnls将输出一个p值，待p＝1时停止步进，并将此时最小的b值记为bd；基于bd～b
max
(b
max
代表最大b值)范围内的s(b)信号，利用非线性最小二乘曲线拟合方法，求解出纤维性水分子的扩散相关参数，即张量fd和d；再以bd为起点，向b＝0的方向回归(以下称为递归)；每次递归后，再次执行nnls，并记录其输出的p值；待p＝2时停止递归，将此时的最小b值记为bc；基于bc～b
max
范围内的s(b)信号以及已求解产生的张量fd和d，利用非线性最小二乘曲线拟合方法计算出毛细血管内水分子的伪扩散相关参数，即张量fc和d
c*
；最后，基于全b值范围的s(b)信号以及已求解产生的张量fd、d、fc、d
c*
，利用非线性最小二乘曲线拟合，求解出小动脉和小静脉的伪扩散相关参数，即标量f
av
、和σ
av
。
[0100]
五、采用围绕nnls的步进递归策略求解出超融合ivim张量模型参数后，利用这些参数表征活体心脏单位体素内心肌的复杂组织学特征。随后，通过对活体心脏的每个体素建立超融合ivim模型并求得其相应的模型参数，从而完成对整个活体心脏的体素内不相干运动张量磁共振成像。
[0101]
本实施方式步骤一中张量(二阶矩阵)及其衍生参数，即局部各向异性指数(fractional anisotropy，fa)、相对各向异性指数(relative anisotropy，ra)、平均扩散率(mean diffusivity，md)以及westin系数，能够有效反映各向异性组织的结构和功能特征。
[0102]
本实施方式步骤二：本实施方式与具体实施方式一不同的是步骤一中描述各向异性毛细血管的指数表达式为：
[0103][0104]
式中j——虚数因子；
[0105]
u0——毛细血管内平均血流速度；
[0106]
μ——毛细血管的俯仰角；
[0107]
v——毛细血管的方位角；
[0108]
δ——两个扩散敏感梯度的时间间隔；
[0109]
——倒数位移矢量，如式(9)所示。
[0110][0111]
式中δ——扩散敏感梯度宽度；
[0112]
——扩散敏感梯度矢量；
[0113]
具体实施方式三：本实施方式与具体实施方式一或二不同的是步骤二中令ρ
av
(d
′
)为狄拉克函数。
[0114]
具体实施方式四：本实施方式与具体实施方式一至三之一不同的是步骤四利用莱温伯格-马夸尔特(levenberg-marquardt)算法拟合非线性衰减的扩散加权信号。
[0115]
实施例一：采集活体人类心脏多方向、多b值扩散加权图像。为了尽可能降低心脏和呼吸运动对扩散加权信号和图像质量的影响，采用前瞻性心电门控采集技术，并将触发延迟设置在心脏舒张末期；同时，在数据采集期间，全程开启实时片跟踪和重定位模式，减少呼吸运动伪影。具体成像参数如下：tr 2beats，te 66ms，flip angle 81
±1°
，slice 5，slice thickness 10mm，fov 260
×
200mm，matrix size 224
×
224，number of excitations 1，sense 2.5，b-values(0，5，10，25，50，75，100，200，300，400s/mm2)，number of diffusion directions 6。针对每个扩散敏感梯度方向执行围绕nnls的步进递归策略。基于在所有梯度方向上解析出的扩散和伪扩散相关参数，利用等式(10)重构fd、f
p
、d和dc*矩阵变量以及和σ
av
参数。
[0116][0117]
式中(g
xi
,g
yi
,g
zi
)——第i个扩散敏感梯度方向；
[0118]
y——扩散或伪扩散二阶张量元素；
[0119]
xi——第i个扩散敏感梯度方向上求解出的扩散或伪扩散相关参数；
[0120]
n——扩散敏感梯度方向的数量。
[0121]
如图2-图8所示，分别为计算生成的活体心脏纤维性水分子扩散张量d的三维可视图、活体心脏心肌纤维间隙及纤维内水分子容积张量fd的三维可视图、活体心脏毛细血管内水分子伪扩散张量d
c*
的三维可视图、活体心脏毛细血管内水分子容积张量fc的三维可视图、活体心脏小动脉和小静脉血管内水分子容积分数f
av
的伪彩图、活体心脏小动脉和小静脉流平均伪扩散系数的伪彩图、活体心脏小动脉和小静脉流伪扩散系数方差σ
av
的伪彩图，从这些图中可看出心肌细胞完整性，细胞膜透性，心肌纤维走向，以及微循环(毛细血管、小动脉和小静脉)网络拓扑结构、血液流速和容积分数等重要的心肌组织学特征。