针对混沌系统的科学机器学习为神经群体学习降阶方程
《PRX IntelligenceNEW》:Scientific Machine Learning of Chaotic Systems Learns Reduced-Order Equations for Neural Populations
【字体:
大
中
小
】
时间:2026年10月28日
来源:PRX IntelligenceNEW
编辑推荐:
从复杂动力系统中提取可解释的数学模型是困难的,尤其是对于伴随噪声实验数据观测到的混沌动力学而言。我们提出PEM-UDE,一种将预测误差方法与通用微分方程相结合的方法,用于从有限的、受噪声污染的观测中发现控制方程。预测误差反馈平滑了混沌优化问题;对于模型类内部生成的无噪声数据,它保
从复杂动力系统中提取可解释的数学模型是困难的,尤其是对于伴随噪声实验数据观测到的混沌动力学而言。我们提出PEM-UDE,一种将预测误差方法与通用微分方程相结合的方法,用于从有限的、受噪声污染的观测中发现控制方程。预测误差反馈平滑了混沌优化问题;对于模型类内部生成的无噪声数据,它保留了数据一致的零损失集合,而噪声和模型设定误差引入了依赖于增益的稳定性-偏差权衡。零损失集合的保留并不能保证唯一的结构可辨识性。我们在两个基准混沌系统——罗斯勒吸引子和一个真实电路——上测试了该方法,即使在一个观测维度的噪声是信号幅值的5倍时,也能恢复正确的函数形式。该方法还能接受系统的先验知识作为初始函数形式,我们利用这一点来学习神经电路方程,该方程考虑了稀疏连接性,这是传统神经质量模型所缺乏的特征。应用于Izhikevich神经元群体时,PEM-UDE生成了一个多尺度神经质量模型,将单神经元参数与宏观网络动力学联系起来,并预测了连接密度、主导振荡频率和同步性之间的关系。我们针对大鼠和人皮层的三个颅内记录数据集测试了这些预测。对于神经科学应用,学习到的方程是特定模拟Izhikevich网络族的降阶闭合;实验记录提供了预测频率和同步趋势的间接一致性检验,而不是将方程直接拟合到记录上。I. 引言发现控制复杂系统的方程在过去一个世纪的大部分时间里一直是动力学研究的目标,通常需要应用领域的专家知识和统计力学的一些理解[1–5]。这种方法已经在许多领域取得了成果,从气候模型[6,7]到化学动力学[8,9]、人类社会互动[10]和神经计算[11]。最近的工作已经产生了自动化方法,以协助通常由专家驱动的用微分方程组描述数据的过程[12–15]。符号回归技术,包括稀疏识别方法(SINDy;[12])和遗传算法[16],直接从数据拟合方程,但通常在噪声数据下效果较差[17]。通用微分方程(UDE)结合了机器学习和科学计算的思想[18],将人工神经网络嵌入微分方程框架中,直接从数据学习未知动力学。网络的训练受到两种方式的约束:一是对动力学施加已知的函数形式,围绕其学习附加项[18];二是对学习过程本身施加结构约束[19]。因此,这里的贡献是具体的:我们不是将观测器、UDE或符号回归作为独立的思想引入,而是将类似观测器的预测误差反馈与UDE训练相结合,使得混沌轨迹可以在学习到的残差场转换为显式符号表达式之前被拟合。这些方法在混沌系统中遇到困难,而混沌系统在神经科学和其他生物学背景中很常见[20–23]。即使是简单的物理系统,当某些参数超过分岔点时,也能表现出复杂、混沌的解[24],而神经动力学通常表现出混沌振荡器的特征[22,25]。混沌的存在给拟合过程带来了数值问题:混沌系统的轨迹在Lyapunov时间之后无法明确定义。任何数值常微分方程(ODE)求解器都会在每一步引入截断误差,而对输入的敏感依赖性意味着求解器输出在Lyapunov时间之后不久就会在产生的轨迹中出现误差。数值模拟给出有意义轨迹的唯一保证由影子引理提供,该引理指出,在初始条件周围存在一个球,在该球内数值求解器的轨迹与球内其他潜在初始条件之一的轨迹相匹配,有效地描述了吸引子周围解密度[26]。影子引理有两个后果,使得混沌识别变得困难。混沌轨迹的发散特性导致用于计算梯度的所有敏感性分析方法不稳定。通过求解器进行的标准自动微分会产生不准确的导数估计,除非应用专门的修改[27–31]。因此,轨迹校准的标准方法(拟合由ODE求解器数值解定义的损失函数)在混沌系统背景下不是一个收敛算法。噪声观测数据也会产生一个在训练期间无法区分相邻轨迹的损失函数,导致拟合动力学的不准确性[26]。由于解对初始参数具有指数依赖性,数据的无穷小变化可能使数据点在时间序列轨迹上指数级地更远。由于这些性质,混沌系统的长期时间序列只能基于统计或遍历性质进行概率解释[32–34]。因此,通过范数拟合时间序列的传统技术在混沌背景下是不良定义的。在实际模型中,我们还只有部分可观测性,这使得任何动力系统的参数拟合都很困难。这在神经科学中尤其如此,因为神经相互作用形成跨多个尺度运作的复杂系统:从神经递质活动等分子过程到可以使用脑电图(EEG)观察到的大脑范围电场[21,35,36]。除了描述这些多尺度相互作用所固有的复杂性外,脑动力学通常接近临界点,允许根据任务需求、代谢资源和神经电路的整体健康状态在不同动态状态之间切换[37–39]。接近临界点的运作在信息处理中提供了灵活性,但也使神经动力学对微小扰动敏感[22,40]。传统的神经元群体建模方法,称为神经质量模型,专注于捕捉特定的宏观属性,而没有系统性地考虑离子通道电导等微观生物物理参数[41–43]。下一代神经质量模型(NGNMMs)此后提供了群体活动的解析表达式,既保留了单个神经元中存在的微观细节,也捕捉了时变同步性等涌现的宏观属性[44–47]。这些解析方法通常依赖于Ott-Antonsen拟设[48],以避免对平均场动力学进行Fokker-Planck推导[49,50]。一个注意事项是,这些推导捕捉了许多理想的多尺度动力学,但依赖于非生物物理假设,特别是神经元之间的全连接[44]和系统中的电流守恒[50]。最近的工作将这些方法扩展到考虑稀疏连接和二次积分放电神经元网络中波动引起的涨落[51]。这样的扩展仍然局限于简化的神经元模型,缺乏生物物理详细的适应和恢复动力学。将解析平均场推导扩展到生物物理详细的神经元模型(例如,带有适应电流的Izhikevich神经元),同时放松连接性假设,仍然是一个未解决的问题,并推动了数据驱动方法的出现。底层神经生物学与这些推导所基于的全连接网络假设之间的脱节是有影响的。这一假设对于深部和皮层下结构大致成立[45],但在皮层中远非生理性的,典型的连接性只有1%–5%[52]。皮层稀疏性不是偶然的:多条证据表明,皮层中稀疏、结构化的连接模式对于支持高级认知过程的结构化信息流至关重要[53,54]。因此,在捕捉群体活动的方程中考虑稀疏性对于准确的大脑模型是必要的。为了解决这些局限性,我们提出了一种训练UDE以拟合混沌系统动力学的方法,并用它来学习跨连接水平的模拟神经群体模型的降阶闭合项。我们描述了这一过程,它使用预测误差方法(PEM)[19,55]在混沌系统数据上训练UDE,然后从UDE中提取学习表达式的符号形式以推广动力学。我们表明,这种方法可以学习混沌系统中的未知状态,准确恢复传统UDE方法无法恢复的罗斯勒系统中的隐藏状态,并且通过从系统中移除对参数选择的敏感依赖性来实现这一点。该方法还在不完整观测下学习混沌电路的动力学,其中观测数据的一个维度被原始信号幅值5倍的噪声所破坏。在这些条件下,直接应用于数据的序贯阈值最小二乘(STLSQ)无法恢复正确的动力系统,而PEM-UDE成功;因此,我们认为该方法比纯符号回归更数据高效。最后,从脉冲网络模拟的放电率和电压观测中,我们学习了现有神经质量模型的降阶修正。符号修正推广到模拟的连接概率,并揭示了连接性变化如何改变该网络族的节律动力学。后来与记录的比较检验这些定性预测是否与区域生理学匹配;它不声称从记录中直接发现方程。这些对照比较在研究的条件下分离了PEM相对于其他相同的UDE和直接STLSQ的贡献;它们不是所有混沌系统识别方法的全面基准。完整流程的图形概述如图1所示。图1. PEM-UDE方法示意图。我们从具有时变活动的物理或生物系统的观测(数据)开始。从这些和先前的工作中,我们形成对被观测系统的初始描述(已知约束)。然后,我们拟合一个UDE来描述该初始估计未涵盖的动力学部分,使用预测误差方法帮助在混沌系统上训练UDE(学习动力学)。然后,我们从UDE描述的动力学中提取符号形式(生成方程),并检查该符号形式是否正确预测超出初始观测的活动(验证)。II. 结果 A. PEM-UDE方法的正式描述我们在这里给出组合PEM-UDE/符号回归方法的正式定义,实现细节见第IV节和随附代码。相关的代码库包括一个使用玩具数据的演示,可在标准计算机上运行,供读者尝试该方法。方法的图形概述如图1所示。我们从导致相关时间序列数据集的物理或生物系统的观测开始。目标是一个描述被观测系统动力学并准确预测未来动力学的微分方程组。从观测和先前对类似系统的研究中,我们对动力学进行初始猜测,利用来自物理系统的先验信息来构建系统行为的第一近似[18]。然后,我们建立一个UDE系统,直接从观测数据学习未知动力学。训练UDE的标准方法难以准确拟合混沌系统,如下面的结果所示。为了解决这个问题,我们开发了一种训练方法,在模拟过程中将PEM动力学与UDE集成,使UDE能够从数据中学习。训练完成后,我们用基于稀疏线性求解器的符号回归(STLSQ)或基于进化算法的符号回归[16]从UDE网络中提取动力学的符号表示,具体取决于具体背景和数据可用性。符号描述使我们既能预测未来时间序列,又能对底层动力学进行解析观察。预测误差项是一个训练设备:在观测器修正被减少或移除后,符号方程被评估为自主动力学,因此成功要求学习到的向量场本身再现观测到的动力学结构。1. UDE公式遵循先前的工作[18],带有一般过程的加性修正的通用微分方程,给定一些可观测状态,为:(1a)(1b) 这里,我们用粗体变量表示向量。是系统的已知(或假设)动力学,参数为,表示这些动力学在初始条件下在参数值处的解,是外部刺激,是由参数化的观测函数,可以显式依赖于时间(例如,由于噪声过程或外部刺激对观测的影响),是基于其参数训练以学习系统未知动力学的通用逼近器网络。为符号简洁,我们省略观测对参数和初始条件的参数依赖,但保留网络参数,因为我们想要优化这些参数而保持其余固定。在这项工作中,我们聚焦于常微分方程系统,但该方法可扩展至随机和延迟微分方程系统,正如先前使用UDE的工作[18]所做的那样。对于观测数据集,训练UDE使用代价函数,该函数最小化每个观测时间点处的误差,使得 (2) 其中,是带有神经网络权重的ODE在测量时间处的观测解。在所有依赖观测数据的方程中,我们将任何观测数据表示为。在本工作中,我们选择为均方误差,但其他代价函数适用于不同场景,取决于测量噪声如何影响测量信号。随后,我们使用预测误差方法在此损失函数上为混沌系统训练UDE。
2. PEM公式
PEM是用于将动力系统参数拟合到观测数据的标准优化技术[19,56]。与方程(1)类似,我们考虑一个一般过程,它产生某个可观测状态,并带有某些测量噪声:(3a) (3b) 在此表述中,是时变系统的状态,带有我们希望拟合的参数。观测状态由过程的测量函数定义,该函数由参数化,并带有测量误差。这里,我们假设到存在一一映射(即,中没有依赖于多个元素组合的项)。关于PEM可应用条件的更一般性讨论,参见参考文献[19,56]。例如,不限于是的线性函数,还可以执行卷积或包含延迟[56]。在我们的所有示例中,的少数分量为零,表明这些分量未被观测。估计误差由函数加权,并在每个时间步添加到原始动力学中,以产生调整后的估计和观测值,使得 (4a) (4b) (4c) 其中是观测数据插值到求解器所需的任意步长处。误差增益因子可以采用不同形式。在本工作中,我们将设为常数,这是先前工作所称的参数化观测器[56]的最简形式。这种观测器形式与基于同步的混沌轨迹稳定密切相关[57];此处它被用于UDE训练循环内部,而非作为独立的状态估计程序。如果有充分理由怀疑观测质量随时间变化(例如,在某些时间区间噪声显著更高),那么时变的——例如,将方程(4)扩展为类似扩展卡尔曼滤波器的形式,其中为卡尔曼增益——可能在拟合时提供更好性能。尽管如此,先前工作已表明,参数化观测器在实践中优于扩展卡尔曼滤波器,因为随时间调整超参数较为困难[56]。为简化符号,在后续章节中我们将使用,而不失一般性:如果有数据驱动的理由要求,可以实现扩展卡尔曼滤波器。
3. 将PEM-UDE结合为单一方法
我们通过引入通用逼近器来学习未知动力学,同时使用预测误差反馈来稳定训练轨迹,从而将UDE和PEM结合为单一系统:(5a) (5b) (5c) 超参数控制优化所依据的损失景观,我们将证明这对于将UDE拟合到混沌系统是必需的。在本工作中,我们仅考虑同时具有初始函数形式和观测数据的系统。对没有初始函数形式概念的未观测潜在状态的估计是一个独立问题,在其他地方讨论[45]。
4. 从PEM训练的UDE中符号恢复各项
在使用组合PEM-UDE方法拟合系统动力学之后,我们从通用逼近器中提取动力学的符号形式。符号形式通常比黑箱神经网络更稳健地泛化到未来情境[3,18],并允许对学习到的函数形式进行直接分析。UDE学习到的项可以通过STLSQ(通过岭回归从库中学习特征[12])或遗传算法符号回归方法[16]恢复。每种方法各有优劣,适当的选择取决于数据结构和已知动力学的结构。当系统函数形式存在复杂约束时,遗传算法是更好的选择[16],而STLSQ在需要大量数据点才能实现精确拟合的系统上效果更好。我们在伴随此方法的代码中提供了两种方法,并在下面各节中讨论我们的选择。由于这些方法并非本工作的创新点,进一步讨论留待第四部分。在此流程中,符号回归并非用于声称神经网络具有内在可解释性;它被用作对训练残差场的稀疏化投影,随后对恢复的方程进行前向模拟。
B. PEM-UDE目标的解析性质
在数值演示之前,阐明预测误差项对优化问题的作用是有益的。设 (6) 表示建模的向量场,设求解标准UDE,设求解具有相同初始条件的PEM-UDE。对于线性观测映射,PEM项可写为,其中为标量增益和直接逐分量观测。这是方程(5)的矩阵形式;对于非线性观测映射,同样的局部论证成立,只需将替换为观测雅可比矩阵。我们比较连续时间中的普通和PEM辅助轨迹损失:(7) (8) 其中方程(2)是采样版本。
为何最小值不移动。首先,考虑隔离优化几何的理想化情况:数据无噪声且由模型类中的某个成员生成。因此,存在使得 (9) 具有相同初始条件和唯一解。在处,预测误差恒为零:(10) 因此,PEM反馈在正确的动力学处恰好关闭:(11) 由唯一性,有,所以。由于损失是平方范数,对所有成立,且对于任何固定增益,是PEM-UDE目标的全局极小值点。
同样的论证也表明,PEM不改变观测状态目标的零损失集:(12) 如果,则残差在所有时刻为零,因此反馈项沿整个PEM轨迹消失,且该轨迹也求解无强迫UDE。反之,任何零损失的标准UDE轨迹具有零残差,不受PEM项影响。因此,反馈重塑了远离数据一致流形的景观,但不移动该流形。在结构可辨识性和充分可观测性条件下,这个共同零集归结为真实动力学;在部分可观测性或模型误设下,不同的向量场仍可能是观测等价的,因此该陈述应理解为数据一致最小值的保持,而非唯一全局可辨识性的保证。
PEM如何使混沌目标可训练。最小值保持本身并不能解释优化器为何能找到最小值。增益之所以重要,是因为它改变了误差动力学的稳定性。设,定义模型误差场 (13) 并设。将PEM-UDE在真实轨迹附近线性化,得到 (14) 对于具有标量增益的完全观测,有,因此PEM项将线性化误差动力学的每个指数偏移。具有指数的非稳定方向在以下条件下变为收缩的:(15) 这为损失景观平滑提供了简单的物理解释:观测器在误差污染损失及其梯度之前,阻尼了邻近轨迹之间指数增长的分离。在R?ssler示例中,李雅普诺夫时间为或,而经验上有效的增益–0.25恰好位于该阈值之上,与补充材料[58]的图S3一致。对于部分观测,同样的阻尼发生在观测到的非稳定子空间中;未观测的非稳定方向必须与观测耦合,反馈才能稳定它们。
设是由生成的转移矩阵。如果观测到的误差动力学是一致收缩的,即对,常数变易公式给出 (16) 其中对具有标量增益的完全观测有。因此,最小值附近的PEM损失由自治向量场误差控制,而非指数发散的轨迹:(17) 因此,将混沌轨迹匹配问题转化为稳定的、观测器正则化的向量场比较,而方程(12)保持精确最小值固定。权衡在相同缩放中可见:增加会抑制敏感依赖性并平滑景观,但它也通过因子压缩邻近向量场之间的差异,从而产生在大增益处看到的平台。
含噪观测与增益选择。对于测量噪声,线性化误差动力学变为 (18) 观测器现在既阻尼混沌发散,也将测量噪声注入模拟状态。设表示与稳定转移算子的卷积。到主导阶,观测残差为 (19) 并且,对于与模型误差无关的零均值噪声,(20) 第一项是不同自治向量场之间的有用对比;第二项是依赖于的噪声基底。小增益使混沌不稳定性保持活跃。非常大的增益使观测器强烈跟踪测量时间序列,以至于模型误差被反馈隐藏,方程(20)中的信息性对比变得浅平。这激发了一个有用的区间:(21) 其中下界稳定轨迹,上标尺由模型误差项与噪声基底相当的点设定。精确的上标尺取决于观测算子和测量噪声的谱,因此并非普适的。在实践中,当李雅普诺夫尺度可从数据估计时,我们选择刚好高于该尺度,并检查所得损失景观;对于神经质量模型,我们则使用延续调度,将从0.3退火到0.1,然后再进行符号提取,以减少观测器对学习到的自治动力学的影响。
神经网络项学到了什么。网络并非旨在学习整个向量场。在缺失物理表述中,它仅学习 (22) 即尚未由已知物理表示的部分。这个较小的目标是理解为何小网络可以工作以及大网络并不自动更好的核心。对于含噪数据,拟合的残差场具有通常的偏差-方差结构:(23) 其中是噪声平均后的学习网络。增大网络规模可减少偏差,但会增加方差项,即跟随特定噪声实现的部分。
这种容量依赖性即标准的偏差-方差权衡[59]。局部线性化可以明确说明。如果是所得特征矩阵,是活跃自由度数,是回归目标中的有效噪声方差,那么使用岭/权重衰减正则化 (24) 因此,噪声拟合项随有效参数数量增加而增长,随信息样本数量增加而减小。最优模型规模是有限的:它是偏差变得相对于方差成本足够小的最小网络。由于已经包含主要物理结构,残差的复杂度远低于完整动力学,因此该最优规模可以很小。这就是我们使用紧凑网络而非大型黑箱神经ODE的原因。
这些论证减轻但并未消除网络吸收归因于已知物理的动力学的风险。在我们的实验中,及其假设参数保持固定,因此网络目标是匹配观测轨迹所需的残差。然而,数据主要在采样状态上约束总场:灵活的网络可以补偿错误设定的已知项或采用不同的轨迹外延。保持较小并使用AdamW权重衰减通过如下形式的目标来限制这些替代方案:(25) 最后,符号回归步骤施加额外的稀疏约束,学习到的符号方程在PEM项减少或移除后得到验证。这些检查减少了但无法排除采样吸引子之外的残差不可辨识性。这些结果也有助于解释为什么STLSQ在PEM-UDE输出上成功,而在表I中的原始含噪数据上却失败:通过训练后的UDE进行积分可产生平滑的生成轨迹,而对含噪数据进行数值微分则会以与功率成正比的比例放大高频噪声。表I.将STLSQ直接符号回归应用于含噪Petrzela-Polak电路数据,在PEM-UDE方法成功的相同数据和噪声条件下,未能恢复正确的动力学。拟合表达式 0.1 0.2 0.3 0.4 0.5 C. PEM-UDE在拟合混沌系统方面优于标准UDE
我们首先展示PEM-UDE方法如何拟合R?ssler系统(图2)。选择该系统的原因是它是一组简单的常微分方程,仅含有一个非线性项,但在适当的参数集下会产生混沌动力学[24]。在本例中,我们假设已知的函数形式以及的函数形式,但对的函数形式一无所知,因此必须从时间序列观测中学习动力学。仅面对这一挑战,传统UDE无法良好拟合[图2(a)]:UDE学习到的动力学通常收敛到时间序列的均值,即目标函数中的局部极小值。相比之下,PEM-UDE方法在训练后成功收敛,即使在整个训练时间的第一个10%之后,性能差异也已明显,而此时标准UDE已接近其局部极小值。
图2. PEM-UDE组合方法能够学习传统UDE无法学习的混沌系统动力学。(a) UDE系统(未能准确捕捉动力学)与PEM训练的UDE(成功)的训练误差比较。在整个训练会话进行到10%时,精度差异已经很大(插图),此时未辅助的UDE不再取得进展。(b) PEM的应用平滑了原本困难的参数损失景观,其中超参数调节陡峭程度。(c) 拟合动力学的均方误差,UDE与PEM训练的UDE对比。(d) 使用STLSQ替换PEM训练的UDE所得到的表达式能够捕捉原始动力学(前200个时间点)并在未来持续准确描述(最后100个时间点)。(e) R?ssler吸引子的表面显示出相空间的混沌填充,在平面附近有一个盆地,并沿轴有偶尔的爆发。(f) PEM训练的UDE在吸引子表面每个点处的相对误差。
性能提升的原因在于PEM方法改变了UDE训练的景观[图2(b)]。未修改的景观极其粗糙,包含多个局部极小值,仅在真实全局极小值处(本例中为)有非常陡峭的下降。当在UDE训练期间将PEM添加到系统动力学中时,景观变得更容易接受标准优化技术。的选择也直接调节了景观的陡峭程度。这应结合考虑PEM将如何在其未应用的状态下平滑其他参数的景观来进行。在无噪声可实现的情况下,PEM保留了与数据一致的零损失集,但它可以沿与未观测状态或错误设定参数相关的方向使景观变平。随着更大引入更高级别的滤波,找到整个系统的真实函数形式更多地依赖于不合适参数估计的准确性(见补充材料[58]的图S2)。因此,对于PEM-UDE组合方法,最佳选择是产生足够平滑的景观以便优化,但仍需足够陡峭,使得全局极小值与周围区域清晰可辨,从而允许对未测量参数的更嘈杂估计。
PEM-UDE方法之所以有效,是因为预测误差修正从被拟合的系统中移除了对参数的敏感依赖性,从而中和了混沌系统的定义特征。如补充材料[58]的图S3(a)所示,当较小时,对初始参数选择的任何偏差都会产生发散[在补充材料[58]的图S3(a)中表示为的偏差],但当足够大时(本例中为),这种敏感依赖性被消除,两条轨迹之间的误差缩小。对于足够大的,残差误差取决于参数偏差的大小,并且可以根据该偏差的估计来调节[补充材料[58]的图S3(b)]。对于在模型类内生成的无噪声观测,反馈在与数据一致的轨迹上消失,并保留零损失集[式(12)]。在存在噪声、部分可观测性或错误设定的情况下,PEM转而以较浅且可能偏置的目标为代价来换取稳定性[式(20)]。在原始混沌系统上不稳定的范数轨迹拟合方法,因此在本文所考察的增益范围内得以稳定。
PEM-UDE学习到的动力学很好地拟合了系统,而单独UDE的初始估计即使在系统的Lyapunov时间内也无法提供准确预测[图2(c)]。对于此处选择的参数,系统的Lyapunov时间为,参数空间中的大多数邻近区域具有相似的特征时间(;详见第四节)。单独UDE的平均误差在模拟时间内持续增长,而PEM-UDE拟合动力学在每个时间步的平均误差比动力学本身小2–3个数量级。
为了将这些动力学推广到原始观测之外,我们通过对PEM-UDE系统动力学每个时间步的UDE输出进行STLSQ拟合,提取了方程的符号形式。缺失状态的函数形式直接由STLSQ恢复(),参数初始拟合接近真实值(,,而真实值为,)。在STLSQ收敛到真实值(精确到三位小数)后,一个额外的优化步骤调整了学习到的符号方程的参数,从而使动力学在时间上推广到原始模拟之外[图2(d)]。因此,该流程提供了动力学的解析形式,能够在时间上向前推广,并可用于对物理系统获得洞察,以通常低于4%的误差恢复系统的原始吸引子[图2(e)和2(f)]。因此,该基准评估了隐藏状态恢复、参数误差、窗口外模拟和吸引子几何结构,而不仅仅是短时域轨迹匹配。
D. PEM-UDE可以在严重观测噪声下学习混沌动力学
接下来我们提出PEM-UDE方法拟合混沌系统的第二个示例,与R?ssler系统有两个重要区别。首先,本例是一个物理实现的电路,电路图如图3(a)所示,物理电路的参数(如参考文献[60]所述)见表II。从先前系统的理想化混沌动力学转向真实世界电路也具有实际意义,因为电路中的混沌本身就有其应用[61]。其次,我们对其中一个电路组件添加了含噪观测[,由图3(a)中的蓝色框表示]。该系统的真实方程如式(32)所示,本例中使用的无量纲形式如式(33)所示。在本节中,我们假设对第一个状态存在未知的错误观测,使得观测周围的噪声为真实信号幅度的[图3(b)]。这给PEM-UDE拟合带来了困难:真实动力学[图3(c)]现在被足以扭曲所拟合吸引子形状的噪声所掩盖[图3(d)]。
图3. PEM-UDE方法能够在存在显著观测噪声的情况下恢复混沌电路的正确形式。(a) 产生本节所考察动力学的模拟Petrzela-Polak电路。观测状态为、和无量纲形式,其中和(绿色框)以少量噪声观测,而(蓝色框)以较大观测噪声观测。(b) 在本示例系统中,我们假设对一个状态()的错误观测,其中真实动力学(顶部)被真实信号幅度的观测噪声所掩盖(底部)。(c) 原始吸引子表面,观测维度的真实形式以颜色显示,观测噪声的误差幅度显示。(d) PEM-UDE方法在第三维存在噪声的情况下拟合的吸引子表面(散乱的颜色表明噪声已显著破坏状态的动力学)。(e) 待恢复状态的真实形式显示在顶部,STLSQ应用于PEM-UDE动力学所恢复的形式显示在底部。函数形式(虽然不是精确参数)从PEM-UDE拟合中恢复(表格比较了STLSQ拟合的参数与原始参数)。
表II. Petrzela-Polak电路的物理和无量纲参数。参数 值 0.5 kΩ 0.08 kΩ 1.0 μF 1.0 μF 1.0 H 1.0 mA 无量纲参数 1.25 0.5 1.0
我们假设在拟合时该噪声是未知的,并使用与上述相同的方法将PEM-UDE拟合到损坏的数据[拟合式(34)所示的函数形式]。即使其中一个状态存在大量观测噪声,应用于PEM-UDE拟合动力学的STLSQ也能恢复动力学的精确函数形式[图3(e)]。该方法不提供原始参数,尽管形式是正确的。该问题实际上已从未知动力学问题转变为简单的参数优化问题,可以应用其他方法(例如基于模拟的推断[62])。这表明,即使观测数据非常嘈杂,PEM-UDE方法也能恢复有关系统真实动力学的有用信息。
该示例还说明了为什么拟合训练后UDE的符号形式优于对原始数据集进行STLSQ拟合。作为对图3(e)所示动力学的比较,我们直接使用SINDy方法[12,63]拟合相同的含噪数据集:我们计算平滑导数估计的数据集,并直接对导数估计执行STLSQ符号回归方法。如表I所示,该方法无法在交替方向乘子法(ADMM)超参数的一定范围内恢复真实符号形式,更稀疏的表示会错误地将含噪状态从表达式中完全丢弃。恢复训练后UDE的符号形式更准确,因为UDE提供了动力学的一种无噪声(或大幅降噪)版本,可以以任意密度采样,从而避开了符号回归方法对含噪数据的敏感性。此比较并非旨在穷尽SINDy的所有可能平滑或弱形式变体;相反,它隔离了基于导数的符号发现中常见的一种失败模式,其中对含噪观测进行微分放大了高频误差,而对学习到的UDE进行积分则产生平滑的生成轨迹,以供后续稀疏回归使用。
E. 学习违反解析推导假设的稀疏网络的新表达式
在展示了该方法可以学习混沌物理系统的动力学之后,我们接下来转向推导神经元群体平均场活动的表达式。在本节中,学习到的方程是模拟微观网络的降阶闭合方程,而不是直接拟合到后续考虑的实验记录的方程。研究人脑所需的数据跨越多个尺度:单个神经元活动产生同步动作,定义了相互交互的功能回路[图4(a)–4(c)]。传统的神经质量模型通过建模特定现象来表示神经元群体的平均活动(例如,Wilson-Cowan模型中的兴奋/抑制平衡[41],Larter-Breakspear模型中兴奋群体的离子梯度动力学[43,64]),但与理解中间神经元动力学所需的许多生物物理参数(例如神经递质)缺乏直接联系。NGNMM通过提供平均场动力学的解析推导来解决这一问题,该推导保留了单个神经元的生物物理细节,同时也捕捉了群体层面的涌现特性[44,45]。该推导及其关键假设在第四节中进一步讨论。
图4.学习能够解释神经连接稀疏性区域异质性的神经质量模型,有助于阐明大脑内不同频带的起源。单个神经元产生脉冲序列[图(a)],每个神经元的每个脉冲被编码为光栅图中的一个点。在群体水平上,这种视图使得群体活动易于读取,同步状态显示出清晰可见的爆发[图(b)]。在行为实验中,可以放置电极阵列来记录来自多个脑区的群体水平活动[图(c)]。脑区内局部连接的程度在大脑各处有所不同[图(d)],较深的区域表现出更高的连接性,而皮层区域在区域内的连接则更为稀疏。为了解释这种稀疏性,我们训练了一代神经质量模型,以学习在连接变得更稀疏时捕获涌现群体动力学的项[图(e)]。PEM-UDE学习到的动力学表现出与真实脑区观察到的频率趋势定性相似,随着网络连接日益稀疏,主导频率从θ波向α波转变[图(f)]。更高的稀疏性也限制了网络内可能实现的同步程度[图(g)]。图由BioRender.com创建。
NGNMM提供的平均场动力学解析形式相对于传统神经质量模型有几个优势。首先,单神经元参数与神经群体动力学之间存在直接联系,这一特征仅在某些经过调校以对这些参数敏感的传统神经质量群体中可用[43,64]。其次,NGNMM捕捉了传统神经质量模型未涵盖的神经群体涌现特性。例如,反弹现象——即被抑制的频带活动在受抑制后功率增加(这种现象常在运动皮层研究中观察到)——被观察为仅在α频带振荡的耦合NGNMM的一个特征,而没有任何功能约束迫使这种活动发生[46]。
NGNMM方法也存在局限性,主要涉及Ott-Antonsen拟设[48]所需的假设,这是NGNMM推导的关键组成部分。从生物物理角度看,最关键的后果是假设神经元群体是全连接的。这在具有显著区域内侧支的深层区域(例如海马、内嗅皮层)中近似正确[45]。然而,这些区域的邻近区域往往侧支较少,例如丘脑纹状体连接模式[65]。全连接的假设在皮层中显然不成立[52],皮层通常具有非常稀疏的连接[约为所有可能连接的1%–5%;图4(d)]。因此,在NGNMM中解释稀疏性是一个必须解决的问题,才能使其在新背景下有用。
关于调整NGNMM推导以适应稀疏性的正确解析方法存在分歧[50,51,66,67],因此我们绕开了解析推导的问题,而是使用PEM-UDE方法从模拟的微观网络数据中学习稀疏网络的形式。遵循Izhikevich神经元的解析方法[45],我们从由方程(35)描述的大量神经元群体出发,其突触动力学由方程(36)给出。在全连接情况下,描述该群体活动的解析NGNMM由方程(37)给出。当我们开始修剪图中的连接以引入稀疏性时,这种形式不再准确表示群体活动。因此,我们学习给定降低的连接概率 (其中 )下的平均放电率 (即平均放电率)和电压 的修正项。尽管恢复电流 和突触动力学 还有其他状态,但我们只从放电率和电压状态学习,并且仅基于它们的数据进行训练,因为这些通常是实验中可获得的数据。将训练限制在这两个宏观状态并不会使后续的实验记录成为训练集;相反,它反映了群体记录所强加的可观测性约束。
PEM-UDE方法以良好的精度捕捉了不同稀疏程度下网络的平均场动力学[图4(e)]。这种精度不包括初始条件;这将在第三部分(讨论)的局限性部分进一步讨论。该方法确实捕捉到了网络振荡主导频率随稀疏性变化的偏移[图4(f)]。密集连接区域通常在α频带有一个主频峰,而最稀疏的网络在θ频带表现出主峰。这是学习动力学的一个自然特征,并与实验数据定性一致。海马节律通常在θ频段,并且源于大脑中密集连接的区域[68],而α频段通常源于稀疏连接的皮层区域,然后传播到整个大脑[69]。因此,PEM-UDE学习到的动力学通过准确捕捉系统的遍历特性,扩展了NGNMM在捕捉平均场尺度上涌现特性(如同步性)的能力。学习到的动力学也能在训练数据之外泛化:当仅在 –0.9上训练时,PEM-UDE仍能恢复低至5%连接性的动力学(见补充材料[58]中的图S4)。
与前两个系统一样,我们通过使用STLSQ将动力学的一种符号形式拟合到UDE网络,扩展了方程(39)中描述的PEM-UDE学习到的动力学,从而给出解的解析形式[方程(40),拟合参数见表IV]。Izhikevich及其相关的二次积分激发型NGNMM在其他参数选择和耦合结构下可能经历倍周期和混沌集体状态[25,70,71]。然而,此处使用的参数集处于稳定的群体振荡状态。因此,神经示例在部分观测的振荡系统中测试了相同的训练框架,并与两个混沌基准系统不同。由于NGNMM基于耦合振子群体的假设推导,群体内同步程度有一个精确表达式,由Kuramoto序参量 给出: (26a) (26b) 对于PEM-UDE学习到的网络动力学,我们观察到一个涌现预测:更密集连接的区域倾向于具有更高的同步程度作为基线状态,而更稀疏连接的区域在没有外部驱动作用时倾向于表现出更低的同步程度[图4(g)]。这是从不同稀疏性网络的降阶方程学习得到的一个意外预测。为了检验这一点,我们检查了来自不同稀疏程度的三个颅内电极记录数据集(图5)。由于这些数据集在物种、解剖区域、电极配置、任务背景和行为状态方面不同,这种比较不能将局部连接概率隔离为因果变量。它反而是一个外部一致性检查:那些已知解剖学暗示不同典型局部连接的区域,表现出由模拟降阶模型预测的静息同步的相同排序。这些数据集包括大鼠内侧内嗅皮层(MEC; ),这是靠近海马的深层区域,具有高区域内连接性[72,73];大鼠前额叶皮层( ),其连接性显著更稀疏[74,75];以及人类眶额皮层( ),它表现出类似的高稀疏度[76,77]。与学习到的降阶方程一致,更深且更密集连接的结构表现出比任一皮层数据集更高的同步性,并且这种模式在个体间保持一致,尽管受试者间存在显著变异。这些结果支持了模型定性预测在模拟集之外的合理性,而将学习项的因果验证留给未来使用直接测量连接性的实验。
图5. 颅内EEG显示各区域间的同步性变化与符号化PEM-UDE衍生的神经质量模型一致。如图4(g)所示,连接性较高区域的典型同步性高于连接性较低的区域。在生物大脑中出现了类似的排序:我们检查了三个数据集。在连接性最高的区域(大鼠内侧内嗅皮层;左列),在实验期间表现出比前额叶皮层(中列)更高的同步性(通过相位锁定值PLV测量)。与静息状态的人类眶额皮层(右列)相比,差异更大,后者在自发活动中整体同步性往往很低。由于所有累积分布(上排)都基于来自许多个体的同步值,因此还显示了每个受试者的同步性分布(下排)以说明个体间的变异性。图由BioRender.com创建。
三、讨论
在这里,我们介绍一种在混沌系统中拟合通用微分方程(UDE)的方法,并用它来学习模拟神经群体活动的降阶方程。该方法使用预测误差法来训练UDE系统,然后通过符号回归[12]提取学习动力学的符号形式。我们证明该方法可以学习经典混沌系统——R?ssler系统[24]的动力学。单独的常规UDE无法有效学习该系统,但组合的PEM-UDE方法能准确捕捉动力学,表明该方法已经识别出动力学的真实函数形式,而非在常规统计边界内进行投影。先前在UDE方面的工作集中于刚性系统[8]、在随机噪声存在下拟合UDE[78],以及从UDE网络中提取可解释函数[18]。另外,先前使用PEM进行参数优化的工作已显示其在拟合动力学系统参数方面的效用[19,56],甚至能优于神经网络的动力学拟合[79,80]。我们的方法贡献在于在混沌轨迹的UDE训练中使用预测误差反馈,随后恢复并对显式残差场进行自主验证;观测器、UDE和符号回归本身都是既定方法。该方法在常规UDE方法失败的地方取得成功,因为它平滑了UDE必须优化的参数景观。在第二节B部分所述的噪声自由可实现性和唯一性假设下,这种变形保持了数据一致的零损失集;在存在噪声、模型误设或不完全可观测性的情况下,分析反而识别出限制该陈述的条件和权衡。我们还将该方法扩展到来自电路[60,81]的第二个混沌系统,并表明即使在存在显著观测噪声的情况下,PEM-UDE方法也能恢复未知状态的真实函数形式。通过采用从PEM-UDE动力学学习到的函数形式,并应用更先进的优化技术来充分利用在不同状态下模拟系统的能力(例如,基于模拟的推断[62]),可以获得这些恢复动力学参数的更准确估计。
第二节B部分的解析结果解释了这些经验观察,同时界定了它们的适用范围。对于噪声自由可实现数据,预测误差反馈使零损失集不变,因为反馈在数据一致的流形上消失[方程(12)]。对于具有标量增益的完全观测,它将线性化误差指数偏移了 [方程(14)],因此 是收缩的局部充分条件。部分观测额外要求将未观测的不稳定方向耦合到测量状态。对于噪声或模型误设,同样的分析表明为什么增益必须保持适中:强反馈可能掩盖自主模型误差。方程(23)和(24)的偏差-方差论证启发了此处使用的紧凑网络,但本身并不保证在采样吸引子之外的不可辨识性。
观测器项相对于UDE求解和微分循环的成本较小:它在求解器时间步长上增加了测量状态的插值和对右侧的向量加法,而主要计算仍然由ODE积分、灵敏度计算和网络优化构成。因此,我们将PEM-UDE视为常规UDE训练的稳定化方法,而非竞争性的状态估计方法。为混沌辨识开发的其他方法,包括多重打靶、弱形式目标、影子化损失和滤波辅助符号回归,是自然的补充基线,并且通常可以与同一个观测器正则化UDE结合使用。最近的SINDy类型工作也通过结合线性多步近似、集成自适应库和内积稀疏回归[82],提高了对噪声数据的鲁棒性。我们的比较建立了相对于相同UDE和直接STLSQ在测试条件下的性能;它们并未确立在所有混沌辨识方法中的最先进性能。对更高维、刚性和多尺度系统的系统性跨方法基准仍然是未来的工作。
在证明PEM-UDE方法可以学习混沌系统动力学之后,我们接着考虑如何在生物物理可行约束下学习神经元系统中的降阶控制方程。这些约束有两类:一类是对可用数据的实验性约束,要求只使用能够在体内可行获取的数据源(放电率、膜电压);另一类是神经生物学约束,其中中间神经元连接的稀疏性必须被精确建模。实验性约束值得牢记,因为它们使方法更直接可应用,但方法论的进步往往能随时间放宽这些约束。生物约束,尤其是连接密度,对于精确建模更为重要,因为它们决定了大脑内可能的信息处理类型[53,54]。我们采用了一种“缺失物理”方法,其中群体活动的初始形式由NGNMM给出[44,45]。NGNMM捕捉了神经元群体的涌现特性,并保留了通常在体外被操纵的微观参数[46,83],但它们依赖于全连接假设,这对大多数脑区来说是非生理性的[44]。已有多种尝试绕开这一限制,通过重新加权突触电导程度或突触电流分布来校正稀疏性,但成功程度不一[50,51,66,67,84]。每种方法都带来了新的假设和自身的局限性。在这里,我们直接从网络本身的观测中,为一个由Izhikevich神经元组成的模拟群体学习NGNMM的降阶闭合模型。以这种方式学到的闭合方程能预测在生物学中直接观察到的涌现特性,例如当网络变得更稀疏时,从θ频段到α频段的转变,这在从密集连接区域(如海马体)到稀疏连接区域(如前额叶皮层)的大脑中是可见的特征[68,69]。这些方程还预测,更密集连接的区域在静息状态下表现出比稀疏网络更高的区域内同步性。为验证这一假设,我们计算了来自三个颅内电极记录数据集(大鼠内嗅皮层、大鼠前额叶皮层和人眶额皮层)的区域内同步性[72,74,76]。这些数据集中观察到的排序与模拟模型的预测一致,但由于连接性未被直接测量,这应被解释为生物学一致性检查,而非因果验证。神经群体结果应在降阶闭合的层面上进行解释。微观模拟使用了Izhikevich神经元,并采用了Chen和Campbell[45]的分析性NGNMM构建中继承的参数值,这些值旨在代表全连接极限下的海马样动力学。随后,学到的项在连接概率变化时校正该平均场系统。这不同于常见的神经记录统计拟合,例如脉冲序列的广义线性或点过程模型[85,86]:那些模型通常更简单,且更适用于编码或解码实验,而PEM-UDE旨在保留一个生物物理层面的群体状态空间,能够保留微观参数和涌现同步性。增加结构的代价是,学到的稀疏性校正不应被解读为独特的生物机制,除非它们与匹配结构连接性或前瞻性变化的微观参数记录进行测试。决定性的生物学测试将是在共同制备中,将记录与测量到的结构连接性或受控扰动配对。该方法也有局限性,我们将在未来工作中解决。非常嘈杂的系统可能难以学习,特别是混沌系统在试验中可能难以平均。还应注意:非物理解可能掩盖物理解[26]。通过PEM在训练期间纳入信息,以及提取动力学的符号形式,有助于缓解这些困难。在包含已知物理的模型中,如果参数不完全准确,PEM可能导致接近但不精确的最小值。这限制了该方法在已知物理模型参数足够准确的情景中的适用性。原则上,可以通过拟合这些参数范围内的模拟来克服,但训练时间会显著增加。该方法最适合于已知数据结构(即某些驱动输入已知)的情况。神经实验通常提供这种精度,但对静息态数据也有很大兴趣[87]。因此,将方法适应于拟合无已知扰动的平稳时间序列,是我们认为对神经科学应用重要的另一扩展。对于神经示例,频率随稀疏性的变化可能取决于所选的神经元参数、突触缩放、图系综、有限尺寸效应、延迟、异质性、外部驱动。因此,我们将θ到α的转变解释为这里研究的模拟Izhikevich网络族的预测,而非所有神经群体的普遍规律。与大鼠内嗅皮层、大鼠前额叶皮层和人眶额皮层的比较同样受到物种、区域、电极、行为、临床差异的限制;未来验证理想情况下应使用共同制备中与测量结构连接性配对的记录。虽然本文展示了该方法在神经科学中的实用性,但仍有几个领域需要进一步开发。未来工作可以纳入包含初始条件优化的损失目标,这有助于推广到现实系统[79,88]。未来工作也可以专注于学习动力学的“黑盒”表示,其中不应用符号回归方法来学习动力学。当复杂性惩罚(即学习表达式中的强制稀疏性)偏向远离正确动力学时,这可能有用。例如,当许多小输入基于其时间模式导致复杂相互作用时,神经网络中呈现的高维表示可能更可取。视觉皮层到前额叶皮层投射的时间序列是这种结构化高维输入空间的一个例子[52]。在这种情况下,不是在初始联合运行后符号化拟合动力学,而是在没有PEM的后续优化运行中,可能导致更可泛化的结果(因为这迫使PEM的贡献为0,而不是让其在优化问题中自然最小化)。包括一种直接在动力学中估计随机过程的方法,如同某些标准UDE方法中所做的那样[78],将为该方法增加另一层泛化。这对于有限或部分观察的神经群体尤其相关:NGNMM中的淬火输入分布通过和出现,但未观察的神经元也能产生随机有限尺寸波动,这些波动不能由这里使用的确定性闭合表示。随机PEM-UDE将允许这些波动作为过程噪声被学习,而不是强迫确定性残差项吸收它们。将PEM-UDE方法扩展到偏微分方程组,也将允许在不同连接密度下直接拟合神经场模型,随着神经科学试图理解大脑中行波的作用[21,89]以及密集与稀疏连接区域之间的关系如何塑造认知过程的流动[53],这越来越受关注。四、方法A. Julia生态系统中的科学机器学习Julia生态系统具有坚实的科学机器学习基础设施[18],我们利用它来构建这里描述的PEM-UDE方法。由于在Julia中设置和训练UDE已在其他地方讨论过(参考文献[8,18]是很好的例子),我们在此仅提供关于求解器选择、网络结构、拟合过程和符号回归的简短技术说明,并请感兴趣的读者参阅附带的代码以获取实现细节。在本工作中,我们默认使用非刚性求解器(通常是Vern7),主要为了速度。只要UDE网络足够大,能够对未知动力学进行合理近似,且系统本身不是非常刚性,非刚性自适应求解器应能提供足够准确的拟合。如果系统已知或强烈怀疑是刚性的,则建议切换为刚性求解器(例如Rodas4),这在提供的代码中很容易实现,因为Julia生态系统提供多种求解器。这种切换可能是必要的,原因在于混沌系统在某些区域往往是刚性的,或者即使系统本身不是刚性的(这里讨论的神经系统在正确的参数空间中可能如此),包含一个通用逼近器通常也会导致刚性问题的出现,至少在初始参数权重调整之前是这样。我们建议在默认使用刚性求解器之前增加通用逼近器网络的规模,因为这对性能有很大影响。我们选择径向基函数(RBF;符号见参考文献[90])作为通用逼近器网络的激活函数,由下式给出: (27)选择RBF既有理论原因也有实践原因。大量文献表明,RBF在通用逼近的背景下很有用[90–92]。通用逼近结果激励了这一函数类,但并未使有限网络规模任意化。对于噪声缺失物理问题,较大的网络并不自动更优:一旦残差场可以被表示,额外容量增加了拟合可能跟随特定噪声实现的方向数量。因此,我们使用紧凑的架构,并将容量视为正则化选择。对于本文展示的网络,我们使用四层深度,对于输入维度和输出维度,每个网络的大小为。由于训练中的PEM部分需要将原始训练数据插值到求解器时间步,我们使用DataInterpolations.jl进行观测之间的简单线性插值[93]。因此,PEM引入的额外数值工作是插值调用和右侧反馈添加;它不添加单独的滤波传递或第二个优化问题。尽管有许多优化器可用于训练UDE,但我们在本文展示的所有结果中默认使用AdamW方法(解耦权重衰减和梯度更新)[94]。项的符号恢复使用STLSQ,通过ADMM[95]求解,如DataDrivenDiffEq.jl[63]中实现。根据参考文献[12],是特征库,权重使得,稀疏性通过LASSO强制执行[96]: (28)由于STLSQ对于大数据集( )比遗传算法方法更快[16],我们默认使用此方法。原则上,训练好的UDE可以生成任意密度的子样本,因此两种方法都可用。B. R?ssler系统的背景R?ssler系统是一组非线性微分方程,最初作为研究混沌动力学的简化模型引入。它只包含一个非线性项,但仍产生复杂的混沌行为。动力学由三维常微分方程组描述: (29a) (29b) (29c) 对于特定的参数值(在本文中为、,),系统产生漏斗形R?ssler吸引子,表现出对初始条件的敏感性、分形结构和倍周期混沌活动。我们选择这作为PEM-UDE组合的初始测试,因为尽管R?ssler系统简单,但具有混沌特性,包括正李雅普诺夫指数、拓扑熵和各种分岔序列,这些已被用于研究动力系统中有序与混沌状态之间的转变。R?ssler系统也被扩展以解析描述某些物理问题(例如机械振动、大气变率)。在本文展示的结果中,我们使用通用逼近器来学习整个第二状态的动力学,并借助观测数据通过PEM技术辅助。因此,学习期间模拟的系统由下式给出: (30a) (30b) (30c) 对于正文中展示的结果,我们使用。这一选择是基于第二节讨论的参数景观的陡峭程度;如果选择合适的求解器,该超参数的其他合理值不会影响所学动力学的准确性。本文的结果基于使用AdamW优化器[94]进行100000个训练轮次。用于训练UDE网络的损失函数是每个时间步的预测状态与真实时间序列之间的均方误差(在训练期间添加少量噪声以增强鲁棒性): (31)对于我们分析中选定的参数,我们通过系统的李雅普诺夫时间来表征混沌程度,为。大多数周围空间中的参数组合具有相似的值( 有 ;有关最大李雅普诺夫指数经验累积分布函数(CDF)的计算,请参见补充材料 [58])。李雅普诺夫指数是使用 ChaosTools.jl 中提供的标准函数计算得出的 [97,98]。C. Petrzela-Polak 电路的背景
Petrzela 和 Polak [81] 引入了一族表现出混沌动力学的模拟电路,其满足三个约束条件:规范形式、结构稳定的吸引子和新颖的动力学。我们使用一个特定的实现,其混沌动力学可以在简单电路中物理实现,并且吸引子可在示波器上观察到(参见文献 [60],我们从中采用了此处讨论的参数和无量纲形式)。对于图 3(a) 所示的电路,描述三个观测状态的方程为:
(32a)
(32b)
(32c)
为了数值模拟的方便,我们使用系统的无量纲形式(即用最少的参数表示相同动力学,以替代物理常数):
(33a)
(33b)
(33c)
有量纲和无量纲系统的所有参数均取自文献 [60],并提供于表 II 中。
PEM-UDE 方法学习得到的方程形式为:
(34a)
(34b)
(34c)
用于训练此 PEM-UDE 组合的损失函数与方程 (31) 相同,以测试系统对观测不良状态的鲁棒性。如果先验知道某个状态的观测值损坏到这种程度,则可能值得修改方程 (31) 中的损失函数,以在训练期间对该状态的误差进行较小加权,或完全舍弃这些观测值。如前一节所示,即使某些状态完全隐藏,PEM-UDE 训练方法也能恢复动力学,前提是对其动力学的初始猜测足够准确。
D. 尖峰神经网络和下一代神经质量模型
下一代神经质量模型依赖于 Ott-Antonsen (OA) 假设,以推导神经元群体的平均场活动 [44,48]。OA 假设的引入是为了避免为耦合振子群体的平均场活动建立完整的 Fokker-Planck 方程,只要满足某些假设条件即可。OA 假设 [48] 假设耦合振子群体存在于热力学极限中(即 ),并且可观测量的概率分布(通常为电压)在相空间中收敛为洛伦兹分布。该推导可以以合理的精度近似其他分布 [50]。OA 假设方法还假设输入到每个神经元的噪声是一个异质的、淬火噪声输入 ,其从一般分布中抽取。电流不能是随机的;否则洛伦兹假设将不再成立。OA 假设进一步假设宏观变量的条件概率是可分离的。这一假设及其局限性已被详细检验 [84],并且已为神经质量模型推导出二阶圆累积近似 [66]。OA 假设还作出了我们在本文中要解决的假设:被近似的网络具有全连接结构。如上所述,这对于大脑的多数区域是非生理学的假设,去除该假设对于开发更准确的神经质量模型是必要的。
对于此处的应用,除了作为训练 UDE 所用群体活动函数形式的初始假设外,我们不需要 OA 假设的完整推导,因此省略大部分推导,仅给出主要要点。此处展示的方程改编自给出完整推导的早期工作 [45]。参数集继承自该工作,旨在表示全连接极限下的海马样 Izhikevich 群体动力学;我们未将这些微观参数调至实验数据集。
遵循该工作,我们模拟一个 Izhikevich 神经元网络,其形式为:
(35a)
(35b)
其中 和 分别是神经元 的电压和恢复电流,网络包含 个神经元。这些神经元的电压活动取决于内部神经元动力学、分布式噪声源 、任何外部驱动电流 以及来自其他神经元连接的突触电流 。突触活动由下式给出:
(36a)
(36b)
因此,每当突触前神经元 放电时,每个神经元都会收到一个脉冲 。如果网络结构满足上述假设,则整个群体的平均场活动可以使用 OA 假设推导,形式为:
(37a)
(37b)
(37c)
(37d)
其中 是群体放电率, 是平均电压, 是平均恢复电流, 是平均突触活动。如早期工作所述,这对于全连接网络的平均场活动是一个准确的表示 [45]。虽然此类模型在其他参数选择下允许混沌区域,但表 III 中的参数产生了此处研究的稳定群体振荡。
表 III. Izhikevich 神经元群体和 NGNMM 的参数,来自文献 [45]。参数 值 解释
0.6215 无量纲 Izhikevich 电压参数
1.2308 突触电导
0 外部电流
1 突触逆转电位
0.0077 Izhikevich 恢复电流参数
?0.0062 Izhikevich 恢复电流参数
0.0189 每次放电的恢复电流更新
2.6 突触时间常数
1.2308 每次放电的突触活动更新
200 的截断值
0.12 用于施加噪声的洛伦兹分布的均值
0.02 用于施加噪声的洛伦兹分布的宽度
对于本文展示的所有模拟,我们模拟了一个包含 1000 个 Izhikevich 神经元的网络,并具有不同程度的网络稀疏性,其中两个神经元形成突触的概率由 给出,如图 4 所示的稀疏性。在一个标准的神经元网络中,随着连接稀疏性的变化,突触电流必须按中位入度 的比例重新缩放,以维持平衡(同步)活动。先前的工作表明,将突触电导按 重新缩放,并对背景电流进行类似缩放,可以在广泛的网络结构(具有洛伦兹、高斯和 Erd?s-Rényi 分布的入度)上提供准确的缩放 [50]。遵循该工作,我们将突触电导按 重新缩放,使得:
(38)
其中 。
与前面几节类似,我们设置 UDE 来学习我们感兴趣状态中的未知项,但有两个不同之处。在此推导中,我们学习第一个和第二个状态的项,因此将 PEM 应用于这两个状态。我们还修改方程 (31) 以仅使用前两个状态的数据进行训练,因为这两个状态通常是实验数据中唯一可观测的。该系统的初始参数见表 III,PEM-UDE 系统的正式描述为:
(39a)
(39b)
(39c)
(39d)
通过 PEM-UDE 方法学习动力学并通过 STLSQ 恢复 UDE 表达式的符号形式后,带有稀疏性修正因子的新 NGNMM 由表 IV 中的参数值和函数形式给出:
表 IV. 在 下通过符号 PEM-UDE 学习到的 NGNMM 参数。参数 值 解释
?0.50 突触对放电率贡献的减少
?0.12 的重新缩放
0.88 电压依赖的突触活动
?5.01 放电率对电压贡献的缩放
9.80 非线性相互作用:放电率-电压
?0.05 非线性相互作用:放电率-突触
?0.27 突触对电压的二次贡献
0.12 的重新缩放
0.06 的重新缩放
(40a)
(40b)
(40c)
(40d)
(40e)
(40f)
此形式仅适用于 50% 连接性时的稀疏性(参见补充材料 [58] 的图 S5),并且需要对不同稀疏性进行额外调整。补充材料 [58] 现在给出了完整的神经拟合序列——训练数据和可观测值、优化器和迭代次数、增益延拓、符号提取、外推测试以及所学项的解释——连同对机制解释的限制。
E. 实验数据采集和预处理
三个实验数据集的采集和预处理在其各自的出版物 [72,74,76] 中有详细描述;我们在此仅简要介绍。大鼠内嗅皮层数据集 [73] 从六只雄性 Long-Evans 大鼠中获取,电极记录是在圆形轨道上长时间自由移动期间进行的。大鼠前额叶皮层数据集 [75] 从四只雄性 Long-Evans 大鼠中获取,经历十次成功学习行为任务的试验,包括任务记录前后的休息期。人类眶额皮层数据集 [77] 是在十名正在接受癫痫评估的神经外科患者进行脑电皮层记录期间获取的。患者进行了一项决策游戏,用于本分析的记录是脑电数据中试验间休息部分。
对于所有数据集,区域内电极间的同步性被计算为相位锁定值 (PLV),其基于来自每个电极的希尔伯特变换信号的角度( )[99,100]。遵循先前的工作 [52],对于记录尖峰时间而非 LFP 的数据集,通过加窗尖峰计数将尖峰数据转换为 LFP 代理,并对此进行 PLV 计算。PLV 由下式给出:
(41)
对所有离散的 和所有电极对 有效。