利用交互式患者特异性模型为深部脑刺激疗法定位神经纤维束

摘要

本项目的目标是开发一个交互式、针对患者个体的建模流程，以近实时地模拟深部脑刺激的效果，并提供有意义的反馈，说明这些设备如何影响大脑中的神经活动。

摘要

脑深部电刺激（DBS）通过植入电极对特定脑区进行刺激，是治疗运动障碍的成熟疗法，并正逐步应用于越来越多的其他疾病。计算建模已被成功用于预测DBS的临床效应；然而，随着DBS设备复杂性的不断增加，亟需开发新的建模技术。这些模型还需能够快速而准确地生成预测结果。本项目的目标是建立一个图像处理流程，将结构磁共振成像（MRI）和弥散加权成像（DWI）整合到一个可交互的、针对个体患者的模型中，以模拟DBS的作用效果。可在患者模型内放置虚拟DBS电极导线，并设定激活触点和刺激参数；当导线位置或方向发生变化时，系统可在近实时条件下（约10秒内）生成新的有限元网格并求解生物电场问题。该系统还支持模拟多个邻近导线，通过调节不同导线上阳极和阴极的配置实现电流导向。本文所介绍的技术降低了构建和使用计算模型的难度，同时为非建模专业背景的研究人员或临床医生提供了关于电极位置、电极设计及刺激配置影响的有意义反馈。

引言

脑深部电刺激（DBS）是一种已确立的治疗运动障碍的方法，适用于原发性震颤¹和帕金森病²。该疗法目前也正在被研究用于治疗越来越多的其他疾病，包括创伤性脑损伤³、抽动秽语综合征⁴和抑郁症⁵。DBS 系统需要通过手术植入电极导线，以在特定脑区施加电刺激，从而调节正在进行的神经活动⁶。电极的位置和刺激参数都会影响提供治疗效果的神经回路调节。电极位置的微小差异可能会影响治疗窗口，甚至在达到治疗效果之前增加不良副作用发生的风险^7,8,9。在实际临床中，通常难以预测刺激对神经活动的影响；因此，这一治疗窗口需由临床医生根据每位患者的具体情况，在编程调试刺激设备的过程中逐一确定^8,9。随着新一代 DBS 设备的不断出现，这一过程正变得愈加复杂。例如，新型电极导线设计正在引入更多电触点^10,11,12，在某些情况下，多个电极导线被植入彼此邻近的区域¹³。因此，有必要能够探索并预测 DBS 在日益庞大和复杂的参数空间中的作用效果。

计算建模与分析可用于在个体化基础上预测深部脑刺激（DBS）对生理和临床效果的影响。这些模型采用有限元建模（FEM）方法，构建脑组织以及植入电极生物物理特性的三维表征。FEM生物电场模型已成功用于预测DBS的效果¹⁴，但迄今为止，这些模型的生成过程耗时较长且计算成本高昂。随着DBS设备复杂性的不断增加，亟需开发新型建模技术。这些个体化患者模型应能够在改变电极位置或刺激参数时，提供接近实时的视觉反馈。用户可在数秒内获得关于电极位置和刺激设置的反馈信息，从而在几分钟内持续优化电极的放置。个体化特性的实现依赖于在构建FEM时整合患者的解剖结构、大脑的形态与尺寸，并应用其脑组织的生物物理特性，例如各向异性组织电导率。各向异性电导率描述了电流在不同脑区中的扩散方式，可通过非侵入性方法对全脑进行测量，其过程类似于常规磁共振成像（MRI）。

不使用患者特异性信息的脑深部电刺激（DBS）建模方法可以快速预测刺激效应，但由于采用了通用的脑组织几何结构和电导率值，其预测准确性较低。在此类方法中，所有患者均使用同一个有限元模型（FEM），且可预先计算预测的神经活动。而患者特异性模型由于需为每位个体单独构建新的有限元模型，因此无法通用化或预先计算。这类模型构建更为耗时，但准确性更高。目前有多个因素限制了这些模型的构建和使用速度：1）在建模流程初期修改参数（如电极位置）需要人工手动更新后续所有步骤；2）建模流程中的各个步骤之间难以有效集成，需在多个软件平台之间传递数据。通常情况下，我们希望评估多种不同情形，例如电极位置、刺激参数或电极设计。为了对这些变化对患者治疗效果的影响提供有意义的反馈，相关结果必须既准确又能够快速生成。

