经颅直流电刺激治疗时颅内电场数值仿真

刘习杰，谭云华

1.宜春学院电子通信系，江西宜春 336000；2.北京大学应用电子学研究所，北京 100871

前言

经颅直流电刺激（transcranial Direct Current Stimulation,tDCS）是一种神经调控技术，作为非侵入式物理治疗，越来越多地用于精神疾病和疼痛管理。美国食品药品监督管理局（FDA）报告认为tDCS具有操作简单，使用成本低廉的优点，是一种相对新颖，较少设备依赖的脑刺激形式，对患者的终生和慢性治疗更可行［1］。

Knotkova等［2］研究论文指出，连续5 d tDCS临床治疗慢性疼痛后，患者自我评估减轻疼痛有效。神经性疼痛组中有76%患者表示有效，包含面部疼痛、疱疹神经痛和复杂的局部疼痛综合征。受伤性疼痛治疗组中77%患者表示有效缓解疼痛，包括腰痛。头痛治疗组中，83%患者表示治疗有效。Knotkova临床试验结果还表明tDCS对多种其它慢性疼痛也具有镇痛效果。tDCS 还能治疗多种精神健康疾病。Nitsche 等［3］、Bennabi 等［4］研究表明tDCS 对抑郁症有改善作用。Herrera-Melendez［5］运用tDCS 对患有精神分裂症以及药物依赖上瘾患者进行治疗，试验具有疗效结果。除此之外，Fregni［6］临床试验验证tDCS可以增强记忆。

tDCS 通过置于头皮的一对或多对电极，以弱直流电无创地刺激大脑皮层的兴奋性。tDCS电极上常施加小于2 mA 的电流。有阳极刺激（anodal stimulation）和阴极刺激（cathodal stimulation）两种工作模式，见图1［7］。阳极刺激产生内向电流，膜电位朝着去极化变化（兴奋性刺激），而阴极刺激引起外向电流，膜电位朝着超极化（抑制性刺激）进行［5］。神经元对刺激活动的响应与治疗电极设定极性相关。颅内极性根本上由颅内电场的空间分布决定，与神经元所在空间位置相关；所以，对颅内电场分布进行数值量化分析，有助于深入理解tDCS治疗的有效性。

有研究指出，大脑个体之间的差异导致tDCS治疗结果变化，使得不同临床研究结果的可重复性普遍较低［8-9］。如果对颅内电场分布进行数值量化分析，可以直接对比个体治疗靶区的电场强度是否变化，是否高于0.7 V/m的有效刺激阈值，可更准确解释tDCS治疗有效性机制。笔者注意到国内关于tDCS的研究较少讨论治疗靶区电场强度，比方胡昔权教授［10］在讨论tDCS结合虚拟现实训练治疗脑卒中患者上肢功能时，就没涉及治疗靶区的有效电场刺激强度。

将tDCS 具体治疗实施方案比喻成治疗处方，那么讨论tDCS 治疗处方的功效，会涉及到电流强度和刺激持续时间等刺激剂量参数。临床试验处方多遵循安全电流2 mA 的限制，采用25 和35 cm2两种规格的电极片，治疗持续时间为5～30 min。处方的不同电流强度和治疗持续时间设定，会影响治疗有效性［11］。tDCS治疗处方还关联电极放置位置和馈电极性。电极位置借用10-20脑电图电极系统标记，如图2所示。作为示例，在表1中，给出了6 种治疗处方。其中重度抑郁症（MDD）处方，阳极连F3，阴极连FP2、F8；精神分裂症治疗处方，阳极连F3 和FP1，阴极连TP3。表1是文献中用35 cm2电极进行治疗的部分处方［6,12-16］。

表1 用35 cm2电极进行治疗的部分处方Table 1 Partial prescriptions for treatment with 35 cm2 electrodes

1 脑组织分割和电场计算模型的网格生成方法

分析计算颅内电场数值，先得获取3D头颅数据。头颅数据可由Philips Achieva机器采集成年男子头颅MRI 数据。Achieva MRI 机器的空间分辨率为（1×1×1）mm3，数据输出格式为DICOM。该机器可扫描256层，每层视野范围（256×256）mm2；层厚1 mm，无间隙。

取得MRI 数据后，接着通过脑组织分割标定不同脑组织相对区域位置。脑组织分割目的是完成脑内电磁场不同介质边界的定义。脑组织分割可由ScanIP软件进行处理。ScanIP是一款数字3D图像建模工具软件，能完成MRI 图像数据导入，组织分割处理等一系列过程。ScanIP还内嵌网格生成模块，输出结果可以直接导入COMSOL，并在COMSOL 中进行有限元计算。此外，ScanIP 还导入事先由Solidworks设计的tDCS电极设计文件。

