数据驱动的爆发形态分析用于神经元培养体系的功能表型表征

《iScience》：Data-driven burst shape analysis for functional phenotyping of neuronal cultures

　　

**亮点**

• 改进了原代/分化神经元培养物中群体群发活动的描述
• 用于量化群发时间特征曲线（群发形状）的数据驱动框架
• 结合群发形状和传统特征持续改善表型可分性

**摘要**

体外神经元培养物是研究正常和病理条件下网络动力学的有力工具。成熟网络通常表现出网络群发活动，传统上通过简化特征（如群发间隔和群发持续时间）进行量化。这些特征增进了我们对发育、疾病表型和药物效应的理解。然而，它们忽略了群发内部的时间结构，而这一结构对网络变化高度敏感，因此可以揭示额外的生理或病理效应。在此，我们开发了一个全面的框架来量化群发形状——即网络在群发期间的放电时间进程。在四个数据集上，包括啮齿动物和人类多能干细胞衍生培养物，我们表明群发形状包含关于底层网络动力学的丰富信息。我们通过使用传统和群发形状特征对条件（遗传性疾病、药物）进行分类来量化这些信息，发现加入形状特征显著提高了准确率。我们提供了一个开放的群发形状特征化流程，引入群发形状作为互补的功能表型特征，扩展了疾病建模和药物筛选的分析工具包。

**图解摘要**

**引言**

同步网络群发是体内和体外发育神经网络的基本活动模式，源于兴奋性、连接性和动态神经元特性的相互作用。

神经细胞培养物通常发展出此类群发活动，通常表现为持续数十毫秒到数秒的同步放电，随后是长而不规则的静默间隔。这种活动模式在体外皮层、海马、纹状体和脊髓培养物中以及人类诱导多能干细胞（hiPSC）衍生神经元培养物中均有观察到。

体外群发动力学统计量偏离健康对照是神经发育障碍的强有力指标。

传统上，神经元培养物中的群体活动通过群发事件的汇总统计量进行表征，如持续时间、群发间隔和平均放电率。

随着监控网络活动工具的不断涌现，这些度量已成为评估疾病相关基因变异、以及单个神经元和网络慢性及急性扰动效应的常用工具。

然而，群发形状——即群发内部放电率的详细时间进程——尽管理论、模拟和实验均指出网络特性与群发内部动力学之间存在密切联系，但受到的关注相对较少。

先前的工作尝试了对群发形状进行参数化描述，针对个别实验进行了定制，包括上升、峰值和衰减时间等特征。虽然信息量丰富，但这些简单度量在两个关键方面存在局限性。首先，它们过于简化了复杂的活动模式，如波动的放电率、反跳放电和振荡。其次，将群发形状总结为每个记录的几个数字会掩盖每种培养物内群发形状的变异性，可能掩盖具有生物学意义的异质性。

为了克服这些局限性，我们引入了一种数据驱动框架，以标准化的方式对群发形状进行标准化和比较。在四个已发表的啮齿动物和hiPSC衍生神经培养物数据集上，该框架表明群发形状信息显著改善了实验条件的分类，揭示了群发内部动力学与底层网络特性之间的密切联系。我们的方法能够对群发形状进行可视化，并提取紧凑的、可解释的描述符，其性能接近完整数据驱动表示的水平。此外，我们发现遗传相似的培养物表现出相似的群发形状模式，证明群发形状携带了关于网络表型的互补、生物学相关信息。

**结果**

我们的目标是最大化从不同实验条件下（包括药物操作和疾病相关突变）获取的多电极阵列（MEA）记录中提取的信息（图1A）。传统上，群体活动通过一组简单的统计特征进行描述，如群发持续时间、群发间隔和放电率（图1B）。由于这些度量忽略了群发内部的时间结构，我们提出用群发形状特征对其进行增强。下文将描述如何表示和比较群发形状，以及这些表示如何被纳入训练用于推断实验条件的统计模型。在该框架内，分类准确率的改善量化了群发形状所携带的额外信息（图1C）。

图1 评估群发形状作为网络表型中传统特征的互补方法的示意图

显示完整图注 图查看器

(A) 神经细胞培养物，包括啮齿动物和hiPSC衍生细胞，铺板并通过多电极阵列（MEA）记录。
(B) 群发活动的传统描述符（群发间隔、群发持续时间、随机脉冲百分比、总放电率和群发期间放电率）由群发形状特征进行补充。
(C) 使用推断实验条件的模型评估传统和群发形状特征的贡献。

**数据驱动的群发形状量化框架**

我们分析了四个不同实验条件下MEA记录的网络活动数据集（数据集简称见方括号，详见STAR方法）：使用GABAA受体拮抗剂比库利阻断抑制的鼠皮层神经元（阻断抑制）、具有CACNA1A单倍体不足的hiPSC衍生神经元（CACNA1A+/?）、克勒夫斯特拉综合征患者和对照的hiPSC衍生神经元（Kleefstra综合征）、以及铺板批次的大鼠皮层神经元（发育中大鼠）。我们在970次记录和604小时活动中检测到31,949次网络群发（表S1）。为确保在此规模下的可靠群发检测，我们开发了一种交互工具，用于详细审查和参数调优（图S1，检测参数见表S2，替代群发检测见表S4）。