我们的目标是介绍一种构建患者特异性模型的新技术，该技术充分利用通用模型所获得的速度优势，并自动化处理流程中的多个步骤，从而创建一个交互式建模环境，能够提供关于深部脑刺激（DBS）效果的近实时可视化反馈。交互式仿真允许用户快速测试预测结果并获得输出，而无需专注于模型构建的细节。当需要探索的参数空间较大且这些参数对仿真结果的影响尚不明确时，这种方法尤为有益。我们将描述从磁共振成像（MRI）数据生成交互式、患者特异性有限元法（FEM）模型的处理流程中的各个步骤。利用本文所述的工具与技术，可显著降低构建FEM生物电场模型的时间成本，并为非建模专业背景的研究人员和临床医生提供一种便捷途径来使用这些模型。

本方案描述了如何从获取的磁共振成像（MRI）数据构建患者特异性的有限元模型，并模拟深部脑刺激（DBS）电极诱导的电场。构建此类模型的主要步骤包括：1）建立代表患者脑组织及植入电极的有限元模型（FEM）；2）向该有限元模型中添加脑组织的生物物理特性以及DBS设备的刺激参数；3）求解电极在模型中产生的电压分布。为构建用于DBS仿真的患者特异性模型，需要两种成像模态：T1加权MRI用于构建脑表面、脑室及特定核团的分割图像；弥散加权成像（DWI）——一种测量水分子弥散性的技术——用于估算全脑组织内的弥散张量¹⁵。将弥散张量转换为电导率张量，以逐体素方式量化组织的非均质性和各向异性生物物理特性¹⁶。通过求解泊松方程计算电极在全脑诱导的电压分布，该方程经有限元法（FEM）处理后简化为线性方程组 Ax = b，其中 A 为表示网格电导率和几何结构的刚度矩阵，x 为网格各节点处的电压解，b 则根据边界条件和电流源进行调整。

方案

1. 图像处理

1. T1 MRI 分割
   1. 下载并安装 FreeSurfer¹⁷。假设已获取 T1 MRI 图像，并且其格式为 DICOM 或 NIFTI。
   2. 输入以下命令以在 FreeSurfer 中创建患者目录并添加其 T1 图像：recon-all -s patientName -i /Full/Path/To/nii
   3. 输入以下命令以运行 FreeSurfer 的自动分割：recon-all -s patientName -all
   4. 使用 mri_convert 将 "aseg.auto.mgz" 从 FreeSurfer 专有的 .mgz 格式中转换出来。
2. 弥散加权成像（DWI）
   1. 下载并安装 Slicer¹⁸。假设已进行弥散加权成像，并以一系列 DICOM 文件形式提供。
      注意：一个优秀的命令行工具用于执行此重建是 FMRIB 软件库（FSL）¹⁹。本实验方案使用 FSL 在重建前完成平面回波、运动及涡流畸变校正。本分析中使用的 DWI 数据包含 41 个不同弥散敏感梯度方向的 DICOM 文件序列。
   2. 在 Slicer 中使用 DWIConverter 模块将 41 次独立采集重建为单个体积数据。指定 DWI 数据所在的 Input Dicom Data Directory。选择 DicomToNrrd 单选按钮，并提供输出 DWI 体积的名称。
   3. 使用 Diffusion Weighted Volume Masking 模块创建张量掩膜，以消除 DTI 估计中的背景噪声。将阈值参数保持默认值 0.5。
   4. 使用上一步创建的掩膜作为 Diffusion Tensor Mask，通过 DWI to DTI Estimation 将该单个 DWI 体积转换为弥散张量（DTI）体积。将估计参数设置为最小二乘法，并确保勾选 Shift Negative Eigenvalues 选项。
   5. 将新生成的 DTI 体积以 NRRD 格式保存。

2. 有限元模型生成

1. SCIRun
   1. 下载并安装 SCIRun 5 版本（http://www.sci.utah.edu/software/scirun.html）。SCIRun 是一个用于建模、仿真和可视化科学问题的通用问题求解环境。本项目中使用该计算工作台来创建、求解并与所开发的模型进行交互。
      注意：SCIRun 网络通过连接执行特定任务的独立模块序列构建而成。将数据传递通过由简单模块构成的网络，可实现更复杂的处理与仿真。