在ScanIP 输入成年男子头颅MRI 数据，通过设置灰度，ScanIP 自动按灰度将组织定义为不同脑组织。选择脑组织灰度值需要脑解剖知识，并反复尝试。此脑组织自动分割过程不完善，后续要借助ScanIP 软件内嵌的3D 可视化交互编辑工具，对自动分割的结果做人工编辑，例如处理骨骼、大脑皮层折叠、脑脊液（CSF）等的MRI对比度缺失。编辑过程需脑组织解剖知识，才可正确完成组织分割修改。全部头颅组织分割完成后，得到脑组织区域分布。图3展示重建后3D 模型结果，依次为白质、灰质、CSF、颅骨、脂肪、头皮。

图3 ScanIP中脑组织分割后部分结果Figure 3 Some results of brain tissue segmentation in ScanIP

接着建立头颅网格数据模型，准备有限元数值计算。在此之前，需要导入事先设计的tDCS电极文件，此电极文件可由Solidworks创建。电极导入，放置在选定位置，结果见图4。考虑网格生成参数关乎有限元数值解运算迭代次数和计算精度，因此在网格生成时要合理平衡网格生成密度、总节点数、模型复杂性、准确性之间的关系。选择网格密度即精度配置，对体积小的脑组织，或对重要研究区域，选高精度。对于体积大，或对于非重点研究区域，选低精度。精度设置完成后，在网格创建算法窗口界面单击“+FE free”。

图4 两对电极放置在F3(正极)-P08（负极）和F4(正极)-P07(负极)Figure 4 Two pairs of electrodes are placed at F3(anode)-P08(cathode)and F4(anode)-P07(cathode)

完成网格生成的各项配置，单击“Full model”，启动ScanIP 网格划分，开始创建有限元计算网格。使用四核3.1 GHz CPU 和512 GB 内存计算机，生成网格大约需要3 h。本例生成网格总计23 104 974 个四面体单元，其中头皮5 549 178 个，颅骨3 490 965 个，脂肪3 614 410 个。导出网格文件（COMSOL*.mph扩展文件)，导入COMSOL得到图5所示模型图。

图5 导入COMSOL中的头颅网格模型图Figure 5 Head mesh model imported into COMSOL

2 tDCS治疗时颅内电场数值计算实现

实现电场数值计算，可使用开源或商用软件工具。本文讨论tDCS 治疗处方剂量数值分析，运用COMSOL Multiphysics 软件实现数值求解。COMSOL Multiphysics 是多物理学科计算仿真软件，亦能用于生物领域电磁场工程计算。基于安全的规范，tDCS 治疗处方电极最大电流设为2 mA。电极片大小为（5×7）cm2。当给电极接通直流电流时，形成一个电场边界，电场边界数值计算按照以下的步骤实施。

第一步，将已经生成的头部和电极合并网格数据导入COMSOL，如图5所示。在COMSOL中，给不同脑组织和电极的介电常数和导电率赋值，建成高分辨率3D人头电场边界问题数字模型。重建过程中，脑组织电导率引用文献［17］。皮肤电导率为0.465 S/m；脂肪电导率为0.025 S/m；颅骨电导率为0.010 S/m；CSF电导率为1.650 S/m；灰质电导率为0.276 S/m；白质电导率为0.126 S/m；电极电导率为5.99×10^6 S/m；头皮和电极间的导电胶电导率为4.5 S/m。

第二步，在COMSOL中选择遵循的物理学科。在AC/DC菜单下面选择电流交互界面。因为tDCS研究的电场是稳定，不随时间变化的，故选静态模式。继续在COMSOL中为相应边界设置电流密度和电位。本文讨论的tDCS处方，4个电极分别放置在F3（正极）-P08（负极）和F4（正极）-P07（负极）。此时边界，第一组电极中与头皮接触的阳极表面，其外加法向电流密度设为n⇀·J⇀1=2 mA/25 cm2=0.8A/m2，阴极表面接地电位φ1=0；第二组电极中，阳极法向电流密度设为n⇀·J⇀3=0.8A/m2，阴极接地电位φ3=0，颅内脑组织的各个分界边界均是电流密度法向连续且电位连续。

完成上面设置后，单击COMSOL 中计算按钮。四核3.1 GHz CPU 和512 GB 内存计算机在大约2 h后收敛，算出头颅静态电场分布。颅内电场计算主要技术路径如图6所示。

