用于视觉假体研究的视网膜神经元计算建模 - 基本方法

我们总结了一个工作流程，以计算方式模拟视网膜神经元对电刺激的反应。该计算模型是通用的，包括自动化步骤，可用于模拟一系列生理场景和预测未来 体内/体外 研究的结果。

计算建模由于其预测 体内 和 体外 系统行为的能力，已成为神经工程中越来越重要的方法。这具有关键优势，即通过提供通常非常精确的生理结果预测来最大限度地减少给定研究中所需的动物数量。在视觉假体领域，计算建模具有一系列实际应用，包括为植入式电极阵列的设计提供信息，以及预测可能通过从所述阵列传递电脉冲而引发的视觉感知。文献中描述的一些模型结合了三维（3D）形态来计算电场和感兴趣的神经元或神经网络的电缆模型。为了增加这种两步法对那些在计算建模方面经验有限的研究人员的可访问性，我们提供了一个视频，介绍了要采取的基本方法，以便构建计算模型并利用它来预测 通过部署 的刺激方案的生理和心理物理结果。 视觉假体。该指南包括在有限元建模（FEM）软件中构建3D模型的步骤，在多区室神经元计算软件中构建视网膜神经节细胞模型的步骤，然后将两者合并。用于数值求解物理方程的有限元建模软件将用于求解组织电刺激中的电场分布。然后，使用专门的软件来模拟神经细胞或网络的电活动。要遵循本教程，需要熟悉神经假体的工作原理以及神经生理学概念（例如，动作电位机制和对霍奇金-赫胥黎模型的理解）。

视觉神经假体是一组向视觉通路中的神经细胞提供刺激（电、光等）以产生光幻觉或看到光的感觉的装置。这是一种治疗策略，对于由退行性视网膜疾病引起的永久性失明的人，已经在临床上使用近十年。通常，一个完整的系统将包括一个捕获用户周围视觉信息的外部摄像头，一个电源和计算单元，用于处理图像并将其转换为一系列电脉冲，以及一个植入的电极阵列，该电极阵列连接神经组织并将电脉冲传递到神经细胞。工作原理允许视觉神经假体放置在从视网膜到视觉皮层的视觉通路的不同位置，只要它位于受损组织的下游。目前视觉神经假体的大多数研究都集中在提高刺激的功效和提高空间敏锐度以提供更自然的视觉上。

在提高刺激功效的努力中，计算建模一直是验证假体设计和模拟其视觉结果的一种经济和时间有效的方法。自1999年以来，该领域的计算建模越来越受欢迎，因为Greenberg¹模拟了视网膜神经节细胞对细胞外电刺激的反应。从那时起，计算建模一直用于优化电脉冲²，³的参数或电极^4，5的几何设计。尽管复杂性和研究问题各不相同，但这些模型通过确定介质（例如神经组织）中的电压分布并估计附近神经元因电压而产生的电响应来工作。

导体中的电压分布可以通过求解所有位置的泊松方程⁶ 来找到:

其中 E 是电场， V 是电势， J 是电流密度，σ 是电导率。公式中的 表示 梯度运算符。在静止电流的情况下，对模型施加以下边界条件:

其中 n 是表面的法线，Ω 表示边界， I₀ 表示比电流。它们一起在外部边界产生电绝缘，并为选定的边界创建电流源。如果我们假设单极点源位于具有各向同性电导率的均质介质中，则任意位置的细胞外电势可以由⁷计算:

其中 I_e是电流，是电极和测量点之间的距离。当介质不均匀或各向异性，或者电极阵列有多个电极时，使用计算套件来数值求解方程会很方便。有限元建模软件⁶ 将体积导体分解为称为"单元"的小部分。这些元素相互连接，使得一个元素的变化影响其他元素的变化，它求解了用于描述这些元素的物理方程。随着现代计算机计算速度的提高，这个过程可以在几秒钟内完成。一旦计算出电势，就可以估计神经元的电响应。

神经元以电信号的形式发送和接收信息。这种信号有两种形式 - 分级电位和动作电位。分级电位是膜电位的临时变化，其中膜上的电压变得更加正（去极化）或负（超极化）。分级电位通常具有局部效应。在产生它们的细胞中，动作电位是可以沿着轴突长度长距离传播的全有或全无反应。梯度电位和动作电位对电气和化学环境都很敏感。当跨膜电位阈值越过时，各种神经元细胞类型（包括视网膜神经节细胞）可以产生动作电位尖峰。动作电位尖峰和传播然后触发信号到下游神经元的突触传递。神经元可以建模为分为圆柱形段的电缆，其中每个段由于脂质^双层膜8而具有电容和电阻。神经元计算程序⁹ 可以通过将细胞离散为多个区室并求解数学模型¹⁰ 来估计可电兴奋细胞的电活动:

在该方程中，C_m是膜电容，V e，n是节点n处的细胞外电位，V_i，n是节点n处的细胞内电位，R n是节点n处的细胞内（纵向）电阻，I离子是通过节点 _n处离子通道的_离子电流。当刺激处于活动状态时，FEM 模型中的 V 值对于神经元中的所有节点实现为 V_e，n。

来自离子通道的跨膜电流可以使用霍奇金-赫胥黎公式¹¹进行建模:

其中g i是通道的比电导，V _m是跨膜电位（V _i，n - V_e，n），E离子是_离子通道的反转电位。对于电压门控通道，例如Na通道，引入了描述通道打开或关闭概率的无量纲参数m和h:

其中 是特定离子通道的最大膜电导率，参数 m 和 h 的值由微分方程定义:

其中 α_x 和β_x 是定义离子通道速率常数的电压相关函数。它们通常采用以下形式:

这些方程中的参数值，包括最大电导率，以及常数A、B、C和D，通常从经验测量中找到。

使用这些构建基块，可以按照所述步骤构建不同复杂性的模型。当泊松方程无法解析求解时，例如体积导体中的不均匀或各向异性电导或电极阵列的几何形状复杂时，FEM 软件非常有用。求解细胞外电位值后，可以在神经元计算软件中对神经元电缆模型进行数值求解。结合这两个软件，可以将复杂的神经元单元或网络计算成非均匀电场。

将使用上述程序构建脉络膜上刺激下的视网膜神经节细胞的简单两步模型。在这项研究中，视网膜神经节细胞将受到一系列幅度的电流脉冲。细胞相对于刺激的位置也不同，以显示距离-阈值关系。此外，该研究还包括使用不同尺寸的刺激电极¹²对皮质活化阈值的体内研究对计算结果的验证，以及显示电极 - 神经元距离与激活阈值¹³之间关系的体外研究。

1. 建立用于电势计算的有限元模型

1. 确定仿真步骤和模型的复杂性
   注意:第一步的目的是阐明建模的目的，这将确定模型的必要元素和模拟过程。需要考虑的重要一点是需要通过模型显示的神经细胞的行为，以及需要什么测试协议来证明这种行为。这项研究显示了细胞外刺激神经元的距离-阈值关系，以及电极大小-阈值曲线。为此，需要将神经细胞模型划分为对细胞外电压敏感的不同部分（以结合神经元中形态和生物物理参数的变化），并模拟一系列电极尺寸和位置。
   1. 定义研究问题和实验变量。
      1. 定义研究问题和测试协议以指导模型的构建。最好从一个明确的问题开始，并建立一个尽可能简单的模型来回答它。
   2. 确定要包含在完整模型中的必要元素
      注意:在这种建模方法中，细胞被视为浸入导电介质，即生物组织。电刺激发生在这个"体积导体"，即介质上，导致电势分布。
      1. 根据要解决的研究问题和变量，确定是否需要这两个元素（FEM 和神经元电缆模型）。例如，如果建模需要一个可以简化为点源且介质均匀的单个电极，则可能不需要 FEM，并且可以对细胞外电场进行分析计算来代替它。
2. 下载并安装软件
   注意:该研究使用了 材料表中指定的软件应用程序版本（COMSOL，NEURON和Python Anaconda）和硬件。如果使用不同版本的软件/硬件，则步骤或结果可能会略有不同。
   1. 如有必要，请下载适合计算机操作系统的软件并购买许可证。确保下载所有必需的仿真模块并安装所有软件。
3. 收集有关要建模的组织和细胞解剖结构的数据
   注意:对于这种方法，解剖学和生物物理参数取自经验发现。由于数据不可用，计算模型通常会混合在不同物种中测量的参数。为了模拟脉络膜上刺激，模型中需要包括刺激电极和参比电极之间的组织层。
   1. 从组织学研究中收集组织的解剖结构。
      1. 在此模型中，包括脉络膜、视网膜组织和玻璃体结构域，其中每个结构域都建模为矩形棱柱，以便于构建模型。从已发表的组织学数据¹⁴ 中收集平均视网膜组织厚度，以稍后用作每个棱镜的高度。
   2. 从细胞染色或公共神经元数据库中收集单细胞形态数据。
      1. 从 NeuroMorpho.org 等数据库中下载详细的神经元形态，该数据库提供搜索元数据功能，可根据物种、大脑区域、细胞类型等查找相关神经元。在本研究中，通过在物种字段中输入兔子>新西兰白和在大脑区域字段中输入视网膜，找到郭的 OFF RGC 模型 （D23WM13_27_1-OffRGC_msa）¹⁵。单击模型并下载 .swc 文件。
4. 收集建模细胞的生物物理数据
   注意:生物物理参数包括每个组织层的电导率值以及神经膜和离子通道的电参数。
   1. 由于数据的可用性，请使用从兔^子16 获取的电导率值作为组织模型，而离子通道的动力学基于虎蝾螈视网膜¹⁷的Sheasby和Fohlmeister模型。