为了比较持续时间和网络放电率范围广泛的群发形状（从数十毫秒到数秒），我们将所有群发形状标准化为单位持续时间和单位积分（图2A）。平均这些标准化形状可提供群体差异的可解释曲线，从而能够设计手动群发形状特征（图2B）。我们通过构建K最近邻（KNN）图来捕获个体群发之间的相似性，将每个群发形状连接到其K个最相似的群发形状，使用沃瑟斯坦距离（也称为搬土距离）进行测量（图2C，见STAR方法）。从该图中，我们通过汇总个体群发的邻居来获得记录级别的预测（图2D）。然后我们对KNN图应用谱嵌入，在保留其结构的同时将群发形状投影到低维空间（图2E），可以进行交互式探索（图S2）。跨数据集的低维嵌入与传统的群发特征存在相关性，如群发偏度、峰值时间和群发宽度。然而，这些相关性的强度和方向在不同数据集间是变化的（图S14和S15）。将记录中所有群发的二维表示取平均值得到数据驱动的形状汇总特征（图2E）。最后，我们结合这些形状特征与传统特征进行XGBoost分类（图2F），比较不同分类器的性能。这一群发形状的表示为基础评估其对实验条件分类的贡献奠定了基础，详见以下章节。

图2 数据驱动的群发形状量化框架

显示完整图注 图查看器

(A) 通过标准化持续时间和大小计算标准化群发形状。
(B) 通过跨记录平均来理解群体差异。
(C) KNN图将每个节点（群发形状）连接到K个最相似的节点。
(D) 基于个体群发形状的邻居聚合到记录级别进行分类。
(E) 个体群发形状的谱嵌入聚合（平均）到记录级别，其中形状维度定义形状特征。
(F) 结合传统和形状汇总特征使用XGBoost预测记录标识。

**用群发形状特征增强传统特征可改善跨数据集的分类性能**

为了测试哪些特征包含关于培养物表型的有效信息，我们训练了四个分类器：基于个体群发形状的KNN分类（图2D）以及基于三组特征的XGBoost分类（传统、群发形状和组合）（图2F），其中所有模型的表现在随机水平以上（图3A）。我们仔细比较了交叉验证的模型性能与频率主义方法（见STAR方法），并在表S3和图S6中报告了贝叶斯比较。总体而言，包含群发形状信息改善了分类器性能，在所研究的数据集内提供了更具区分力的培养物表型描述。

图3 用群发形状特征增强传统特征可改善跨数据集的分类性能

显示完整图注 图查看器

(A) 四个模型的平均分类准确率：基于个体群发形状的KNN聚类，以及基于三个汇总特征集（传统、形状-2D和组合）的XGBoost模型。数据表示为均值±均值的标准误（n=100次随机分割，80-20比例）；y轴最小值表示随机水平。
(B–D) 不同模型在轴标签所示数据集上的准确率差异。组合特征集优于传统集（B），个体群发特征优于形状汇总特征（C），在某些情况下优于传统特征（D）。数据表示为均值±均值的标准误（n=100次随机分割，80-20比例）。星号表示Wald检验的显著性（*p<0.05，**p<0.01，***p<0.001）；荟萃分析（假设声明在每个图像上方）使用随机效应模型。
(E) 七个组合特征的Shapley值分析表明，两个形状特征对预测做出了重大贡献。
(F) 阻断抑制数据集中4×6孔板布局的KNN分类预测（留一法）。饼图份额代表记录的加权邻居；实线边框表示真实类别。
(G) CACNA1A数据集的KNN分类预测。虚线分隔真实标签（上：对照，下：CACNA1A+/?）。

组合特征集（群发形状和传统描述符）比仅使用传统特征更有信息量，如分类器准确率的显著改善所示（荟萃分析，p=3.4×10?3，图3B；对于沃瑟斯坦指标最强，见STAR方法；图S10）。为了确定七个组合特征中哪些对这一改善贡献最大，我们计算了Shapley值（图S7）。群发形状特征贡献显著，在阻断抑制和CACNA1A数据集中占重要性的50%以上（图3E），在其他两个数据集中约占30%。

基于KNN图的分类器直接利用了不同记录中个体群发之间的相似性。在某些药物和疾病数据集中，这种相似性提供了良好的培养物表型表征，在阻断抑制（98%）和CACNA1A数据（100%）中实现了接近完美的准确率，并显示出高模型置信度（图3F和3G，其他数据集见图S4和S5）。与基于汇总群发形状特征的表示相比，个体群发形状的相似性包含的关于实验条件的信息显著更多（图3C，p=2.3×10?2）。

综上所述，两种方法解决了不同的分析目标。基于KNN的分类器最适合于同一数据集内已知条件之间的标签迁移和分类，因为它直接利用了个体群发之间的相似性。相比之下，群发形状嵌入提供了一种灵活的表示，可以与传统特征结合，实现对培养物表型更丰富的表征，并支持跨数据集的探索性和比较性分析。

**平均标准化群发形状揭示改变的起始和终止动态**

虽然数据驱动的谱嵌入形状特征紧凑且在分类中表现良好，但其抽象维度不易于解释，仅为研究潜在机制提供有限指导。为此，我们设计了基于平均群发形状的特征，这些特征捕获了起始和终止动态的生物学有意义方面，同时保持了接近数据驱动特征预测性能。