图6 颅内电场计算的技术路径Figure 6 Technical path of intracranial electric field calculation

3 颅内电场数值计算结果及讨论

用COMSOL 完成数值计算，不同的计算机性能会需要不同的时间，对于内存低配的计算机可能会出现内存不足而中断。图7～图9呈现了颅内电场的可视化图解。图7依次是白质顶视、侧视，灰质顶视、侧视，CSF 顶视、侧视，颅骨顶视、侧视，皮下脂肪，头皮的可视化电场分布。

图7 颅内脑组织电场分布Figure 7 Electric field distribution of intracranial tissues

电场强度在白质和灰质部位为0 V/m（蓝色）～1.8 V/m（红色）。电极片之下区域的电场强度高于0.7 V/m 的有效刺激阈值，表示tDCS能够在这些区域有效调控神经，达到物理治疗的作用。可选择数值-颜色对照的范围，如果最大值选0.7 V/m，就能清楚看到白质和灰质中电场强度高于0.7 V/m 的区域范围（红色），见图8。

图8 颅内灰质和白质高于0.7 V/m的电场分布Figure 8 Electric field distribution of intracranial gray matter and white matter above 0.7 V/m

图9是电场分布横截面图解。图中颅骨和其它脑组织分界面的颜色分界对比清晰，说明颅骨中电场大小明显不同其它脑组织中电场大小，表现为电场在跨颅骨进入到CSF、灰质和白质，经历了快速减弱的过程，最大从18 V/m降为1.2 V/m。

图9 颅内电场在头颅横截面的分布（多层片电场分布,V/m）Figure 9 Distribution of the intracranial electric field at the cross-section of the skull (electric field distribution of the multilayer sheet,V/m)

临床治疗中电极位置偏移导致个体治疗响应存在差异。笔者以tDCS 治疗抑郁症状为例，讨论电极位置偏移对治疗量化剂量的影响。假定治疗的电极放置在F3 和F4，即左背外侧前额叶皮层（DLPFC）和右DLPFC 的正上方。笔者将无偏差的F3 电极位置称之为正位。以正位为中心，假设电极向8个方向偏移：东、西、南、北，以及东南、西南、西北、东北，每次移动为1/4 电极大小（圆盘直径），最大移动一个直径大小（4/4），共32种中央位置偏移方案。完成32个偏移电极放置场景的颅内电场计算，读取每次计算结果中的P 点和Q 点电场数值。结果如图10所示。图10显示在2 mA 激励下，皮质中P 和Q 点电场范围为0.29～0.65 V/m。图11显示此时各位置误差所产生的电场剂量偏差为14.55%～53.85%。

图10 不同电极偏移位置下P点和Q点的电场强度Figure 10 Electric field intensity values at points P and Q with different offsets of the electrode positions

图11 不同电极偏移位置下P点和Q点的电场强度变化率Figure 11 Electric field intensity change rates at points P and Q with different offsets of the electrode positions

针对模型所得数值解的可靠性，即模型计算结果与直接人类颅内测量对照，本文预期有14%的误差。Huang等报道［18］和本文介绍的方法是一致的，他通过10 名癫痫患者头建模计算，对比电场计算数值与测量数值，在头皮位置电场强度相关为r=0.86。

4 结论

颅内电场数值仿真计算可以让我们在从事tDCS治疗时，以及在从事相关仪器研制工作中受益。tDCS优化治疗，即是希望电场高于0.7 V/m 的区域范围和希望调控的区域重合，而在调控外的区域低于0.7 V/m。个体治疗响应会存在差异，临床治疗当中对个体剂量进行优化自然会有帮助。人工经验化的优化思路成本高且盲目。而以头部数字电磁场模型建模，通过精确了解大脑中电流分布，预测了整个组织上产生的电流空间分布是否适合，可以找到tDCS治疗合理优化处方，克服人工优化的缺点。

其次，电极位置误差会成为颅内电场剂量差异的主要来源之一。仿真计算可为tDCS治疗操作选择正确电极位置提供指导。

最后，tDCS 计算模型中显示电场在跨颅骨时经历快速减弱。启示在新的神经调控设备研制，如果tDCS 治疗的电极置于颅骨下，一方面可使电极激励方案具备替代DBS 潜力，另一方面会比目前的tDCS更有效实现目标区域神经调控。有关这两方面的研究，将在后继工作中展开。

致谢：Marom Bikson 教授提供作者使用其脑组织分割完成后模型，特此向Bikson 教授表示感谢。Huang Yu 博士提供了非常重要的意见，使本文得到完善，亦特此表示感谢。