5. 在 FEM 软件中构建组织和电极有限元模型的几何形状
   注意:组织和电极阵列的几何形状都会影响电位分布，进而影响神经细胞的行为。因此，构建细胞所在的培养基以及电极的真实几何形状非常重要。本教程中使用的 FEM 软件具有一个 GUI，可以轻松构建模型几何图形。
   1. 在软件的 GUI 中设置 FEM 模型:
      1. 运行 FEM 软件，然后单击 3D >模型向导。在"选择物理场"列表框中，展开"AC/DC">"电场"和"电流>电流 （ec）"，然后单击"添加"。单击"研究"并在"通识研究"选项下添加"固定研究"，然后单击"完成"（补充图1）。
   2. 设置电极的单位和几何参数。
      1. 在模型生成器树上，单击 参数 1。在表中，在"名称"字段中键入"elec_rad"，在"表达式"字段中键入"50"，以创建半径为 50 个单位的电极。然后，单击 几何形状，并将长度单位更改为μm，因为典型视网膜神经节细胞的体细胞直径约为10μm（补充图2）。
   3. 使用块域创建组织层
      注意:为了构建模型几何形状，使用了三个块来表示眼睛中的不同结构。方框 1 代表脉络膜，块 2 代表视网膜组织，框 3 代表玻璃体。
      1. 右键单击 几何 1 >块 以创建块域。再重复此步骤两次，总共创建三个块。对于所有块，将 深度 和 宽度 都设置为 5，000 μm，并将 基础 选项（ 在位置下）更改为 居中。为每个块分配以下 高度 （在 大小 和 形状下）和 z （ 在位置下）值:
         模块 1:高度 = 112 μm，z = 0 μm
         块 2:高度 = 151 μm，z = 131.5 μm
         模块 3:高度 = 5，000 μm，z = 2，707 μm
   4. 创建工作平面以向模型添加电极
      1. 在"模型树"中的几何 1 上单击鼠标右键，然后选择"工作平面"。单击工作平面 1 并将平面类型更改为面平行，单击平面类型下方的激活选择按钮，然后选择块 1 的底面（大牌 1 > 1）。
   5. 在工作平面上绘制圆盘电极
      1. 单击"工作平面 1"下的"平面几何图形"，然后单击主工具栏中的"草图"。选择"圆"，单击"图形"选项卡中矩形内的任意位置，然后拖动以创建圆盘电极。将半径更改为"elec_rad"μm，xw 和 yw 更改为 0 μm，然后单击"全部构建"。
   6. 为每个域指定材料属性
      注意:通过按照构建几何体的步骤，模型将被分成几个"域"，这些"域"是构成完整几何体的单个 3D 零件。应为每个域分配一个电导率值，以计算整个模型中的电场分布。
      1. 在"模型树"中，右键单击" 材料">"空白材料"，然后单击 "材料 1 "并将 "选择" 更改为 "手动"。
      2. 单击图形窗口中的属性域，以便仅选择属性域 1。选择 材料属性 > 基本属性 > 电导率，单击添加到 材料 按钮，然后将电导率值更改为 0.043 S/m¹⁵。
      3. 对域2和3重复这些步骤，电导率值分别为0.7 16和1.55¹⁶ S/m（补充图3）。
   7. 对 3D 模型进行网格划分:要对模型进行网格划分，请转到"模型树"，然后右键单击 网格 1 >自由四面体。单击 自由四面体 1 ，然后选择 全部构建。
      注意:网格划分过程将整个几何体分解为更小的"元素"（元素是模型几何体的虚拟段，其中物理方程以数值求解）。从理论上讲，使用较小元素进行网格划分可以提高近似的准确性，但在计算上是详尽的。一种常见的做法是使用稀疏网格启动模型并记录仿真结果，然后每次使用较小的网格元素连续重复仿真并比较结果。当计算结果与后续细化步骤没有显著差异时，细化过程可以停止。
      1. 评估网格质量:右键单击 网格 1 并选择 统计数据 以显示元素质量的直方图。按照下面的网格细化步骤来提高元素的质量。
         注意: 使用默认网格划分可能会产生许多低质量的元素，这反过来又会导致计算不准确。在大多数情况下，需要一定程度的网格细化。
      2. 细化电极周边的网格
         注意:电场可能突然变化的区域通常需要更精细的网格。在这里，使用边缘分布特征在电极周边添加了更密集的网格。
         1. 首先，删除现有的自由四面体 1 网格。然后，右键单击 网格 1 >分布，单击 分布 1，将 几何实体级别 更改为 边缘 ，然后选择 边 19-22 （圆盘电极的周长）。
         2. 将分布类型设置为固定元素数，并将元素数字段更改为 elec_rad*3/10 以使元素合理小。
      3. 细化脉络膜和视网膜组织的网格
         1. 在"模型树"中，右键单击 网格 1 >扫描。单击扫描 1。选择 域 1 和 2。然后，右键单击 网格 1 >自由四面体，将 几何实体级别 设置为 剩余，然后单击 全部构建。（可选:再次检查元素质量直方图，以确保低质量元素已按比例减小）。
6. 将物理场应用于有限元模型
   注意:FEM 软件中的"物理场"是需要分配给模型的一组数学方程和边界条件。它是在运行模型的过程中对联立方程组的解的计算。选择应用于几何的物理场取决于所模拟的物理现象。例如，该模型中使用的电流物理场观察电场分布而忽略了磁（感应）现象。如果要解决其他物理问题（例如，温度分布，机械应力等），则可以将其他物理应用于几何。
   1. 物理场的选择和边界条件的应用
      注意: 如果要应用恒压脉冲，则应将浮动电位边界条件替换为电势边界条件。
      1. 展开"模型树"中的电流 1，并检查是否列出了电流守恒 1、电绝缘 1 和初始值 1。然后，右键单击电流 1 > 接地（这将 0 V 分配给远平面，模拟远距离参比电极），并将其应用于离电极最远的表面（表面 10）。
      2. 接下来，右键单击分配给圆盘电极（表面 14）的 电流 1 >浮动电位（模拟恒定电流的电流源），并将 I ₀ 值更改为 1[μA] 以施加单一电流。
   2. 使用参数化扫描运行模拟。
      注意:此步骤将运行模拟并添加参数化扫描，其中在每次模拟中更改参数值的情况下进行了多次模拟。在这里，扫描电极半径参数，并将每个仿真的电势分布存储在模型文件中。运行仿真后，模型树中的结果分支填充了电势 （ec） 多层图。
      1. 在"模型树"中，右键单击" 参数化扫描">"算例 1"。单击" 参数化扫描"，在 "研究设置 "表中，单击" 添加"，然后选择" elec_rad "作为 "参数名称"。
      2. 为参数 值列表 键入"50、150、350、500"，为 参数单元键入"μm"，然后单击" 计算 "以运行研究（补充图4）。