从原始（未标准化）网络群发的代表性示例中可以看到组间群发形状的差异（图4A）。然而，原始群发形状的平均值（图4B和S3A–S3C）被大型、长时间群发所主导，掩盖了形状特征。相比之下，标准化群发形状的平均值具有非常紧凑的组均值（图4C，阴影表示±3×均值的标准误）。这凸显了将群发形状标准化为单位持续时间和单位积分作为分析关键步骤的重要性。

改变的起始和结束动态是各组条件之间最显著的形状差异，我们通过两个手工设计的形状特征来捕捉它们（图4C中的虚线）。

**图4 平均归一化爆发形状揭示了改变的起始和结束动态**

(A) 原始（未归一化）爆发形状的代表性示例：叠加的每10毫秒 bin 的突触计数和光点图；上行显示对照条件，颜色编码与(C)相同。

(B) 每次记录原始（未归一化）爆发形状的平均值受大型爆发主导，掩盖了归一化后显而易见的特征（以阻断抑制的培养物为例）。

(C) 对照培养物以及具有阻断抑制、CACNA1A单倍体不足和Kleefstra综合征的培养物中归一化爆发的平均值。

(D) 平均峰值时间（C中的彩色虚线）在所有实验中均显著区别于对照，其中阻断抑制和Kleefstra综合征的峰值更早，而CACNA1A+/?的峰值更晚。

(E) 爆发结束附近的振幅（持续时间的80%，C中的黑色虚线）在阻断抑制条件下显著更高，表明突然终止，而在CACNA1A+/?和Kleefstra综合征中则降低，与延长衰减动态一致。

(F) 两个设计的特征——起始(C)和结束(D)——的信息量几乎与数据驱动的特征相当，但随着使用维度的增加，最终被超越。

在(B)和(C)中，阴影表示±3个平均值的标准误（N = 记录数，N = 48，N = 81，N = 390）。在(D)和(E)中，箱线图表示中位数（线）、四分位距（框）、1.5×四分位距（须）和异常值（点）；星号表示t检验的显著性（?p < 0.05，??p < 0.01，???p < 0.001）。在(F)中，数据以均值±平均值标准误表示（n = 100次随机分割，80–20比例）。

我们的第一个设计形状特征通过峰值活动时间来刻画爆发的起始阶段（图4D，图4C中的彩色虚线），类似于对原始爆发形状的传统方法。阻断抑制培养物中的爆发峰值达到时间显著快于对照（抑制关闭 = 0.06，对照 = 0.21，p = 3.1 × 10??），与网络活动募集加速一致。CACNA1A+/?爆发形状显示峰值显著延迟于对照记录（对照 = 0.17，CACNA1A+/? = 0.26，p = 6.4 × 10??），而Kleefstra综合征平均爆发形状峰值早于对照（对照 = 0.28，Kleefstra = 0.15，p = 7.7 × 10?1?）。这些发现与原始爆发形状绝对峰值时间的分析（图S3）及原始研究结果一致。

我们的第二个设计形状特征通过归一化爆发持续时间80%处的活动水平来量化终止动态（图4E，图4C中的黑色虚线），评估爆发是呈现突然结束的持续活动还是延长衰减。阻断抑制培养物在终止附近的活动水平显著高于对照（t = 0.8，抑制关闭 = 0.018，对照 = 0.012，p = 4.0 × 10?12），表明抑制平滑终止爆发，而抑制缺失导致活动突然下降。CACNA1A+/?培养物在爆发终止前显示活动水平显著降低（t = 0.8，CACNA1A+/? = 0.008，对照 = 0.016，p = 1.1 × 10?1?），Kleefstra综合征显示类似但幅度较小且显著的降低（t = 0.8，Kleefstra = 0.012，对照 = 0.015，p = 3.2 × 10??），表明更连续、延长的衰减动态。

精心设计的起始和结束特征包含与数据驱动的爆发形状维度大致相同数量的培养物表型信息（图4F）。我们两个设计特征包含的信息略多于两个数据驱动维度（图S8A）。然而，增加数据驱动维度的数量可捕获更多信息，在各数据集上超越了设计特征（图S8B和S8C）。这表明手工设计的形状特征可以捕获足够的信息，尤其是在简单数据集中，而数据驱动的形状特征可靠且具有可扩展性。

**遗传相关培养物表现出更大的爆发形状相似性**

接下来，我们探讨了遗传因素是否在塑造爆发动态中起重要作用。为此，我们分析了不同遗传相关程度培养物的记录，使用成对Wasserstein距离测量爆发形状相似性，并将其与随机配对的基线距离进行比较（图5）。

**图5 遗传相关培养物表现出更大的爆发形状相似性**

(A) 发育过程中大鼠原代培养物每次记录的平均爆发形状示例，按窝和培养物分组。框架指示(B)中的比较组。

(B) 爆发间的平均Wasserstein距离显示逐渐增加的相似性：同一窝中的配对比不同窝中的配对更接近，同一窝中同一培养物的配对比不同培养物中的配对更接近，而培养物内逐日变异是所有比较中最小的。

(C) hiPSC疾病数据集中平均爆发形状示例。

(D) 同一受试者的记录具有显著更小的爆发形状距离，而疾病组无显著影响。

