周 鹏，贺威忠，林冬冬，葛家怡，王明时

(天津大学精密仪器与光电子工程学院，天津 300072)

1985年，Barker等[1]将平面线圈置于正常人运动区的头皮上，观察到手肌抽动，用表面电极在小指外展肌记录到运动皮质诱发电位，这种方法后来被称为经颅磁刺激(transcranial magnetic stimulation，TMS)．TMS克服了电刺激会带来患者不适的缺点，研究证明，TMS可通过对神经系统Ca2+活动、神经递质和肽类物质代谢以及免疫功能等广泛调节作用治疗抑郁症、老年痴呆、帕金森氏、偏头痛、多发性硬化症和癫痫等疑难病症，避免了药物或手术治疗的痛苦和创伤，使得 TMS得到了迅速商业化，并加速了对它的研究[2-4]，脑磁刺激是在体外经颅使用强脉冲磁场实现的，耦合了某频率的脉冲磁场对生物体的积极影响[5-7]．笔者对 TMS在大脑中产生的电磁场进行分析，了解 TMS下大脑内电磁场的变化规律和大小分布对医学临床应用和探寻磁刺激生物效应具有指导意义．在临床中，脑内特定部位的磁刺激强度需要量化，这就对线圈放置位置和参数提出了要求．

1 经颅磁刺激原理和有限元方法

1.1 经颅磁刺激原理

TMS的脑磁刺激是在脑外经颅施加强脉冲磁场实现的．细胞膜保持一个电位差，静态细胞的跨膜电位差是-70,mV．磁刺激仪电路加载到线圈上产生的脉冲磁场在脑内感应出外加电场，电场叠加到细胞膜两侧可以改变细胞膜电位差，因此外加电场能够极化细胞膜，激活可兴奋组织．当感应电流超过神经组织兴奋域值时，磁刺激就像电刺激一样刺激相应部位的组织．由麦克斯韦方程组时变磁场产生感应电场和感应电流见式(1)，知道磁场大小和分布后就可以推导出感应电场和电流[8]．

式中：σ为电导率；r为矢径．

1.2 有限元方法原理

大脑电磁场的计算问题实际上是一个求解位场边值的问题，即求解同时满足拉普拉斯方程又满足边界条件的位函数．由能量极值定理知道，求解边值问题的解函数等价于求解满足同一边界条件的能量函数的极值函数．在大脑和线圈这样的复杂电磁场的整体场域内求得能量泛函的极值函数是做不到的，而有限元方法以变分理论为基础，通过区域剖分和分片插值，把二次泛函的极值问题转化为一般多元函数的极值问题，而后者就是一组多元线性代数方程的求解．利用计算机求解矩阵即可得出最终的数值解．

2 经颅磁刺激电磁场分析

使用有限元方法[9-10]对线圈产生脉冲磁场进行仿真分析．为尽可能准确真实地分析刺激的全过程，必须解决两个关键问题：①头模型[11]和线圈模型能够对现实原物进行模拟；②能够模拟出经颅磁刺激下线圈内电流变化的过程．

2.1 头部和线圈的建模

2.1.1 真实头模型

随着对大脑电磁场研究的深入，头模型也从最初的简单均匀单层球模型逐渐发展到了现在的 4层同心球模型．为此提出一种真实头模型[12-13]，真实头模型是指对CT和MRI等大脑断层图像进行处理并三维重建出的头模型[14-15]．真实头模型能更加真实地反映大脑结构，因此，使用真实头模型能使分析更加符合现实情况．

采用天津市第一中心医院核磁室的165张256×256真人头部MRI图像作为原始数据，使用医学图像软件对图像进行处理，根据不同灰度阈值对组织边缘进行提取，把头部分成头皮、颅骨、脑脊液和大脑4层组织进行三维重建．为了能将得到的三维几何实体导入有限元分析软件，需要做一些平滑处理，在处理过程中，只对大脑和颅骨的边缘进行三维重建，脑脊液采用布尔运算生成，头皮采用膨胀生成．图1分别显示了大脑、颅骨的三维重建模型和导入有限元软件ANSYS剖分后的有限元大脑模型、有限元头模型．

图1 真实头模型Fig.1 Real head model