2. 构建电极几何结构
   注意：本项目中建模的电极为 Medtronic 3387 DBS 引线²⁰。它是一种圆柱形引线，包含四个触点，每个触点高度为 1.5 mm，直径为 1.27 mm，触点之间相距 1.5 mm。四个触点之间由非导电材料隔开。
   1. 根据所列尺寸，为 DBS 引线中的每个组件创建封闭的表面网格。通过使用两种基本几何体——圆柱体和球体——来构建完整的引线。
      注意：引线几何结构可在多种 3D 建模软件中创建。
   2. 创建一个厚度为 0.5 mm 的包膜层表面网格，包裹整个引线。
3. 创建全脑网格
   1. 加载第 1 部分中创建的大脑表面，作为有限元法（FEM）的外边界。
   2. 创建两个围绕 DBS 引线的同心盒状表面，以控制电极周围的网格密度。
      注意：DBS 引线和盒状表面应沿正 z 轴对齐，电极杆的底端位于坐标原点。这对于后续电极的旋转和平移操作至关重要。
   3. 创建一个点云，在网格的每个不同区域中放置一个点。网格的区域包括：杆段、四个触点、盒状表面以及大脑表面。附加在这些点上的数值代表该区域内每个四面体单元的最大体积约束。通过调整这些约束来控制电极内部及周围的网格密度。
   4. 使用 "InterfaceWithTetGen" 模块创建全脑网格。设置输入端口如下：
      端口 1：大脑表面（来自 1.1）
      端口 2：区域属性点云，体积约束（来自 2.3.3）
      端口 3：无
      ​端口 4：电极模型表面（来自 2.2）
4. 交互式电极放置
   注意：此步骤的目标是能够在大脑内部任意位置轻松移动电极，并使系统自动更新生成新的网格和仿真结果。
   1. 使用两个级联的 CreateGeometricTransform 模块，在 X 和 Y 平面内旋转电极。点击 Rotate，将第一个模块的 Axis Y 设置为 1.00，第二个模块的 Axis X 设置为 1.00。旋转角度通过底部滑块 Rotate Theta (degrees) 设置。
   2. 为每一个需要变换的变换操作和对象连接一个 TransformMeshWithTransform 模块。
   3. 创建一个 GenerateSinglePointProbeFromField 模块，并将其几何端口连接到 ViewScene 模块。在 ViewScene 窗口中按住 Shift 键并单击左键来移动球体。释放鼠标后，点探针模块将更新球体的位置。
      注意：该位置用于确定施加于电极几何结构的平移变换。
   4. 使用一个 TransformMeshWithTransform 模块，根据点探针的位置对电极几何结构的每一部分进行平移。
   5. 将两个旋转变换以及来自点探针控件的平移应用于除大脑表面以外的所有对象。
      注意：这三个变换应在每个对象送入 InterfaceWithTetGen 模块之前依次完成，以便网格生成算法能够知晓需在大脑内部的新位置和方向上对电极几何结构进行网格划分。每次在大脑内移动电极时，都将重建一个新的网格。

3. 生物电场计算

1. 设置电导率参数
   1. 加载第1步中创建的DTI体积数据，并使用MapFieldDataFromSourceToDestination模块的线性（加权）选项，将数据映射到第2步中通过InterfaceWithTetGen生成的四面体网格上。
      注意：此步骤中可实现任何各向异性电导率估计。
   2. 使用CalculateFieldData将电极杆区域的电导率设为1e-6，接触区域的电导率设为1e6。设置完所有电导率值后，将场数据输入至BuildFEMatrix模块。
2. 设置电流/电压源及边界条件
   1. 将第2步网格化过程中使用的脑表面数据输入至SetFieldDataToConstantValue模块，并将其值设为0，从而在整个脑表面创建一个电压汇。然后使用InsertVoltageSource模块将电压值施加到由InterfaceWithTetGen输出的四面体网格上。
   2. 对于单极刺激，通过在其中一个接触面使用ReportFieldInfo模块计算FieldCenter值，在该接触位置中心创建一个单一点。将该点的场数据设为-1伏特，然后使用InsertVoltageSource模块将该点源插入四面体网格中。
      注意：此时存在两个四面体网格副本，一个在脑表面设置了电压值，另一个在接触区域内包含一个点源。
   3. 合并两个InsertVoltageSource模块的输出结果，并将其与BuildFEMatrix的输出一同输入至ApplyFEMVoltageSource模块。