数据以均值和95%置信区间（误差线）表示。星号表示t检验的显著性（?p < 0.05，??p < 0.01，???p < 0.001）。

在发育中的大鼠原代培养物中，我们发现爆发形状在不同窝之间和整个发育过程中保持稳定。来自同一窝（批次）的记录表现出显著更小的爆发形状距离，与随机比较相比（双侧t检验，p = 3.0 × 10?12，图5B）。在所有记录日期间测量到的单个培养物内距离甚至更小（p = 3.1 × 10??），而培养物内日与日之间的变异是所有比较中最小的（p = 7.1 × 10?1?），且随着记录间隔增大，距离逐渐增长（图S12）。这些发现表明每个培养物都发展出特征性的爆发形状，且随时间缓慢变化。强烈的窝内相似性暗示了遗传成分，尽管接种条件等其他因素也可能有所贡献。这种纵向分析仅适用于Wagenaar数据集，因为其他三个数据集每孔仅有一次记录。

在hiPSC衍生的培养物中，我们发现爆发形状强烈受受试者特异性生物身份影响。来自同一个体衍生的培养物的记录显示显著更小的距离，与随机比较相比（p = 4.3 × 10?1?，图5D）。有趣的是，疾病状态对平均距离的影响较小，尽管仍然可以从KNN图和手工设计的爆发形状特征中读取出来（图3A，图4D和4E）。总体而言，我们的发现表明爆发形状在发育过程中保持稳定，并反映了培养神经元所携带的共享生物背景，包括遗传因素以及共享的接种和处理。

**讨论**

神经元培养物和体内发育过程中的网络爆发呈现出各种各样的轮廓和持续时间，使得系统性地量化其爆发形状变得困难。在本工作中，我们提出了一种数据驱动的方法来解决这个问题，产生了一个提取有意义且实用的爆发形状描述框架。我们将爆发归一化到共同的持续时间和大小，使得可变爆发形状能够汇聚，并使用基于距离的度量来捕捉其轮廓中的信息。这种方法通过利用爆发形状信息增强了网络的表型分析，无需提取特定描述符。因此，爆发形状特征包含有意义的信息以区分实验条件，如病理性基因变异或对神经元活性的药理学扰动。此外，我们展示了爆发形状也反映了培养神经元所携带的共享生物背景。

我们在药理学和遗传学扰动下观察到爆发形状的明显变化，最显著的是爆发的起始和结束。我们将起始解释为活动在整个网络中募集的速度，这可能遵循预设的发放序列，受兴奋性、连接性和突触传递的影响。相反，结束反映了终止爆发的过程，包括兴奋-抑制平衡、短期突触抑制和发放率适应。基于此，我们假设兴奋性增加导致在阻断抑制时起始的加速募集，与既往关于抑制阻断的研究一致。然而，我们还揭示了结束动态中的突发终止，这比依赖于绝对平均值从而掩盖爆发形状的替代方法观察得更清楚。这可能是因为网络缺乏平滑终止爆发的抑制。相反，爆发可能通过兴奋耗竭或积累适应而崩溃，此前已有机制建模。CACNA1A单倍体不足导致单神经元兴奋性升高但突触传递降低。虽然增加内在兴奋性预计会加速募集，但我们观察到的较慢起始和更平滑的结束表明突触传递降低在网络水平上占主导地位。延迟峰值在原研究报告中有记载，但结束变化是新颖的，其效应量似乎甚至更强。Kleefstra综合征中爆发形状更早的峰值与N-甲基-D-天冬氨酸（NMDA）受体NR1亚基的上调一致。然而，爆发结束的异常仍然不明确，因为Kleefstra综合征作为一种发育性网络疾病很复杂，可能导致单细胞和网络水平的可塑性适应。

尽管爆发在发育网络中稳健涌现，但个体培养物中爆发行为的变异性已知较大。此外，在刺激下，培养物可能表现出个体特异性的活动模式。Wagenaar等人表明，这部分变异源于遗传差异。具体而言，他们证明了"姐妹培养物"——来自同一窝的培养物——与非姐妹培养物相比具有较低的变异性。在我们的工作中，我们重新分析了这些数据，发现来自同一窝的培养物其爆发形状也更为相似。扩展这些发现，当我们分析来自同一细胞系的hiPSC衍生培养物的爆发时，我们发现具有共享生物背景（包括相同基因组）的培养物表现出更相似的爆发形状。这表明个体遗传差异与网络爆发形状变异性之间存在联系，尽管这种关联的精确机制仍有待确定。

**研究的局限性**

我们的方法扩展至适用于不同场景有四个主要方向。首先，我们的方法依赖于准确的爆发检测，因为定义爆发起始和终点时的错误可能导致形状错位并掩盖平均值。我们确认了使用替代爆发检测算法时关键洞察得以保留（表S4；图S11），与我们基于网络突触间期特性的简单检测算法相比。然而，更复杂的数据集可能需要更精密的检测或对齐。其次，嵌入选择也会影响框架。虽然谱嵌入在此有效，但替代方法（如图嵌入算法如node2vec）可能捕获额外结构。第三，我们将电极上的突触活动合并为单一时间序列，这便于应用于高通量MEA系统，但忽略了活动的空间传播。纳入二维爆发动态，尤其是高分辨率记录，可能提供关于变异性和网络结构的进一步洞察。最后，虽然我们展示了四个多样化数据集上判别能力的提升，但洞察主要适用于这些数据。需要在更多实验室、MEA平台和筛选队列中进行前瞻验证，以确立疾病建模和药物筛选的广泛适用性。我们未分析性别对结果的影响。