真实头模型的 4层由外到内分别代表头皮、颅骨、脑脊液和大脑．各层的电导率分别为：头皮 0.33 S/m，颅骨0.042,S/m，脑脊液1.0,S/m和皮层 0.33,S/m．磁导率都为1,H/m．

2.1.2 线圈模型

实际线圈中的一些参数，如绕向、线径 xj、最大外径wj、材料和额定电流等在有限元分析软件中都不能得到提供．需将这些参数转化为有限元软件可以识别的参数：绕向通过规定电流流向解决；材料通过设定电阻率和磁导率解决；计算结果中，线圈的最大电流如果超过额定电流相应的参数，将被视为无效参数，其他参数转化为

式中：F和 N都是有限元软件中的参数，分别代表填充率和匝数；A是圆形线圈的横截面积．图 2是线圈有限元剖分图．

图2 线圈的有限元模型Fig.2 Finite element model of coil

2.2 TMS脉冲刺激过程分析

磁刺激仪主要由储能电容、线圈和控制电容充放电的电路组成．其原理可等效成 RLC串联电路，其中R、L、C分别是线圈的电阻、电感和储能电容值．对大脑起作用的主要是电容的放电过程．在这个放电过程中，线圈电流是以指数衰减形式变化的，并且对于不同的线圈参数，电流的变化也分为3种情况：非振荡过程、振荡过程和临界阻尼过程，其电流公式分别为

式中 U0指的是电容初始的脉冲电压．图 3显示了在磁刺激分析中两组不同线圈参数所导致的振荡过程(图 3(a))和非振荡过程(图 3(b))的线圈电流变化图(临界阻尼过程电流变化与非振荡过程类似)．可以看出 TMS中电流和磁场的变化规律是复杂的．所以为了更好地了解脉冲刺激的过程和线圈在其中起到的效果，须将电路考虑进去进行电路磁场混合分析(利用ANSYS电路模块进行电路建模)．

图3 线圈电流变化曲线Fig.3 Current curves of coil

考虑电路后，再观察变化的电流导致的磁场变化规律，图4列出了对应非振荡情况脑内几个点的磁场变化图．结果显示，磁场随电流变化而变化，变化曲线、最大值、最小值出现的时间与电流变化过程的相应特征是一致的．实际上，由于大脑与空气的导磁系数相类似，大脑内的任意其他点的磁场也都会有相同的变化规律(振荡或临界阻尼情况类似)．

图4 非振荡情况时各点的磁场变化Fig.4 Magnetic field in several places under non-oscillation condition

大脑各处的磁场的变化规律是一致的，区别只在于大小的不同．考虑到线圈额定电流参数对电流的要求，选取电流最大值出现的时间来观察磁场分布(因为脑内各点磁场的变化规律是相同的，任意时刻下脑内磁场大小的相对分布规律不变，所以只需选择某个时刻来观察)．图 5显示了大脑表面和沿大脑线圈中心连线大脑半剖面的磁场分布．结果显示在大脑表面的磁感应强度以线圈中心在大脑投影处最强，向边缘减弱．大脑内部离线圈最近的地方磁感应强度最大，离线圈越远，进入脑内越深，磁感应强度衰减越快．

图5 磁场分布Fig.5 Distribution of magnetic field

3 线圈参数对电磁场的影响分析

3.1 电路性质参数

线圈连接于电路中并由电路中的储能电容供电，故将与电路有关的参数都归入电路性质参数一类．这样的参数有：储能电容 C、储能电容初始电压U0、线圈电阻R、线圈电感L、最大电流 Imax．这些参数共同决定了线圈电流的大小和变化规律．

从式(4)、(5)、(6)可以看到线圈电流的大小、线圈电流以何种形式变化，电流变化曲线性质完全由U0、C、R、L这 4个参数共同决定．由前面的分析可知，磁场变化规律与电流一致，因此电流变化形式一旦确定，磁场变化形式也随之确定．

3.2 几何参数

线圈内径 r、外径 R、厚度 h、体积 V和横截面积A由于与线圈的结构有关所以被称为几何参数．几何参数主要起着两方面的作用：①通过与性质参数共同作用决定线圈的电路性质参数；②影响电磁场的分布．选取一组线圈几何参数为：内径 0.02,m，横截面宽度(外径减去内径)0.02,m，厚度 0.02,m．将其他参数保持不变，内径以 0.02,m 的差量增加，最终得到 6组线圈参数．图 6(a)为线圈在空气中产生的磁力线截面图，选取图 6(a)中所示观察路径，观察 6种不同大小线圈在观察路径上产生的磁场大小分布，如图6(b)所示．