3. 求解正向问题
   注意：ApplyFEMVoltageSource模块的输出在端口1提供刚度矩阵，在端口2提供右端项矩阵。这两个组分是求解方程组所必需的。
   1. 将ApplyFEMVoltageSource的两个输出连接至SolveLinearSystem模块，选择共轭梯度法作为求解方法，Jacobi作为预条件子。该模块的输出为解矩阵，即四面体网格中每个节点处的电压值。
   2. 使用SetFieldData模块将电压解矩阵赋值到四面体网格，以实现电压分布的三维可视化。
4. 等值面可视化
   1. 将SetFieldData模块的结果输入至ExtractSimpleIsosurface模块。打开用户界面中的Quantity选项卡，将均匀分布的等值数设为10。
   2. 使用ShowField模块并将其连接至ViewScene模块，以透明方式显示等值面，并采用彩虹色图对电压解的最小值和最大值进行缩放显示。

结果

在本方案结束时，将提供构建患者特异性模型所需的全部输入信息：脑表面、电极几何结构以及电导率张量。同时还应已创建一个 SCIRun 网络，用于整合所有输入数据，构建有限元模型并模拟诱导产生的生物电场。该仿真模型支持电极在脑内的位置调整，并可修改诸如激活触点或刺激幅度等参数。

图1 展示了在多个电极位置下模拟和比较纤维激活预测的能力。本方案中未展示基于扩散成像生成纤维束的过程，但该分析可通过任何确定性示踪方法完成。与位置1相比，位置2在所有触点上均能在较低刺激强度下实现对目标纤维束更有效的激活。该分析在脑深部电刺激（DBS）手术规划阶段用于确定能够有效刺激目标纤维通路的电极植入路径。

该系统的一个新特性是能够快速切换模型中的电极几何结构，并能够模拟多个彼此靠近的电极导线。这两种方法均可用于更精确地控制电极周围电场的形态与导向。图2展示了在研究电极导线位置以刺激中央丘脑纤维通路并避开感觉核等邻近区域时，轴对称型与定向型脑深部刺激（DBS）导线在纤维激活预测方面的比较。对于相同的电极路径，可以看出定向型导线能够将纤维激活引导至目标区域，同时避免刺激非目标纤维。该分析可用于判断在何种情况下需要更换患者所接受的DBS治疗装置，以避免刺激邻近区域引发不必要的副作用。该方法在靶向刺激治疗帕金森病症状的丘脑底核时尤为有用，同时可避免刺激邻近的内囊。

图 1：电极位置与目标纤维激活预测的比较。 针对两种不同电极位置，预测在不同刺激强度（-0.5 V 至 -5.0 V）下四个接触点对总纤维束的激活情况。请点击此处查看该图的放大版本。

图2：轴对称与定向深部脑刺激（DBS）电极对丘脑中央部纤维束激活预测的比较。 (A) Medtronic 3387 圆柱形 DBS 电极与 Sapiens 定向 DBS 电极的可视化图示，显示在单一单极刺激幅度下，目标丘脑中央部纤维束的激活函数计算结果。(B)(i) 穿过 DBS 电极与纤维束的横截面切片。(ii) 横截面区域的二维投影，标示出 DBS 电极、目标纤维以及需避免激活的纤维。(C) 随着刺激幅度增加，轴对称与定向 DBS 电极在目标纤维及需规避纤维上的激活扩展情况。特定刺激幅度下被激活的纤维以红色表示，未被激活的纤维以蓝色表示。(D) 汇总面板 C 中所有模拟结果，显示纤维束各区域的预测刺激阈值幅度。请点击此处查看该图的放大版本。

讨论