**结论**

我们的框架提供了额外信息，以准确表征发育中神经元培养物的自发活动。爆发形状分析实现了健康和疾病神经元网络更详细的功能表型分析。此外，我们预期纳入爆发形状可以贡献于药理学治疗的评价，以及机制模型的参数推断。

**资源可用性**

**主要联系人**

有关进一步信息和资源的请求应发送并将由主要联系人Tim J. Sch?fer（[e-mail]）接收和回复。

**材料可用性**

本研究未产生新材料。

**数据和代码可用性**

- 本文分析了现有的公开可用数据：
- 解离小鼠皮层培养物（阻断突触抑制），Vinogradov等人，可在Figshare获取：https://doi.org/10.6084/m9.figshare.27110542
- Kleefstra综合征、MELAS和对照hiPSC衍生培养物，Mossink等人，可在Mendeley Data获取：https://doi.org/10.17632/bvt5swtc5h.1
- 大鼠皮层原代培养物（发育时间序列），Wagenaar等人，可在https://potterlab.bme.gatech.edu/development-data/html/daily.spont.dense.text.html获取
- CACNA1A单倍体不足和对照hiPSC衍生培养物数据，Hommersom等人，可按合理请求从原始出版物作者处获取
- 所有原始代码已存入Zenodo，公开可用：https://doi.org/10.5281/zenodo.21918740

该代码库也在 GitHub 上公开可用，地址为 https://github.com/LevinaLab/burst-shape。

如需重新分析本文报告数据所需的任何其他信息，可通过主要联系人请求获取。

致谢

本研究由 Alexander von Humboldt 基金会的 Sofja Kovalevskaja 奖项资助。T.J.S. 通过 Else Kr?ner Medical Scientist Kolleg "ClinBrAIn: Artificial Intelligence for Clinical Brain Research" 获得资助，并通过 Deutsche Forschungsgemeinschaft（DFG，德国研究基金会）的 SFB 1690/1, A07（项目编号 528760423）以及德国卓越战略（EXC 编号 2064/1，项目编号 390727645）获得资助。T.J.S. 和 E.G. 感谢 Intelligent Systems 国际马克斯·普朗克研究学校（IMPRS-IS）的支持，O.V. 感谢心智功能与功能障碍机制国际马克斯·普朗克研究学校和 Joachim Herz 基金会的支持。我们感谢 BMBF 通过 Tübingen AI Center（F.K.Z.: 01IS18039B）提供的支持。A.L. 是 Machine Learning 卓越集群成员，EXC 编号 2064/1——项目编号 390727645。

我们感谢所有实验室成员，特别是 Aaron Spieler 和 Elias Seiffert，感谢他们的有益讨论和支持。我们还感谢 Farzin Negahbani 的有益讨论。特别感谢 Marina Hommersom 和 Nina Doorn 分享他们的数据并提供有益讨论。

作者贡献

概念化：T.J.S.、E.G.、A.L. 和 O.V.；方法论：T.J.S.、E.G.、A.L. 和 O.V.；调查研究：T.J.S.、E.G. 和 P.S.-B.；撰写——初稿：T.J.S.、E.G.、A.L. 和 O.V.；撰写——审阅与编辑：T.J.S.、E.G.、A.L. 和 O.V.；资金获取：A.L.；资源：A.L. 和 O.V.；监督：A.L. 和 O.V.

利益声明

作者声明不存在利益冲突。

生成式AI及AI辅助技术使用声明

在准备本工作期间，作者使用了 Claude、GitHub Copilot 和 ChatGPT 来辅助编码，并使用 Claude、ChatGPT、Gemini 和 Perplexity 来辅助撰写。在使用这些工具或服务后，作者根据需要审查和编辑了内容，并对出版物内容承担全部责任。

STAR★方法

关键资源表

试剂或资源 | 来源 | 标识符

存入数据：MEA 记录，分离的小鼠皮层培养物（阻断突触抑制） | Vinogradov 等人；Figshare | Figshare: https://doi.org/10.6084/m9.figshare.27110542

存入数据：MEA 记录，CACNA1A 单倍剂量不足及对照 hiPSC 衍生培养物 | Hommersom 等人 | 可通过向原始发表作者提出合理请求获取

存入数据：MEA 记录，Kleefstra 综合征、MELAS 及对照 hiPSC 衍生培养物 | Mossink 等人；Mendeley Data | Mendeley Data: https://doi.org/10.17632/bvt5swtc5h.1

存入数据：MEA 记录，原代大鼠皮层培养物（发育时间过程） | Wagenaar 等人 | https://potterlab.bme.gatech.edu/development-data/html/daily.spont.dense.text.html

软件与算法：本研究的自定义分析代码 | 本文；Zenodo | Zenodo: https://doi.org/10.5281/zenodo.21918740

打开新标签页

实验模型与研究参与者详细信息

本研究是对先前已发表的神经网络培养记录的重新分析。关于物种、品系/系、基因型、年龄、培养及伦理监督的细节，如各数据集所引用的原始出版物中所述，并在"数据集"段落中进行了总结。我们未分析性别因素的影响。