图 1:创建几何图形。 将块几何图形插入到 FEM 模型中以表示组织。 请点击此处查看此图的大图。

图 2:创建电极的几何形状 。 （A） 制作工作平面以绘制圆盘电极。（B）在工作平面上绘制圆以创建圆盘电极。 请点击此处查看此图的大图。

图 3:有限元模型的单元质量直方图。直方图显示了整个模型中元素的质量。如果很大一部分元素位于低质量区域，则需要网格细化。请点击此处查看此图的大图。

图 4:为电极分配电流值。 在 FEM 软件中施加到电极几何形状的单一电流。 请点击此处查看此图的大图。

2. 在神经元计算套件的 GUI 中导入神经细胞的几何形状

1. 构建单元模型的几何图形
   1. 使用单元格生成器功能导入形态。
      1. 从神经元计算套件的安装文件夹运行"nrngui"，单击" 工具">"杂项">"导入 3D"，然后勾选" 选择文件 "框。
      2. 找到下载的 .swc 文件，然后单击 读取。导入几何图形后，单击导出 >单元格生成器 （补充图 5）。
   2. 创建导入的细胞形态的 .hoc 文件
      1. 转到"子集"选项卡并观察模型中预定义的 子集 （例如，体细胞、轴突、基底等）。勾选 "连续创建 "框，转到 "管理>导出"，然后将形态导出为"rgc.hoc"。
   3. 查看细胞形态
      1. 单击工具栏上的 1 个实像元>>模型视图>根 Soma[0] ，右键单击出现的窗口，然后单击轴 类型>视图轴。通过目视检查，该模型的树枝状场直径约为250μm。暂时关闭神经元窗口。

图 5:将神经元模型信息导出为 .hoc 文件。 神经元的几何形状被导出到一个.hoc文件，以允许进一步的修改。 请点击此处查看此图的大图。

图 6:测量神经元的尺寸。 神经元的形态（顶视图）显示在神经元计算套件的GUI中，并叠加了x-y轴。刻度以微米为单位。 请点击此处查看此图的大图。

3. 编程神经元计算模拟

1. 通过用 .hoc 语言编程来调整细胞的形态
   注意:可以通过GUI的单元格生成器功能 调整 单元格的形态。但是，演示了如何通过编辑 .hoc 文件来加快该过程来完成此操作。.hoc文件定义了建模细胞的拓扑（神经元每个部分之间的物理连接），形态（每个神经元部分的长度，直径和位置）和生物物理特性（离子通道参数）。有关 .hoc 编程的完整文档，请参见:https://neuron.yale.edu/neuron/static/new_doc/index.html#，
   1. 使用文本编辑器（例如记事本）打开生成的 .hoc 文件。如Sheasby和Fohlmeister 17中所述，添加40μm长度的轴突初始段和90μm长度的窄轴突段，并改变树突的长度，使树枝状场的大小变为180μm以匹配Rockhill等人中的G1细胞^18。
      1. 创建新的单元格部分并为每个部分定义拓扑连接。
         1. 要为轴突初始片段 （AIS） 和窄轴突段 （NS） 创建新的单元格部分，请将以下行添加到 rgc.hoc 文件的开头:
            创建 AIS、NS // 声明称为 AIS 和 NS 的单元隔室
            然后，将行"连接轴突（0），soma[1]（1）"替换为:
            连接 AIS（0）， 相细胞[1]（1） // 将 AIS 的第一段连接到体瘤的末端[1]
            连接 NS（0）， AIS（1） // 将 NS 的第一段连接到 AIS 的末端
            连接轴突（0）， ns（1） // 将轴突的第一段连接到 NS 的末端
      2. 定义单元格部分的 3D 位置、直径和长度
         1. 通过在"proc shape3d_31（）"括号内编写这些行来定义 AIS 和 NS 隔间的 3D 位置和直径:
            ais { pt3dadd（-2.25， -1.55， 0， 1） // 前三个数字是 xyz 坐标，直径为 1 μm
            pt3dadd（37.75， -1.55， 0， 1）} // 第一个点在 x = -2.25 μm 处，最后一个点在 x = 37.75 μm 处
            ns { pt3dadd（37.75， -1.55， 0， 0.3） // NS 段的 3D 坐标和直径
            pt3dadd（127.75， -1.55， 0， 0.3）}
         2. 在文件末尾，通过键入以下内容来移动轴突的 3D 坐标，使其初始点与 NS 的终点相遇:
            axon {for i=0，n3d（）-1 {pt3dchange（i， x3d（i）+130， y3d（i），z3d（i）-5， diam3d（i））}} //Shift the x-坐标
         3. 在文件末尾，通过键入以下内容将树突状隔室缩短 18%:
            Forsec baseal {L=L*0.82} // 缩放长度以使树枝状场尺寸更小
            define_shape（） // 填写缺失的 3D 信息