本方案展示了一系列技术，可显著降低构建用于脑深部电刺激（DBS）的患者特异性模型所需的时间和计算成本，从而实现仿真结果的近实时反馈。快速反馈使得能够探索较大的参数空间，以更深入地理解参数变化如何影响模型结果。这些参数包括电极位置、激活触点的选择，以及在这些触点上的刺激波形的幅度、脉宽和频率。所提出工具的关键特性包括：1）提供一个简洁的用户界面，可调节模型参数，并对这些参数如何影响仿真结果进行近实时可视化；2）通过少量输入自动完成模型构建：包括脑表面、脑组织电导率张量，以及电极几何结构的表面表示。该自动化过程加快了针对不同脑部几何结构和组织电导率的多个个体构建患者特异性模型的速度，同时也便于评估将不同电极设计植入现有模型所产生的影响。本方案中描述的图像预处理步骤尚未完全自动化，可能需要长达一天的处理时间。然而，一旦完成，这些步骤生成的数据被视为静态数据，即在仿真过程中不会被修改。模型构建的自动化源于系统能够将这些数据自动应用于患者特异性的有限元模型（FEM），而无需人工干预。用于执行全部模型生成、仿真和可视化步骤的 SCIRun 网络只需构建一次，因此，为新患者生成患者特异性模型时，仅需重新执行图像预处理步骤即可。

建模流程在生成结果方面的性能提升，得益于将网格生成、生物电场计算以及解的可视化整合到单一软件环境中。采用现有的建模技术（如自适应网格生成），可在电极周围生成较高密度的网格，而在远离电极的区域生成较低密度的网格，从而缩短有限元法（FEM）模型的构建和求解时间。软件 SCIRun 还支持网格生成和生物电场计算的自动化。用户通过交互式控件自定义电极移动时，会触发根据更新后的电极位置重建新网格的过程，包括为新的电极位置调整边界条件和电导率值。

在将电极位置整合到有限元模型（FEM）之前，电极几何结构被视为脑容积内可自由移动的物体。这种建模方法在网格构建中的一个关键优势是可轻松地将多个电极插入模型中。例如，可将同一电极几何结构的第二个副本放置在距离原位置数毫米处，两个电极都将被包含在有限元模型中。在近期研究中，已有将两个电极近距离植入以治疗多发性硬化症震颤的案例¹³，并在非人灵长类动物实验中用于探索有效的刺激靶点²¹。使用多个电极的优势在于能够在更大区域内更精确地控制组织中产生的电场。通过在两个电极的激活触点上同时进行刺激，可将电流引导至目标区域，并避开可能引起不良副作用的区域。在目标靶点确切位置未知的情况下，对更大区域实现更精细的刺激控制，对于探索不同刺激位置尤为有用，这正是许多新兴脑深部电刺激（DBS）治疗应用所面临的情况。然而，由于参数空间本身已较大，增加电极数量进一步扩大了参数空间，因此确定实现治疗性刺激的参数比单电极情况更具挑战性。

我们设想，这种交互式建模工具可在脑深部电刺激（DBS）植入术的术前规划过程中提供帮助。关于电刺激在脑组织中作用范围的反馈信息，可使外科医生在手术计划中调整电极位置，从而对目标区域提供有效的治疗性刺激。尽管DBS治疗是开发该工具的主要动因，但本文介绍的技术也可应用于其他有限元生物电场模型，适用于不同的刺激或记录范式。例如，用于治疗抑郁症的经颅直流电刺激²²，或用于癫痫治疗的深部电极应用²³，在确定最佳刺激位置以实现治疗效果方面，面临着与DBS相同的挑战。皮层脑电图（electrocorticography）是一种在脑表面布设电极阵列以识别癫痫发作起始区域的记录技术，其挑战在于如何确定电极的放置位置，以有效记录脑内目标区域的信号²⁴。所有这些应用均依赖于电极的位置，同时又需应对电流在脑组织中传导路径的不确定性。本文所介绍的技术降低了构建和使用计算模型的负担，同时为非建模专业背景的研究人员和临床医生在使用这些设备时提供了有意义的反馈信息。

材料

本文使用的材料清单

姓名	公司	目录编号	评论
FreeSurfer	哈佛大学 Athinoula A. Martinos 生物医学成像中心		https://surfer.nmr.mgh.harvard.edu/
3D Slicer	布莱根妇女医院，哈佛大学		https://www.slicer.org/
SCIRun	犹他大学整合生物医学计算中心		http://www.sci.utah.edu/cibc-software/scirun.html