数据集

我们使用了来自神经网络培养物的四个不同的微电极阵列（MEA）数据集，包括大鼠和人类诱导多能干细胞（hiPSC）衍生培养物。其中三个数据集是公开可用的，另一个是向原始出版物作者请求获取的。对于所有数据集，我们使用了处理后的尖峰数据，该数据由每个电极的时间戳尖峰列组成。每条记录对应每个孔的单一时间点，并被视为独立样本；不同体外天数（DIV）的记录未被合并。

阻断抑制

该数据集由 Vinogradov 等人发表，包括分离的小鼠皮层神经元的记录（数据在线可用，地址为 https://figshare.com/articles/dataset/MEA_recordings/27110542）。神经元培养在24孔MEA板上，每个孔包含12个金电极。对于两个批次，自发性网络活动在体外第17天和第18天记录（批次1和批次2）。为阻断突触抑制，在每个批次的24个孔中的10个孔中加入40 μM bicuculline，结果产生N = 20条阻断抑制记录和N = 28条对照记录（总共N = 48条）。经过10分钟稳定期后，网络活动记录了20分钟。虽然原始研究未讨论爆发形状，但由于其对网络活动的系统性调节，该数据集非常适合我们的分析。

CACNA1A 单倍剂量不足

该数据集由 Hommersom 等人提供，包括人类诱导多能干细胞衍生培养物。对照细胞系通过 CRISPR/Cas9 修饰产生 CACNA1A 单倍剂量不足，并验证了多能性标志物和三系分化潜力。原始研究未报告细胞系是否进行了支原体污染检测。培养物在DIV 42时使用48孔MEA板记录了5分钟，每个孔包含16个电极。数据包含四个批次，其中24条记录来自 CACNA1A+/? 培养物，57条记录来自对照（未显示爆发活动的记录被排除在分析之外）。原始研究报告了 CACNA1A 单倍剂量不足爆发相对于对照的峰值时间略有延迟，这表明存在与我们的研究相关的系统性爆发形状差异。

Kleefstra 综合征

该数据集由 Mossink 等人发表，包括从对照和疾病患者皮肤成纤维细胞衍生的 hiPSC 培养物。数据集包含来自四名 Kleefstra 患者、五名健康对照和五个等基因对照（修饰过的患者细胞）的培养物。我们将 N = 287 条记录（每孔一条）的对照组分组，与 N = 58 条 Kleefstra 记录进行比较。该研究还包含三名 MELAS（线粒体脑病、乳酸酸中毒和卒中样发作）综合征患者，我们将其排除在分析之外，因为原始研究报告 MELAS 爆发形状与对照非常相似。相比之下，Kleefstra 爆发被报告具有不同的起始动态，使其适合我们的研究，并将分类简化为 Kleefstra 与对照网络的二分类问题。记录跨越 DIV 27 到 35，即 Mossink 等人确定的稳定爆发参数窗口，每个批次的每个孔贡献一条记录。根据原始研究的补充信息，细胞系已通过鉴定（短串联重复序列分析、单核苷酸多态性芯片分析或核型分析）。未明确说明是否进行了支原体污染检测。

原代大鼠神经元培养物

该数据集由 Wagenaar 等人发表，包括来自 Wistar 大鼠胚胎的前部皮层神经培养物，在妊娠第18天提取。同一窝的细胞合并为一个批次，总共产生八个批次（窝）。因此，每个批次是若干具有共同遗传背景的胚胎的混合物。培养物在无预先稳定期的情况下记录了30分钟。我们选择了体外第7天到第35天之间的数据进行我们的分析，结果为八个批次（1–8）分别有74、137、102、28、18、40、28和24条记录，总共 N = 451 条。虽然原始研究重点关注爆发统计量的发育轨迹，但也包含了爆发内动态的简单特征，如起始和结束持续时间。他们报告爆发持续时间随发育而缩短，特别是爆发起始持续时间，如原始研究图6所示。

方法细节

预处理

爆发识别

我们使用了由 Bakkum 等人先前提出并适用于各种设置的尖峰时间合并与爆发检测方法。首先，我们将所有记录单元（电极）的尖峰时间合并为单个全网络尖峰列。然后我们使用检测算法基于峰间间隔（ISI）找到爆发的开始和结束。爆发被定义为当 ISI 低于预定义阈值 ISIstart 时开始，当随后的 ISI 超过第二个阈值 ISIend 时结束。虽然这些阈值通常相同，但可以独立设置。

为过滤有意义的网络事件，我们应用了附加标准。检测到的爆发仅在其满足以下全部条件时才被接受为爆发：(1) 最小爆发持续时间 BDmin，(2) 最小爆发间间隔 IBImin，以及 (3) 最小爆发大小（尖峰数量）BSmin。我们还应用了可选的最小爆发长度和最小发放率过滤器。在某些情况下，还应用了 4 ms 平滑核以减轻极短爆发中波动的影响。所选参数如表 S2 所示。

为测试鲁棒性，我们实现并运行了使用附加爆发检测算法的分析。该算法首先基于 ISI 在单个电极上检测爆发事件。我们根据电极数量 N 选择了从我们的爆发检测到该算法的简单确定性参数映射：ISI 阈值 ?SIstart/end = ISIstart/end · √N，最小爆发大小 B?Smin = BSmin / N（所有其他参数相同）。然后，为检测网络事件，我们对同时活跃的电极进行阈值化，选择 20% 的阈值。选择此阈值的理由是我们认为高活动时间窗口对爆发形状很重要。