2. 定义每个部分的段数
   注意:神经元的每个部分都可以在空间上离散化，就像 FEM 模型中的网格划分过程一样。空间离散化将神经元实际上分成更小的部分，在那里进行计算。对于段数"nseg"，请确保使用奇数以确保在单元格部分的中点有一个内部节点，并尝试将 nseg 数增加三倍，直到计算产生一致的结果⁹。更多的段将产生更准确的数值近似，但也增加了计算负载。
   1. 为了举例说明离散化过程，请在 rgc.hoc 文件中添加以下行，将体细胞和轴突子集中的神经元部分划分为几个片段:
      Forsec somatic {nseg=31}
      Forsec axonal {nseg=301}
      模型中的其他部分也需要通过键入这些行进行离散化，但根据需要更改子集名称（在"forsec"之后）和段数（在"nseg"之后）。
3. 插入定制的离子通道机制
   1. 将自定义的离子通道机制编写为 .mod 文件:要应用离子通道机制，请创建 .mod 文件并按照以下步骤将文件插入 .hoc 文件的生物物理部分部分 3.3.1-3.3.3..mod文件包含每个离子通道要求解的变量和微分方程。
      注意:正确的离子通道机制定义和实现对于准确的神经元模拟至关重要。写入 .mod 文件时，请检查单位是否正确（提供的"modlunit"实用程序，可以运行以检查单位一致性）以及方程是否正确编写。为了测试离子通道机制是否正确，可以绘制细胞内或细胞外刺激期间每个离子通道的电流并与经验发现进行比较。
      1. 电压门控离子通道
         注意:用于创建电压门控离子通道的.mod文件通常包括一个 DERIVATIVE 模块，其中包含要求解的微分方程，一个 BREAKPOINT 模块，其中包含使用所选数值近似方法求解微分方程的命令，以及告诉程序计算门控参数（例如，本例中的mt，ht和d）的 PROCEDURE 模块。此代码将计算每个时间步长通过通道的离子电流值。
         1. 为了举例说明该过程，请创建一个电压相关的Ca通道，该通道具有一阶微分方程来求解门控变量。
         2. 在文本编辑器中打开一个新文件，然后在 补充材料 - 定义电压相关的 Cat 通道中键入行。将此文件保存为 Cat.mod 与 .hoc 文件位于同一文件夹中。需要对模型神经元包含的其他离子通道重复此过程。
      2. 电压和浓度依赖性离子通道
         1. 某些离子通道的动力学，例如视网膜神经节细胞中的钙激活钾通道，除了跨膜电压¹⁹外，还取决于细胞内钙浓度。要对这种类型的机制进行建模，请创建一个名为KCa.mod的文件，并键入补充 材料-电压和浓度依赖性离子通道中所示的线。在这个.mod文件中，计算了变量"cai"，它被定义为Ca离子的内部浓度，然后在方程中使用该变量来计算iKCa电流。
   2. 编译 .mod 文件
      1. 通过从安装文件夹运行神经元计算套件的 mknrndll 实用程序来编译所有 .mod 文件。找到包含 .mod 文件的文件夹，然后单击编译以创建 O 和 C 文件。在此之后，可以将机制插入到该细胞模型中。
   3. 在主神经元模型文件中应用 .mod 文件。
      注意:除了插入离子通道外，最大Na电导仅针对"体细胞"子集定义。如有必要，我们可以单独调整不同神经元段的最大膜电导率。
      1. 为简洁起见，将所有离子通道机制合并到一个 .mod 文件（补充材料-完整 .mod 文件）中。通过在rgc.hoc文件的"biophys"程序中键入以下行，将包含所有离子通道和被动泄漏通道的组合.mod文件插入"体细胞"子集中的所有片段:
         forsec somatic {insert rgcSpike
         插入 PAS // 被动泄漏通道
         gnabar_rgcSpike = 80e-3
         g_pas = 0.008e-3 // 泄漏膜电导率}
   4. 设置轴质电阻率
      1. 细胞具有轴质电阻率，可以按隔室改变。对于该模型，所有段都具有相同的电阻率，为110 Ω·cm。更改 rgc.hoc 文件中的轴质电阻率:
         forall {Ra = 110}
4. 插入细胞外机制并定义脉冲波形
   1. 将细胞外机制插入细胞模型中
      1. 为了使细胞模型响应细胞外电压，请通过键入 rgc.hoc 文件底部的行将细胞外机制插入所有片段中:
         全部{插入细胞外}
   2. 产生双相脉冲
      注意:在本演示中，通过在 .hoc 文件中创建一个过程，用户可以在脉冲宽度、相间间隙和重复次数方面进行调节的恒流双相脉冲。对于更结构化的程序，请使用 rgc.hoc 文件作为文件来创建细胞模型，而刺激过程应用于单独的 .hoc 文件，该文件将加载应用刺激的细胞模型。
      1. 创建一个名为 stimulation.hoc 的新文本文件，并通过加载单元格模型文件来启动代码;然后，通过定义补充 材料 - 在神经元模拟中创建双相脉冲中所示的程序来制作双相脉冲。
         注意:此步骤创建一个恒流阴极优先双相脉冲，其中激励参数由用户在运行仿真时声明。目前，阳极脉冲和阴极脉冲的幅度为±1μA，但该幅度需要根据圆盘电极提供的刺激电流而变化。

4. 运行和自动化多个模拟

1. 组合模型
   1. 提取神经元细胞模型中节点的坐标
      注意:组合模拟的目的是获取与细胞模型的每个节点相对应的细胞外电位值。但是，两个模型的坐标必须对齐。在本例中，体细胞的中心段（soma（0.5））对齐，位于视网膜组织的水平中平面（对应于视网膜神经节细胞层），体细胞的中心节点位于圆盘电极的中心正上方。
      1. 打开 FEM 模型，并注意参考点的坐标（例如，视网膜组织的水平中平面，在圆盘电极中心上方），在这种情况下，它是 [0， 0， 131.5] μm。
      2. 在神经元计算套件中，创建一个名为 calculateCoord.hoc 的文件来提取每个段的质心的坐标并移动每个截面，以便 soma 的中心段与 FEM 模型中的参考点具有相同的坐标（补充材料 - 计算每个节点的坐标）。
   2. 将坐标点保存到文本文件
      1. 运行 calculateCoord.hoc 文件（通过双击文件资源管理器或打开神经元计算套件的 GUI;然后单击工具栏中的 文件>加载 hoc ）。将要评估的细胞外电压值的坐标保存到名为"坐标.dat"的文本文件中。
   3. 运行仿真并将电压数据保存到文本文件中
      注意:在此步骤中，我们从FEM模型中提取计算出的细胞外值，但我们只会保存与每个细胞段中心重合的相关坐标中的数据。当需要大量电位导出时，请遵循步骤 4.1.6.2。
      1. 在有限元软件中打开组织模型文件;转到"模型树"中的"结果"标题，然后单击">数据导出>数据 1"。确保数据集设置为研究 1/参数解 1，然后在表达式列中键入"V"，在单位列中键入"mV"。
      2. 在"输出"下，将"文件名"更改为"单元外.dat"，然后选择"要评估的点:从文件"。加载坐标.dat坐标文件字段，然后单击导出。
   4. 将双相脉冲应用于细胞模型
      注意:在此阶段，每个细胞段在一个时间点（电流为1μA）的细胞外电压值可用。由于该研究打算使细胞受到双相脉冲，因此使用"vector.play"方法使细胞经历的细胞外电压值随时间变化。
      1. 添加补充 材料-双相脉冲在刺激中的应用.hoc中显示的线条。
   5. 运行组合模拟
      注意:需要定义数值近似值的时间间隔"dt"才能运行模拟。与nseg类似，较短的dt可以提高计算精度，但也会增加计算成本。
      1. 将 补充材料 - 执行神经元模拟中显示的行添加到刺激.hoc的末尾。然后，双击 stimulation.hoc 文件以加载脚本并自动运行模拟。感兴趣片段的跨膜电位可以显示在神经元计算套件的GUI中（步骤4.2.1）或保存到文本文件中以便在其他程序中读取（步骤4.1.6.1.2）。如果需要重复计算和大量膜电位，请执行步骤4.1.6.1和4.1.6.2。
   6. 额外:自动化模拟
      注意:要找到阈值幅度，请每次使用不同的电流幅度多次循环仿真。可能需要另一种自动化来找到相对于刺激电极位于不同位置的神经元的阈值。自动化步骤可以在神经元计算套件中使用程序完成，也可以在 FEM 软件中使用称为"方法"的脚本完成。
      1. 自动执行神经元模拟以查找阈值幅度
         注意:一批神经元模拟可以自动完成。在神经元模拟程序中实现以下步骤，以查找不同刺激参数下神经元的阈值振幅。
         1. 创建一个过程以在神经元模拟程序中重复模拟:在 stimulation.hoc 中，创建一个矢量以包含要测试的电流幅度范围。然后，创建一个过程来应用电流幅度并记录任何尖峰的存在（从负跨膜电压到正跨膜电压的正变化），阈值幅度定义为导致发生尖峰的最低电流幅度。为此，请在 stimulation.hoc 文件末尾定义一个名为 findTh（）（在一系列电流振幅上循环的补充材料）的过程
         2. 将阈值处的响应保存到文本文件:将以下行添加到 stimulation.hoc 中的 findTh（） 过程中，以将每个时间步的所有神经元区室的跨膜电压值存储在文本文件中:
            sprint（save文件名， "Response_%d.dat"， th） // 存储阈值
            saveFile.wopen（save文件名）
            for i=0，（responseVector.size（）-1）{
            saveFile.printf（"%g， "， responseVector.x[i]）
            if（i==responseVector.size（）-1） {saveFile.printf（"%g\n"， responseVector.x[i]）
            saveFile.close（save文件名）
            }}
      2. FEM 软件中的自动化功能可查找不同位置神经元的电压值
         注意:另一个可以完成的自动化是自动获取不同位置神经元的细胞外电压值。FEM 软件中的"应用程序构建器"菜单提供了一种定义"方法"的方法，或一种脚本来自动执行软件执行计算所需的步骤。为了证明，细胞在x方向上的位置以100μm的步长移动了5倍（补充图6）。
         1. 编写用于 FEM 模拟自动化的代码。
            1. 转到应用程序构建器，右键单击应用程序构建器树中的方法，选择新建方法，然后单击确定。转到文件>首选项>方法，勾选查看所有代码框，然后单击确定。
            2. 编写一个 .hoc 脚本，该脚本将加载坐标文件，移动值以匹配所需位置，并通过键入 补充材料 - 定义自动化 FEM 仿真方法中显示的代码来保存包含电池新位置的电压值的文本文件。
         2. 在 FEM 软件中运行自动化步骤:切换到 模型生成器、开发人员>运行方法>方法 1。这将生成具有适当电压值的.dat文件，命名为 extracellular_1.dat、extracellular_2.dat 等。
      3. 在通用编程语言中循环模拟
         注意:要循环模拟，每次都需要将适当的文本文件加载到神经元计算套件的模拟中，并且可以轻松加载和操作文本文件的编程语言²⁰ 可以方便地执行此步骤。任何方便的集成开发环境 （IDE）²¹ 都可用于此步骤。
         1. 打开选定的 IDE，单击" 新建文件 "以创建新脚本。此示例中使用了一个.py文件。键入 补充材料 - 使用通用编程语言运行模拟中显示的行。
         2. 最后，单击"运行"或按 F5 运行脚本，这也将打开 GUI（补充图 7）。
2. 显示模拟数据
   注意:按照上述所有步骤，仿真结果应存储在文本文件中，其中包含阈值和阈值处的跨膜电位。但是，用户可以选择在使用 NEURON 的 GUI 运行模拟时显示模拟结果。
   1. 在神经元计算套件的 GUI 中绘制神经元模型对细胞外刺激的反应。为此，请运行 stimulation.hoc，单击工具栏中的">电压轴图形"，然后在图形窗口中右键单击任意位置并选择"绘制内容"。
   2. 在" 变量到图 "字段中输入"axon.v（1）"，这意味着它将绘制每个时间步长轴突最后一段的跨膜电位。

图 7:显示 FEM 计算结果并将其导出到文本文件。 图形窗口显示电势的多切片图（以 V 为单位）。"数据导出设置"中的选项允许将计算变量导出到文本文件中。 请点击此处查看此图的大图。

图 8:使用电压图显示跨膜电位图。神经元跨膜电位显示在神经元计算套件的GUI中。x 轴是以毫秒为单位的时间，而 y 轴是以 mV 为单位的所选神经元段的跨膜电位。请点击此处查看此图的大图。

我们进行了两种仿真协议来演示该模型的使用。第一个协议涉及改变电极尺寸，同时保持神经元的位置和电脉冲参数相同。第二个方案涉及以100μm步长在x方向上移动神经元，而电极的大小保持不变。对于两种方案，使用的脉冲都是宽度为0.25 ms的单阴极优先双相脉冲，相间间隙为0.05 ms。对于第一个方案，电极的半径变化为50，150，350和500μm，而对于第二个方案，电极的半径保持在恒定的50μm。

这里描述的模型表明，在0.25 ms脉冲宽度下增加脉络膜上电极尺寸会增加模型神经元的激活阈值（图9A）。该结果反映了Liang等人的 体内 发现¹²，他们表明皮质激活阈值随着该脉冲宽度下电极尺寸的增加而增加。

由于几个因素，模型激活阈值的大小与经验发现不同。首先，该模型仅涉及特定类型的单个RGC，该RGC可能不存在于 体内 研究中被激活的细胞组中。接下来，该模型不包括视网膜网络，该网络可能通过双极细胞的兴奋性输入促进RGC的激活。造成差异的另一个可能原因是电极-视网膜距离。 由于 解剖学变异性或手术，体内研究中的电极-视网膜距离可能低于该模型。因此，我们高估了电极-视网膜距离，从而高估了激活阈值。同样重要的是要注意，尽管在我们的结果中没有证明这一点，但对单细胞阈值进行建模通常会低估 体内 皮质阈值。这是由于皮质测量的技术限制（主要与信噪比有关），皮质活动通常仅在多个视网膜神经节细胞被激活后才能检测到。因此，视网膜和皮质激活阈值的大小存在差异是可以预期的。尽管存在这些差异，但由于电极尺寸的增加，该模型成功地显示了活化阈值的增加趋势。这是由于当电极尺寸增加时，与周围环境相比没有高电场区域，这不利于神经激活²²。

接下来，我们观察动作电位特征以验证此处描述的模型。潜伏期或刺激开始与动作电位峰值峰值之间的时间范围为1-2.2毫秒（图9B）。这对应于由于非网络介导的视网膜激活引起的短延迟峰值²³。该模型的尖峰宽度为1 ms，这与 体外²⁴中测量的兔RGC的尖峰宽度处于同一范围内。

在第二个刺激方案中，只有神经元在x轴上（沿着轴突的长度）相对于电极的位置是变化的。在距离 0 处，体细胞截面的质心紧邻圆盘电极的中心。负距离意味着圆盘电极的位置更靠近轴突侧，而正距离意味着圆盘电极的位置更靠近树突侧。该模型显示，当轴突的窄段紧邻圆盘电极上方时，达到了最低阈值，并且随着x距离变大而增加（图9C）。与将电极向树突移动相比，将电极进一步向远端轴突移动会产生较低的阈值，因为存在轴突初始段和钠通道更普遍的狭窄段。这一结果与Jensen等人的 体外 发现一致¹³，其中兔RGCs被超细微电极刺激，当电极靠近树突时，活化阈值最高。

图 9:建模方法的结果 。 （A）位于圆盘电极上方的视网膜神经节细胞的活化阈值。电极半径变化（50、150、350和500 μm），阈值随电极尺寸的增加而增加。（B）神经元模型在0.25 ms脉冲宽度下的动作电位形状。不同电极尺寸的阈值动作电位具有相同的1 ms尖峰宽度，但潜伏期随着电极尺寸的增加而增加。刺激开始时间为1 ms，阴极相在膜处引起去极化，但不足以引起动作电位。（C）神经元沿x轴移动，激活阈值表明，其体细胞位于电极中心正上方的神经元达到了最低阈值。电极的半径为50μm。 请点击此处查看此图的大图。

补充图1:有限元模型的初始化。 "研究"和"物理"类型决定了模型中求解的方程列表。这些是在初始创建 FEM 模型文件期间设置的，但也可以在创建模型后进行修改/添加。 请点击此处下载此文件。

补充图 2:更改长度单位。 长度单位和角度单位确定几何定义过程中使用的单位。 请点击此处下载此文件。

补充图 3:插入材料属性。 材料属性是为 3D 模型中的每个域定义的。可用的材料属性列在"材料设置"窗口的"材料属性"中。对于电势计算，仅定义了电导率属性。 请点击此处下载此文件。

补充图 4:创建参数算例以循环参数值列表。 参数研究允许 FEM 软件自动重复计算并更改每次重复的电极半径值。存储每次重复的计算结果。 请点击此处下载此文件。

补充图 5:从 SWC 文件导入神经元形态。 神经元计算套件能够读取从神经元追踪中获得的SWC文件。导入的文件包含有关每个神经元段的形态和拓扑的信息。 请点击此处下载此文件。

补充图 6:通过定义方法自动执行 FEM 操作。 通过编写脚本来定义一种方法，以在 FEM 软件中自动执行流程，这是通过定义参数算例无法完成的。 请点击此处下载此文件。

补充图 7:使用通用编程语言集成模型并自动执行仿真。 通用编程语言用于循环神经元模拟，同时更改用作输入的细胞外电压文件和用作循环中每个步骤的输出的神经反应电压文件。 请点击此处下载此文件。

补充材料:用于 （1） 定义电压相关 Cat 通道的命令行。（2） 电压和浓度相关的离子通道。（3） 完成 .mod 文件。（4）在神经元模拟中创建双相脉冲。（5）计算每个节点的坐标。（6）双相脉冲的应用。（7）执行神经元模拟。（8） 在一定范围的电流幅度上循环。（9） 定义一种自动化有限元模拟的方法。（10） 使用通用编程语言运行模拟。 请点击此处下载此文件。

在本文中，我们展示了一种结合了有限元和生物物理神经元建模的建模工作流程。该模型非常灵活，因为它的复杂性可以修改以适应不同的目的，并且它提供了一种根据经验发现验证结果的方法。我们还演示了如何参数化模型以实现自动化。

两步建模方法结合了使用FEM和神经元计算套件的优势，在细胞外刺激存在的情况下求解神经元的电缆方程。FEM 可用于精确计算体积导体上的细胞外场，在复杂几何形状或电导率不均匀的情况下，解析求解通常是不切实际的。该模型的计算成本也相对较低，因为假设是静态条件。

虽然所描述的建模方法在易用性和灵活性方面具有优势，但这种建模工作流程存在局限性。首先，这种方法不允许在计算电场时存在神经膜。Joucla等人²⁵ 将两步法与整个FEM方法进行了比较，其中神经几何和膜特性包含在FEM模型中。他们表明，当几何结构中包含较大的细胞结构（例如细胞体）时，将神经元包括在电场计算中将改变跨膜电位计算。具体来说，两步法中神经元几何形状的简化意味着隔室中任何点的跨膜电位由隔室节点或中心点处的跨膜电位表示。相比之下，Joucla提出的全FEM模型包括神经元3D几何形状的显式表示，这使得能够单独评估隔室内任何点的跨膜电位。因此，如果需要跨膜电位的确切形状和位置，全有限元模型可能更合适。但是，此方法在计算上比两步法更昂贵。

建模方法的第二个限制涉及形态学和离子动力学数据的可用性。这里使用的模型基于虎蝾螈数据，该数据已用于模拟其他物种的RGC，但存在的离子通道类型可能存在差异，尚未阐明。因此，在某些情况下可能需要进行 体外 工作以调整离子通道参数。

第三，FEM软件的成本可能是一个制约因素。在这种情况下，具有内置泊松方程求解器的开源 FEM 程序²⁶ 可能是一种替代方法。除了使用的 FEM 软件外，此工作流程中使用的软件是免费的。虽然所使用的FEM软件提供了直观的GUI和即用型电流建模，但可以在通用编程软件中执行细胞外值计算。然而，这将需要手动定义物理方程和数值方法来解决方程²⁷。此外，当使用复杂的组织或电极阵列几何形状时，这种方法可能很乏味。

作者声明没有竞争利益。

这项研究由国家卫生和医学研究委员会项目资助（资助号1109056）。