图6 路径磁场分布Fig.6 Distribution of magnetic field along path

从图 6中可以看到，内径在 0.02,m 以下的线圈产生的磁场在中心处为最强，而向远处衰减．内径为0.04,m及以上的线圈产生的磁场则在对应内径的距离处为最强．磁力线围绕电流，线圈横截面附近磁力线最密集．小线圈由于内径较小，其磁力线在中心处汇聚而显得较大线圈中心处更密集，从而能够实现中心聚焦功能．

3.3 性质参数

线圈的性质参数是指选用材料的性质和制作方法，包括匝数 N、导线线径 xj、导线外径 wj、导线材料m、填充率f、绕向及额定电流Icur．

一旦性质参数和形状参数确定，那么线圈的电阻和电感也被确定，进而电流的大小和变化形式也被确定．所以性质参数主要通过电路性质参数来对电磁场产生影响，例如电阻与匝数平方成正比，电感与匝数成正比，增大匝数可使线圈电流的变化形式从图3(a)中的振荡形式过渡到图3(b)的非振荡形式．

4 电磁场分析系统

4.1 分析系统的组成和功能

以对线圈参数在电磁场中的作用的分析为理论基础，便可设计出电磁场分析系统．图 7为系统的流程图．系统包括两个子系统：电磁场自动分析系统和线圈参数获取系统．

大型有限元分析软件一般安装于大型服务器上，对客户终端机发来的指令进行响应并后台计算．电磁场自动分析系统便借鉴了这种思想，输入线圈和电路参数便可以实现经颅磁刺激电磁场分析和结果观看．

系统利用数据库技术实现了参数的反向自动获取(根据所需的磁场大小，给出符合要求的线圈参数)，并可以根据系统提供的优化准则选择最优参数．将磁刺激中常用线圈相对大脑放置位置，线圈参数以及与之相对应的所有结果数据作为一组数据库数据单元，数据库里则存储了总计 10,380组数据单元的信息(数据库利用自动分析系统提前制作完成)．

图7 系统流程Fig.7 System flow chart

4.2 系统应用

经研究表明，通过将耦合了特定频率并达到刺激强度要求的脉冲磁场施加于人脑中缝核能对睡眠产生显著影响．笔者所在实验室研发的睡眠仪便是以此为理论基础研制而成．

根据睡眠仪需要的在中缝核上的磁场大小，经过参数自动获取系统反向获取所有符合要求的数据单元(包含了线圈大脑相对摆放位置、线圈参数及其结果等信息)．再通过系统提供的优化准则如刺激深度最深准则，从中筛选出线圈中心与中缝核距离最远的数据单元．将这组数据单元中的线圈参数输入到自动分析系统便可以观察这组参数下脑内磁场的分布．

5 线圈磁场检验

将线圈置于空气中，分析其在空气中产生的电磁场．首先从数据库中取出一组实际线圈参数制作一款线圈，其实测电阻为 8,Ω，电感为 0.22,H．然后仿真计算线圈模型的电阻和电感，将实际参数按照前面提到的方法转化为软件参数后，计算得出的结果为：电阻为 8.027,Ω，电感为 0.227,548,H．结果说明模型建立准确．再将计算出的磁场值和现实检测得到的磁场值相比较，检测仪器为 PF-035,B型数字特斯拉计．表 1是线圈中垂线上不同点两组磁场值的比较(施加30,V、3,Hz正弦电压)．

表1 磁场对照Tab.1 Magnetic field comparison

从表中可知两种结果变化趋势是一致的，数据平均误差在 15%内．考虑到数字特斯拉计的精度、探针的形状和实验室所处的电磁环境，理论分析结果准确可靠并能够指导实际线圈设计．

6 结 语

通过使用真实头模型对大脑进行模拟，结合有限元方法，经颅磁刺激电磁场分析系统实现了对 TMS在大脑中产生的电磁场的自动分析．同时根据临床需要，分析系统能够自动为使用者提供符合刺激磁场要求的关于线圈参数和摆放位置的参考意见．实际使用表明：经颅磁刺激电磁场分析系统使用方便，分析准确．