表 S4 显示，这些算法在爆发检测上基本一致。即使对于一致性最低的数据集（83.9% 重叠时间），平均爆发形状看起来也非常相似（图 S11A）。对所有数据集运行结果显示，关键结果得以保留，即通过爆发形状特征增强传统特征可以显著捕获更多关于网络表型的信息（图 S11B 和 S11C）。

我们还测试了爆发检测对主参数扰动——峰间间隔阈值 ISIstart/end 的敏感性。我们系统地将 ISIstart/end 从 -50% 变化到 +50%，发现对于大多数扰动，检测保持稳定（表 S5）。对于高达 20% 的小扰动，存在几乎一一对应的匹配。仅在 -50% 或 +50% 的极端扰动下，爆发存在一致性（F1 值）和重叠时间（Dice 指数）才下降超过 10%，且仅在少数数据集中出现。我们估计在 F1 值和 Dice 指数 10% 变化范围内是安全的，因为替代爆发检测方法的偏差更大，但仍保留了我们的关键洞察（表 S4；图 S11）。

爆发审核

在运行自动检测算法后，我们使用自定义交互工具（图 S1）手动审核了爆发检测。该工具使我们能够高效检查大量记录的爆发检测质量。该工具将检测到的爆发叠加在发放率轨迹和光栅图上，可选择放大以进行更仔细的检查。如果初始参数未能准确捕获网络爆发，我们可以调整爆发检测参数，直到达到可接受的水平。平均而言，88.8% 到 98.9% 的电极参与了这些爆发，具体取决于数据集（图 S13），确认了它们是网络事件。

基于 Wasserstein 距离的光谱嵌入

归一化

为了使用 Wasserstein 距离比较爆发，我们在时间和幅度上对爆发进行了归一化，从而隔离出形状特征。我们将持续时间归一化为长度 1，将其离散化为 50 个时间箱，并根据尖峰时间计算每个箱的发放率。最后，我们通过其积分（所有箱的总和）对每个爆发形状进行归一化，以确保总积分为 1。这些步骤产生了具有单位持续时间和单位积分的标准化爆发形状，随后用于进一步分析。

Wasserstein 距离

为量化爆发形状之间的相似性，我们选择了 Wasserstein 距离度量（也称为 Earth Mover's distance）。选择该度量是因为它对微小偏移具有鲁棒性，为相似的爆发形状赋予小距离，为不相似的爆发形状赋予大距离。对于两个爆发 b1 和 b2，一维 Wasserstein 距离 d(., .) 可以从它们的累积和计算得出：

d(b1, b2) = Σj |Σi=0^j bi1 - Σi=0^j bi2|，（公式 1）

其中 bi? 是爆发 n 在箱 i 中的值。

我们比较了替代距离度量（KL 散度、欧几里得距离、相关性），发现 Wasserstein 距离提供了最一致的改进，尽管 KL 散度对于特定数据集和爆发形状可能是可行的替代方案（图 S9 和 S10）。

K 近邻（KNN）图

为捕获单个爆发形状之间的关系，我们构建了 KNN 图，将每个爆发形状与其 K 个最相似的爆发形状连接起来。所有爆发形状之间的成对距离使用 Wasserstein 距离计算。最近邻数量 K 的选择使得每个节点连接到约 1% 的所有节点（类似于 Kobak 和 Berens 对 t-SNE 困惑度的建议），在我们的数据集中 K 在 21 到 150 之间。

该KNN图作为KNN聚类和谱嵌入的基础。

**谱嵌入**

为从KNN图中获得低维嵌入，我们进行了谱嵌入。首先，我们从KNN图的非对称邻接矩阵A构造了对称邻接矩阵 $A_{\text{sym}} = \frac{1}{2}(A + A^T)$。然后，我们计算了图拉普拉斯矩阵 $L = D - A_{\text{sym}}$，其中D为度矩阵，使得L的每一行（和列）之和为零。接着，我们使用了L对应于最小非零特征值的特征向量作为低维坐标，因为它们能最好地保留底层图的连接结构。

**记录的分类**

为了了解哪些爆发特征包含有生物学意义的信息，我们使用各种模型和特征集对记录进行了分类。我们采用了两种主要分类策略。首先，我们使用基于KNN的分类方法，以利用单个爆发形状中的丰富信息，并将这些信息聚合到记录层级。其次，我们应用XGBoost基于汇总特征对记录进行分类，包括传统特征、形状特征和组合特征集。对于所有分类任务，我们使用平衡准确率和重复分层随机拆分（100次，80-20训练-测试划分）来评估模型性能，以实现稳健的交叉验证。

**基于单个爆发形状的KNN预测**

为了利用单个爆发形状中存在的丰富信息，我们进行了基于KNN的分类。为防止数据泄漏，我们首先确保测试集记录中的单个爆发未被包含在训练集中。然后，对于这些"测试爆发"中的每一个，我们在"训练爆发"中识别其K个最近邻。我们使用基于这些邻居的投票方案预测"测试爆发"的类别，这些邻居属于具有已知类别标签的训练记录。为考虑类别不平衡和记录特异性偏差，我们对每个邻居的投票进行了加权。权重计算为：该爆发所在记录中爆发数量的倒数 $N_{r_j}$ 和该类中记录数量的倒数 $N_{c_j}$。

