一种基于流场分布的神经纤维续跟踪方法

摘要

一种基于流场分布的神经纤维续跟踪方法，针对现有物理设备成像分辨率有限、有限的分辨率获得的mask存在的锯齿形边界导致纤维跟踪过早终止等不足，将纤维束看作流束，建立纤维成像中体素的流场分布模型，由建立的流场分布模型扩展流线型跟踪算法，基于流场分布模型计算纤维跟踪过程中的方向信息，通过该方向信息，计算纤维停止点的连续性，判断在当前停止点的基础上是否需要继续跟踪。本发明计算纤维跟踪停止点的连续性判断是否进行纤维续跟踪，有效避免了纤维跟踪过早终止这一问题，提升了纤维跟踪结果的准确性与纤维分布的全面性。

说明书

技术领域

本发明涉及医学图像处理，尤其是一种神经纤维续跟踪方法。

背景技术

基于磁共振成像技术的扩散磁共振成像是目前唯一的一种非侵入式的神经纤维成像方法，在实际中得到了广泛应用。然而，目前的磁共振仪器成像的体素大小在毫米级，人脑白质神经纤维的直径却是微米级，并且因为硬件仪器的限制，成像体素的空间分辨率无法无限制地增加，导致成像的分辨率有限。而在纤维跟踪的过程中，通常需要通过计算成像体素信息，得到判断纤维跟踪终止条件之一的“掩膜(mask)”。由有限的分辨率获得的mask存在的锯齿形边界，会造成纤维在跟踪中过早地到达边界而终止跟踪，无法到达实际的区域。

针对这个问题，一方面，有研究者提出解剖学约束的纤维束成像算法，对磁共振过程中的数据进行匹配、分割等操作，从而得到更加准确的mask，同时加入解剖学的相关先验知识作为约束去剔除不符合解剖学经验的纤维。但这个算法需要和其他的跟踪算法结合使用，且其主要是对最终的跟踪结果的筛选，不能使过早停止跟踪的纤维继续跟踪到实际区域；另一方面，随着近些年来机器学习发展迅速，有学者提出基于机器学习的纤维跟踪算法，对样本数据进行训练，以给定的100个连续的采样方向下各自的概率值作为纤维微结构的方向分布，从而进行纤维跟踪。然而该方法的训练样本来自其他跟踪算法，只是基于现有的跟踪算法结果的后处理，纤维跟踪的结果受训练标本影响极大。

发明内容

为了克服现有物理设备成像分辨率有限、有限的分辨率获得的mask存在的锯齿形边界导致纤维跟踪过早终止的不足，本发明提出一种基于流场分布的纤维续跟踪方法，有效避免了纤维跟踪过早终止这一问题，提升了纤维跟踪结果的准确性与纤维分布的全面性。

本发明解决其技术问题所采用的技术方案是：

一种基于流场分布的纤维续跟踪方法，所述方法包括以下步骤：

步骤一：建立纤维成像中体素的流场分布模型，过程如下：

将纤维束看作流束，空间内点的流线切线即流场，计算流场系数矩阵A得到体素的流场分布；

步骤二：扩展流线型跟踪算法，过程如下：

由建立的流场分布模型扩展流线型跟踪算法，基于流场分布模型计算纤维跟踪过程中的方向信息。

步骤三：计算纤维停止点的连续性，判断在当前停止点的基础上是否需要继续跟踪，过程如下：

基于流场分布模型的方向信息，在纤维跟踪的停止点计算当前停止点的连续性，判断是否继续进行跟踪。