［1］Barker A T，Jalinous R，Freeston I L. Non-invasive magnetic stimulation of human motor cortex[J]. Lancet，1985，325(8437)：1106-1107.

［2］Kahkonen S，IImoniemi R J. Transcranial magnetic stimulation：Applications for neuropsychopharmacology[J]. Journal of Psychopharmacology，2004，18(2)：257-261.

［3］Bailey C J，Karhu J，Llmonieml R J. Transcranial magnetic stimulation as a tool for cognitive studies[J].Scandinavian Journal of Psychology，2001，42(3)：297-306.

［4］Munchau A，Bloem B R，Irlbacher M R，et al. Functional connectivity of human premotor and motor cortex explored with repetitive transcranial magnetic stimulation[J].The Journal of Neuroscience，2002，22(2)：554-561.

［5］赵梅兰，马 强，曹晓哲，等. 电磁脉冲对大鼠学习记忆能力和长时程增强的影响[J]. 中华物理医学与康复杂志，2001，23(2)：69-71．

Zhao Meilan，Ma Qiang，Cao Xiaozhe，et al. Effect of electromagnetic pulse on learning-memory abilities and long-time potentiation in rats[J].Chinese Journal of Physical Medicine and Rehabilitation，2001，23(2)：69-71(in Chinese)．

［6］郭明霞，王明时，王学民. 脉冲磁场对大鼠记忆能力和海马神经递质的影响[J]. 中国生物医学工程学报，2002，21(2)：179-181．

Guo Mingxia，Wang Mingshi，Wang Xuemin. The influence of the pulse magnetic field brings to the rats’memory capacity and neurotransmitters of the hippocampus[J].Chinese Journal of Biomedical Engineering，2002，21(2)：179-181 (in Chinese).

［7］路会生，汪 曣，王明时，等. 经颅磁刺激诱导人愉快状态的研究[J]. 天津大学学报，2007，40(5)：616-622．

Lu Huisheng，Wang Yan，Wang Mingshi，et al. Human pleasure induced by transcranial magnetic stimulation[J].Journal of Tianjin University，2007，40(5)：616-622(in Chinese)．

［8］Roth B J，Saypol J M，Hallett M，et al. A theoretical calculation of electric field induced in the cortex during magnetic stimulation[J].Electroencephalography and Clinical Neurophysiology，1991，81(1)：46-56．

［9］Wang Weiping，Eisenberg S R. A three-dimensional finite element method for computingmagnetically induced currents in tissues[J].IEEE Transactions on Magnetics，1994，30(6)：5015-5023．

［10］Wagner T A，Zahn M，Grodzinsky A J，et al. Threedimensional head model simulation of transcranial magnetic stimulation[J].IEEE Transactions on Biomedical Engineering，2004，51(9)：1586-1598．

［11］Cerri G，De Leo R，Moglie F，et al. An accurate 3D model for magnetic stimulation of the brain cortex[J].Journal of Medical Engineering and Technology，1995，19(1)：7-16．

［12］Mouchawar G A，Nyenhuis J A，Bourland J D，et al.Magnetic stimulation of excitable tissue：Calculation of inducededdy-currents with a three-dimensional finiteelement model[J].IEEE Transactions on Magnetics，1993，29(6)：3355-3357.

［13］郑建斌，霍小林．经颅磁刺激中大鼠真实头模型感应电场分布的研究[J]．北京生物医学工程，2006，25(5)：490-492.

Zheng Jianbin，Huo Xiaolin. Analysis of the distribution of induced electric fields during transcranial magnetic stimulation using the real rat′s head model[J]．Beijing Biomedical Engineering，2006，25(5)：490-492(in Chinese).

［14］郑旭媛，万柏坤. 脑电逆问题的研究进展[J]. 国外医学科学：生物医学工程分册，2001，24(4)：179-183．

Zheng Xuyuan，Wan Baikun. The progress of the inverse problem for EEG[J].Biomedical Engineering Foreign Medical Science，2001，24(4)：179-183(in Chinese)．

［15］Taschereau R，Chow P L，Cho J S，et al. A micro CT X-ray head model for spectra generation with Monte Carlo simulations[J].Nuclear Instruments and Methods in Physics Research Section A：Accelerators，Spectrome-ters，Detectors and Associated Equipment，2006，569(2)：373-377．