单个爆发的类别投票向量 $v_i$ 计算为：

$$v_i = \sum_{j=1}^{K} \frac{1}{N_{r_j}} \cdot \frac{1}{N_{c_j}} \cdot b_{ij} \quad \text{(公式2)}$$

其中 $b_{ij}$ 是一个独热向量，表示测试爆发i的第j个邻居的类别标签。为获得一个测试记录的预测，我们对该记录中所有爆发组成的投票向量 $v_i$ 取平均值。平均投票最高的类别被分配为该记录的预测标签。

**基于汇总统计量的XGBoost预测**

为评估形状特征相对于传统爆发特征的预测能力，我们基于三组记录级汇总统计量训练了XGBoost：形状特征、传统特征和组合特征集。

作为预测准确率基线的传统爆发特征，是从常用的爆发指标中导出的。

- 平均放电率（MFR）：培养物放电率的平均值
- 爆发尖峰率（BSR）：爆发放电率的平均值
- 网络爆发率（NBR）：爆发次数除以记录持续时间
- 网络爆发持续时间（NBD）：爆发持续时间的平均值
- 随机尖峰百分比（PRS）：未分配给爆发的尖峰百分比

我们排除了van Hugte等人使用的一些特征（平均爆发率（MBR）、爆发持续时间（BD）和高频爆发（HFB）），因为它们依赖于单通道爆发和高频爆发，这些在我们的分析中未被计算。

我们的爆发形状特征基于谱嵌入，即KNN图的低维投影。对于每个记录，我们对二维中的单个爆发形状位置取平均。这为每个记录产生了两个形状汇总特征（谱嵌入维度1和维度2）。虽然可以使用多于或少于两个的维度，但我们始终选择两个维度，与可视化中使用的一致。组合特征集通过将五个传统特征和两个爆发形状特征组合获得。

对于所有XGBoost模型，我们采用了嵌套交叉验证方法，包括外循环（测试）和内循环（验证）。外循环由重复分层随机拆分组成（100次，80-20划分）用于测试（与KNN分类相同）。内循环通过5折交叉验证中的网格搜索优化XGBoost超参数。我们报告测试集上的平衡准确率来评估每个特征集的预测能力。

**通过Shapley值衡量特征重要性**

为量化传统特征和形状特征的相对贡献，我们计算了在组合特征集上训练的XGBoost模型的Shapley值（特征重要性）。Shapley值是来自合作博弈论的一种稳健方法，用于在"玩家"（此处为特征）之间公平分配"收益"（即模型的预测）。该方法用于在考虑复杂交互作用的同时量化特征的重要性。Shapley值的一个特别有用的性质是可加性，这使得我们能够计算传统特征集的总重要性，并与形状特征进行直接比较。

**定量分析与统计分析**

统计细节分别在相应的图表说明和正文中报告。除非另有说明，N为记录数量，数值为均值±均值标准误。样本量由每个数据集中可用的记录数量固定，因此不是由先验的效能分析确定的。显著性定义为 p < 0.05，我们使用星号表示显著性水平：*：p < 0.05，**：p < 0.01，***：p < 0.001。为确保可重复性，所有随机化过程（例如重复的训练-验证-测试划分）均使用固定随机种子。

**排除标准**

未检测到任何爆发的记录被排除。在hiPSC疾病数据集中，三名MELAS患者被排除，因为原始研究报告称MELAS爆发形状与对照组无法区分。在大鼠发育数据集中，分析仅限于体外培养7至35天之间的记录。未进行进一步的记录级排除。

**模型准确率比较**

对于模型准确率比较（图3B–3D），我们认识到标准方法（如t检验）会导致显著性膨胀，因为交叉验证的准确率由于数据的重复使用而相互依存。因此，为稳健地比较模型准确率，我们使用了对均值标准误的修正 $\hat{s}_d^2 = \hat{\sigma}^2\left(\frac{1}{n} + \frac{1}{n_{te}} + \frac{1}{n_{tr}}\right)$，其中n为试验次数，$n_{te}$ 为测试集大小，$n_{tr}$ 为训练集大小。根据我们的交叉验证参数（n = 100，$n_{te} = n_{tr} = 14$），我们计算了配对双尾t检验。为测试跨数据集的总体结论，我们使用随机效应模型进行了荟萃分析，并应用了Wald检验。

**特征比较**

手动设计特征之间的差异（图4D和4E）使用非配对t检验进行比较。

**距离比较**

我们使用多步骤过程分析记录之间的距离（图5B和5D）。首先，我们通过对其组成爆发形状之间的距离取平均来计算记录之间的成对距离。然后，我们根据特定标准过滤成对距离，例如包含所有距离（"随机"）或仅包含来自同一受试者记录之间的距离（"受试者"）。为平衡数据集并避免人为夸大样本量，我们通过对所有向外距离取平均为每个记录计算一个平均距离。最后，我们使用双尾t检验确定这些平均距离是否存在显著差异。

**补充信息（1）**

PDF（5.02 MB）文档S1。图S1–S15和表S1–S5