进一步，所述步骤一中，把纤维束作为一种流束考虑，流束可以认为是一系列流线的集合S＝{s

通过一组三元二次多项式来近似表示纤维体素的流场分布：

υ(x,y,z)＝AX(x,y,z) (2)

其中，系数矩阵A：

X(x,y,z)＝[x

在同一束纤维束中，假设水分子的扩散产生的位移具有连续性，根据连续不可压缩流体理论，通过流体发散概念来描述纤维的空间连续性：

当某一神经纤维束的起始区域或终止区域不在当前流场时，认为：

divΩ＝0 (5)

考虑到使流场模型最接近纤维的方向分布FOD，将系数矩阵A优化为如下代价函数计算：

其中，Φ(v(x,y,z))是在点(x,y,z)时FOD的概率分布，取其中心体素的26邻域体素中FOD的峰值P＝[p

计算得到系数A后，体素内的任意一点的纤维方向信息用AX表示。

再进一步，所述步骤二中，扩展流线型跟踪算法：

其中s(t)表示纤维弧线，t表示体素内轨迹上的某一点，v(t)表示该点处的切向量；

把流场分布AX作为纤维的局部方向信息，将公式(8)改为：

其中t

更进一步，所述步骤三中，在某一体素P2中，因为缺失方向信息或跟踪超出mask区域，纤维跟踪在点t停止；选择点t的纤维走向来作为点t所在体素即P2的方向信息，通过公式(7)计算该体素的流场分布模型。如果该处流场不存在，那么纤维跟踪在点t终止，设点t为终止点；如果该处流场存在，根据该点的流场分布，继续通过公式(9)进行纤维跟踪。纤维继续跟踪到下一个体素P3，如果体素P3为空即不含流场，则纤维不执行跟踪，纤维跟踪在点t终止，设点t为终止点；若体素P3不为空，即存在流场，计算跟踪得到的纤维区间段(P1、P2、P3)纤维分布的连续性，用该连续性表示点t的连续性；

通过点t所在体素P2为中心体素，计算区间段纤维分布的连续性：

其中，A为点t所在体素P2的流场系数，X为A对应坐标系，A

如果纤维分布的连续性高，即Con值小，可以认为在该点t满足纤维续跟踪条件，根据公式(9)继续作跟踪迭代；如果纤维分布的连续性低，即Con值大，可以认为在点t不满足纤维续跟踪的条件，则纤维跟踪在点t终止，设点t为终止点。

本发明的有益效果表现在：有效避免了纤维跟踪过早终止这一问题，提升了纤维跟踪结果的准确性与纤维分布的全面性。

附图说明

图1是本方法判断是否在停止点做续跟踪的示意图。

具体实施方案

下面对本发明作进一步说明。

参照图1，一种基于流场分布的神经纤维续跟踪方法，包括以下步骤：

步骤一：建立纤维流场分布模型，过程如下：

把纤维束作为一种流束考虑，流束可以认为是一系列流线的集合S＝{s

可以通过一组三元二次多项式来近似表示纤维体素的流场分布：

υ(x,y,z)＝AX(x,y,z) (2)

其中，系数矩阵A：

X(x,y,z)＝[x

在同一束纤维束中，假设水分子的扩散产生的位移具有连续性，根据连续不可压缩流体理论，通过流体发散概念来描述纤维的空间连续性：

当某一神经纤维束的起始区域或终止区域不在当前流场时，可以认为：

divΩ＝0 (5)

考虑到使流场模型最接近纤维的方向分布(fiber orientation distribution,FOD)，将系数矩阵A优化为如下代价函数计算：

其中，Φ(v(x,y,z))是在点(x,y,z)时FOD的概率分布，取其中心体素的26邻域体素中FOD的峰值P＝[p

计算得到系数A后，体素内的任意一点的纤维方向信息可以用AX表示；

步骤二：扩展流线型跟踪算法，过程如下：

扩展流线型跟踪算法：

其中s(t)表示纤维弧线，t表示体素内轨迹上的某一点，v(t)表示该点处的切向量；把流场分布AX作为纤维的局部方向信息，将公式(8)改为：

其中t

步骤三：计算纤维停止点的连续性，判断在当前停止点的基础上是否需要继续跟踪，过程如下：

如图1左图，在某一体素P2中，因为缺失方向信息或跟踪超出mask区域(图1右图白色区域即mask区域)，纤维跟踪在点t停止。选择点t的纤维走向来作为点t所在体素即P2的方向信息，通过公式(7)计算该体素的流场分布模型。如果该处流场不存在，那么纤维跟踪在点t终止，设点t为终止点；如果该处流场存在，根据该点的流场分布，继续通过公式(9)进行纤维跟踪。纤维继续跟踪到下一个体素P3，如果体素P3为空即不含流场，则纤维不执行跟踪，纤维跟踪在点t终止，设点t为终止点；若体素P3不为空，即存在流场，计算跟踪得到的纤维区间段(P1、P2、P3)纤维分布的连续性，用该连续性表示点t的连续性；

通过点t所在体素P2为中心体素，计算区间段纤维分布的连续性：

其中，A为点t所在体素P2的流场系数，X为A对应坐标系，A

如果纤维分布的连续性高，即Con值小，可以认为在该点t满足纤维续跟踪条件，根据公式(9)继续作跟踪迭代；如果纤维分布的连续性低，即Con值大，可以认为在点t不满足纤维续跟踪的条件，则纤维跟踪在点t终止，设点t为终止点。

本说明书的实施例所述的内容仅仅是对发明构思的实现形式的列举，仅作说明用途。本发明的保护范围不应当被视为仅限于本实施例所陈述的具体形式，本发明的保护范围也及于本领域的普通技术人员根据本发明构思所能想到的等同技术手段。