跳转至

A Novel Schur-Decomposition-Based Weight Projection Method for Stable State-Space Neural-Network Architectures

会议: NeurIPS 2026
arXiv: 2605.14489
OpenReview: https://openreview.net/forum?id=WbOem86LmL
代码: https://codeberg.org/sergiovaneg/SchurSS
领域: 优化/理论
关键词: 状态空间模型, 系统辨识, 舒尔分解, 渐近稳定性, 权重投影

一句话总结

针对非线性动力系统辨识中状态空间神经网络层稳定性保障难题,本文提出基于实舒尔分解的截断式准三角矩阵投影与预分解正交投影方法,在不引入冗余超参数与过度参数化的前提下,以反向传播兼容方式实现严格的舒尔稳定性保证。

研究背景与动机

状态空间模型(State-Space Models, SSM)作为序列建模与动态系统辨识的核心范式,因其天然契合物理系统状态演化规律而在控制理论与机器学习领域重获广泛关注。在构建堆叠式神经网络架构(如结合静态非线性模块的 Hammerstein-Wiener 结构)拟合连续动力学系统时,离散时间状态转移矩阵 \(A\) 的严格渐近稳定性(即舒尔稳定性,要求所有特征值模长均落入复平面单位圆内,即 \(|\lambda_i| \le 1\))是防止模型在长程递归仿真中出现输出发散或数值崩溃的绝对先决条件。

然而,现存的稳定状态矩阵参数化与约束策略面临着显著的理论折中与实用瓶颈。传统的子空间辨识方法虽具备稳定性理论保证,但受限于线性可逆性假设,难以与端到端反向传播优化无缝融合;软惩罚正则化方法极度依赖惩罚系数的精细调优,超参数稍小会导致系统失稳,超参数过大又会使系统状态衰减为零矩阵,且训练收敛极其迟缓;基于线性矩阵不等式(LMI)的硬约束参数化(如 SIMBa)虽能覆盖完整稳定空间,却需要将参数量过度放大至 \(5n_x^2\)(其中 \(n_x\) 为状态维度),引发严重的梯度扇出效应(gradient fan-out),拖慢深层非线性架构的收敛效率;而直接求解黎曼流形上的最近 \(\Omega\)-稳定矩阵双层优化问题单步耗时高达数千秒,根本无法作为训练内循环的投影算子。

本文的核心切入点在于打破“流形全量双层优化代价过高”与“启发式正则化无法严格保稳”的两难困境,提出将实舒尔分解与分块对角最近稳定投影巧妙结合。核心 idea:将状态矩阵在正交基下解耦为准上三角舒尔因子,仅在训练批次后对独立的标量与 \(2 \times 2\) 对角块进行闭式/代数根极值投影,以 \(\mathcal{O}(n_x)\) 的极低截断开销实现严格的舒尔稳定性,并进一步提出避免前向分解计算的预分解正交参数化架构。

方法详解

整体框架

本文针对离散时间线性时不变(LTI)状态转移方程 \(x[k+1] = A x[k] + B u[k]\) 中的状态矩阵 \(A \in \mathbb{R}^{n_x \times n_x}\),设计了两种兼顾计算效率与严格稳定性的权重投影机制。第一种为动态分解投影模式(Schur Proj.),在每个训练批次梯度更新后,计算当前状态矩阵的实舒尔分解 \(A = Z T Z^\top\)(其中 \(Z \in O(n_x)\) 为正交矩阵,\(T\) 为准上三角矩阵),随后固定正交基 \(Z\),仅对准三角矩阵 \(T\) 的对角块实施最近舒尔稳定投影得到 \(\hat{T}\),重构出稳定矩阵 \(\hat{A} = Z \hat{T} Z^\top\);第二种为预分解参数化模式(Schur Built),将状态矩阵直接解耦为可训练的正交因子矩阵 \(Z\) 与准上三角因子矩阵 \(T\),分别通过奇异值分解(SVD)正交投影与分块舒尔投影进行独立更新,避免了反向传播各步昂贵的舒尔分解计算开销。

%%{init: {'flowchart': {'rankSpacing': 24, 'nodeSpacing': 28, 'padding': 6, 'wrappingWidth': 400}}}%%
flowchart TD
    In["输入状态矩阵 A 或可训练因子 (Z, T)"] --> Dec{"参数化模式选择"}
    Dec -->|动态投影 Schur Proj.| SF["实舒尔分解<br/>A = Z T Z^T"]
    Dec -->|预分解 Schur Built| SO["正交因子SVD投影<br/>Z_hat = U V*"]
    SF --> QP["截断式准三角块稳定投影<br/>1x1 截断与 2x2 极值集合检索"]
    SO --> QP
    QP --> Recon["重构稳定状态矩阵<br/>A_hat = Z_hat T_hat Z_hat^T"]
    Recon --> HW["集成至非线性状态空间层<br/>Hammerstein-Wiener 架构反向传播"]

关键设计

1. 截断式准三角块稳定投影:固定正交基下的局部闭式投影

在理论最近 \(\Omega\)-稳定矩阵辨识中,对任意实矩阵 \(A\) 寻找 Frobenius 范数下最近的舒尔稳定矩阵需要在正交流形与准三角流形间进行双层交替迭代优化,单次收敛需耗费成千上万秒。本文发现,若在更新时冻结当前的实舒尔正交基 \(Z\),仅在准上三角因子 \(T\) 上寻找最近稳定准三角矩阵,该优化问题可完全沿对角线解耦为独立的标量子块与 \(2 \times 2\) 块投影问题,即算法 1 所定义的准三角舒尔因子投影。对于对角线上的标量元素 \(T_{i,i} \in \mathbb{R}^{1 \times 1}\),其特征值即自身数值,投影操作简化为饱和截断:

\[\hat{T}_{i,i} = \frac{T_{i,i}}{\max\{1, |T_{i,i}|\}}\]

这一截断式投影在保持 \(\mathcal{O}(n_x)\) 线性复杂度的同时,完全消除了全局流形搜索的死循环隐患,为将稳定投影嵌入批次反向传播提供了工程可行性。

2. 2x2 实分块代数根检索机制:解析边界穷举判别

对于对应复共轭特征值对的 \(2 \times 2\) 对角分块 \(T_{i,i} \in \mathbb{R}^{2 \times 2}\),其舒尔稳定性等价于特征多项式根落于闭单位圆盘内。根据 Routh-Hurwitz 准则及其 Cayley 变换,实矩阵稳定等价于三个线性不等式约束构成的紧凸交集:\(\det(X) \le 1\) 且 \(\pm \mathrm{tr}(X) \le 1 + \det(X)\)。当局部块失稳时,其最近稳定解必位于上述约束边界或其交线上。本文将候选解集显式收敛为包含至多 15 个代数矩阵的离散集合 \(\hat{\mathcal{S}}(T_{i,i})\):

\[\hat{\mathcal{S}}(A) \triangleq \{A, \hat{A}_+, \hat{A}_-\} \cup \hat{\mathcal{A}}_0 \cup \hat{\mathcal{A}}_+ \cup \hat{\mathcal{A}}_- \cup \hat{\mathcal{A}}_*\]

其中,集合 \(\hat{\mathcal{A}}_0\) 与 \(\hat{\mathcal{A}}_*\) 分别对应将矩阵投影到行列式为 1 或迹为 0 的流形,其极值点由四阶多项式求导得到的驻点方程闭式确定:

\[t^4 - \sigma_1 t^3 + \sigma_2 t - 1 = 0\]

投影算法依次遍历 \(\hat{\mathcal{S}}(T_{i,i})\) 中至多 15 个解析解,直接筛选出严格满足模长约束 \(|\lambda_{1,2}| \le 1\) 且距离原矩阵 Frobenius 范数最近的合法分块 \(\hat{T}_{i,i}\)。这种代数级检索无需任何局部梯度迭代,计算极度精确且天然支持 JIT 即时编译。

3. 预分解正交因子参数化:解耦状态权重的零分解更新

虽然动态投影模式将单次批次投影时间缩短了 5 个数量级以上,但其在每个训练 step 依然需要调用一次完整的矩阵实舒尔分解,这在现代自动微分框架中存在算子支持有限(如 JAX 早期缺乏原生 GPU Schur 算子)或底层调度耗时的问题。为此,本文设计了预分解正交架构(Schur Built)。该方案直接分配一个全权重矩阵 \(Z \in \mathbb{R}^{n_x \times n_x}\) 和一个准上三角矩阵 \(T \in \mathbb{R}^{n_x \times n_x}\),总参数量仅为 \(2n_x^2\)(显著低于 SIMBa 的 \(5n_x^2\))。在训练前向或参数同步时,\(T\) 仅需直接应用上述算法 1 投影,而正交因子 \(Z\) 则采用 Golub & Van Loan 的 SVD 闭式投影直接映射至最近正交流形:

\[\hat{Z} = U V^*, \quad \text{其中 } Z = U \Sigma V^\top\]

该设计不仅完全剥离了运行时矩阵特征值分解和 Householder 约化步骤,还将正交矩阵的数值误差严格控制在单精度浮点极限之内(误差量级 \(\approx 10^{-15}\)),实现了极佳的训练速度与梯度平滑度。

损失函数 / 训练策略

模型辨识优化目标采用归一化均方误差(NMSE):

\[\ell(\hat{\Theta}; y, x[0]) = \frac{1}{N} \sum_{k=0}^{N-1} \mathcal{L}\left(y[k], \hat{y}[k; \hat{\Theta}, x[0]]\right)\]

针对非线性架构(Hammerstein-Wiener),端到端优化由 AdamW 算法驱动(初始学习率 \(\eta_0 = 10^{-3}\) 或 \(10^{-2}\))。为防止模型在训练初期陷入边界数值敏感区,初始化时设定谱半径模长上界阈值为 \(0.95\) 或 \(0.90\)。在测试与验证阶段,结合扩展卡尔曼滤波(EKF)与 Rauch-Tung-Striebel 平滑器(RTSs)重构非零初始状态 \(x[0]\),消除非零初始条件对输出响应指标带来的系统性偏差。

实验关键数据

主实验

论文首先在 5 种合成 LTI 系统基准(包含高斯测量噪声 \(\sigma=0.01\)、非零共轭虚根)及 5 个真实物理非线性基准数据集(Silverbox、CED 耦合电驱动、EMPS 机电定位系统、KUKA KR300 工业机器人、CubeSpec 精细瞄准镜 FSM)上进行了系统评测,对比了动态投影法(Schur Proj.)、预分解法(Schur Built)、LMI 稳定参数化(SIMBa)以及谱范数惩罚法(Regularized)。

在真实物理非线性基准测试中,最佳模型测试集归一化均方误差(Test NMSE)与收敛轮数(Epoch Count)对比数据如下(误差后跟标准差):

数据集 状态维数 \(n_x\) 指标 Schur (Proj.) Schur (Built) SIMBa Regularized
Silverbox 2 Test NMSE
收敛 Epochs
\((1.52 \pm 7.86)\times 10^{-2}\)
\(8.36\times 10^3\)
\((1.51 \pm 7.87)\times 10^{-2}\)
\(\mathbf{6.89\times 10^3}\)
\(\mathbf{(1.50 \pm 7.97)\times 10^{-2}}\)
\(1.40\times 10^4\)
\((1.53 \pm 7.49)\times 10^{-2}\)
\(1.05\times 10^4\)
CED 3 Test NMSE
收敛 Epochs
\(\mathbf{(3.30 \pm 7.40)\times 10^{-2}}\)
\(1.14\times 10^4\)
\((3.70 \pm 8.79)\times 10^{-2}\)
\(\mathbf{7.83\times 10^3}\)
\((4.12 \pm 9.92)\times 10^{-2}\)
\(2.65\times 10^4\)
\((4.08 \pm 8.49)\times 10^{-2}\)
\(1.75\times 10^4\)
EMPS 4 Test NMSE
收敛 Epochs
\((3.21 \pm 3.24)\times 10^{-3}\)
\(3.16\times 10^4\)
\(\mathbf{(3.07 \pm 4.51)\times 10^{-3}}\)
\(1.33\times 10^4\)
\((4.86 \pm 5.78)\times 10^{-3}\)
\(\mathbf{1.31\times 10^4}\)
\((2.50 \pm 5.18)\times 10^{1}\)
\(4.78\times 10^4\)
Robot 12 Test NMSE
收敛 Epochs
\((6.63 \pm 12.3)\times 10^{-1}\)
\(1.25\times 10^4\)
\(\mathbf{(5.76 \pm 10.9)\times 10^{-1}}\)
\(\mathbf{1.20\times 10^4}\)
\((7.25 \pm 17.5)\times 10^{-1}\)
\(2.13\times 10^4\)
\((7.40 \pm 17.5)\times 10^{-1}\)
\(2.95\times 10^4\)
FSM 28 Test NMSE
收敛 Epochs
\((4.76 \pm 9.31)\times 10^{-2}\)
\(2.91\times 10^4\)
\((3.96 \pm 7.76)\times 10^{-2}\)
\(\mathbf{2.33\times 10^4}\)
\(\mathbf{(3.62 \pm 6.68)\times 10^{-2}}\)
\(2.68\times 10^4\)
\((3.99 \pm 7.82)\times 10^{-2}\)
\(4.54\times 10^4\)

消融与微观基准分析

为验证截断式准三角投影(Algorithm 1)相比全局双层流形优化(Noferini & Poloni [9])的近似保真度与计算加速比,论文评测了不同矩阵维度 \(n\) 下的归一化谱半径误差(NSSR)、归一化 Frobenius 范数误差(NSFE)与计算耗时:

矩阵阶数 \(n\) 投影方法 NSFE(矩阵范数偏差) NSSR(特征值匹配偏差) 违例率 MSVR 单次投影时间 [s]
\(n=10\) 双层流形优化 (Bi-level)
截断块投影 (Algorithm 1)
\(4.08\times 10^{-1}\)
\(6.67\times 10^{-1}\)
\(5.57\times 10^{-1}\)
\(\mathbf{4.06\times 10^{-1}}\)
\(2.18\times 10^{-11}\)
\(2.43\times 10^{-6}\)
\(1.02\times 10^3\)
\(\mathbf{2.43\times 10^{-6}}\)
\(n=20\) 双层流形优化 (Bi-level)
截断块投影 (Algorithm 1)
\(4.55\times 10^{-1}\)
\(\mathbf{3.42\times 10^{-1}}\)
\(6.83\times 10^{-1}\)
\(\mathbf{5.54\times 10^{-1}}\)
\(9.10\times 10^{-9}\)
\(1.36\times 10^{-6}\)
\(1.24\times 10^3\)
\(\mathbf{9.49\times 10^{-5}}\)
\(n=50\) 双层流形优化 (Bi-level)
截断块投影 (Algorithm 1)
\(4.42\times 10^{-1}\)
\(\mathbf{3.24\times 10^{-1}}\)
\(7.48\times 10^{-1}\)
\(\mathbf{6.13\times 10^{-1}}\)
\(1.36\times 10^{-3}\)
\(5.66\times 10^{-2}\)
\(8.20\times 10^3\)
\(\mathbf{4.76\times 10^{-5}}\)
\(n=100\) 双层流形优化 (Bi-level)
截断块投影 (Algorithm 1)
\(4.85\times 10^{-1}\)
\(\mathbf{3.73\times 10^{-1}}\)
\(7.61\times 10^{-1}\)
\(\mathbf{6.28\times 10^{-1}}\)
\(2.65\times 10^{-2}\)
\(4.17\times 10^{-1}\)
\(7.05\times 10^4\)
\(\mathbf{7.03\times 10^{-5}}\)

关键发现

  • 计算加速比极其显著:截断式准三角投影耗时在各维度下均在 \(10^{-5}\) 秒量级,相比全量流形优化的 \(10^3 \sim 7\times 10^4\) 秒实现了高达 8 甚至 9 个数量级的加速,且特征值匹配距离(NSSR)反而优于流形优化提前截断所带来的病态解。
  • 收敛速度大幅领先:在所有真实物理非线性基准上,Schur (Built) 达到最佳验证损失所需的训练 Epoch 数均显著低于 SIMBa 与 Regularized 方法(例如在 CED 任务上 Schur Built 仅需 7,830 轮,而 SIMBa 耗费 26,500 轮,Regularized 耗费 17,500 轮),直接验证了减少状态参数量(\(2n_x^2\) 对比 \(5n_x^2\))能有效缓解梯度扇出效应。
  • 极端动力学下的鲁棒性断崖:在含有纯积分器与复杂非线性摩擦力的 EMPS 数据集上,由于输入包含高频剧烈扰动,Regularized 软约束彻底失效,测试 NMSE 恶化至 \(25.0\)(发散失败),而基于舒尔投影的方法均稳定收敛至 \(3\times 10^{-3}\)。

亮点与洞察

  • 正交基冻结诱导子空间局部凸化:通过将最近矩阵投影约束在当前正交分解的子空间内,原本高度非凸的舒尔稳定空间瞬间解耦为一系列低维相互独立的解析多项式根求解问题,以微乎其微的理论最优性折中换取了工业级反向传播可用性。
  • 参数紧凑性赋能深层优化:相比 SIMBa 为保障稳定性构建的 \(5n_x^2\) 参数过度扩张,预分解 Schur 仅用 \(2n_x^2\) 权重,原生消除了额外中间变量的链式法则放大效应,使得深层非线性系统辨识训练兼具线性层的紧凑与非线性层的灵敏。
  • 普适迁移于序列状态空间架构:该准三角对角分块投影算法完全独立于特定输入输出维度,可即插即用作为 S4、Mamba 或各类递归神经网络(RNN)内部稠密/准稠密循环转移矩阵的投影层。

局限与展望

  • 任意特征值模长阈值约束泛化受限:当前闭式代数投影严格依赖单位圆边界(\(|\lambda| \le 1\))导出的 Routh-Hurwitz 根多项式判据;若需施加严格小于 1 的阻尼紧缩界(如 \(|\lambda| \le \rho < 1\)),需重新推导代数流形驻点方程。
  • GPU 原生 Schur 分解生态支持滞后:主流深度学习框架(PyTorch/JAX)中 GPU 端实舒尔分解算子优化不充分,高维状态在小批量场景下 CPU 吞吐反而更高,阻碍了超大规模语言建模中的吞吐扩展。
  • 高阶刚性系统的数值稳定性界限:在极端高阶系统(\(n_x > 100\))且伴随病态条件数时,单精度浮点累计舍入可能导致微小的复特征值模长溢出,需引入自适应收缩安全裕度 \(\epsilon\)。

相关工作与启发

  • vs SIMBa (Di Natale et al., 2024):SIMBa 利用 LMI 参数化消除显式约束,支持任意阻尼上界,但代价是 \(5n_x^2\) 的权重膨胀与收敛轮数倍增;本文方法仅需 \(n_x^2\) 或 \(2n_x^2\) 参数,训练收敛速度提升 1.5 到 3 倍。
  • vs 谱正则化 (Bemporad, 2025):正则化方法单步计算轻量但无法严格保证渐近稳定,极度依赖惩罚系数超参微调,且在强非线性物理系统(如 EMPS)中易彻底发散;本文方法提供硬性保证且无超参数。
  • vs 最近稳定矩阵黎曼优化 (Noferini & Poloni, 2021):原算法依赖交替流形优化求解全局最优,单次耗时千秒;本文证明冻结正交因子单步投影已足够为梯度下降提供精确的稳态引导。

评分

  • 新颖性: ⭐⭐⭐⭐☆ 首次将实舒尔分解最近稳定投影简化并重塑为端到端反向传播兼容的状态空间稳定算子。
  • 实验充分度: ⭐⭐⭐⭐⭐ 涵盖 5 组不同状态维度的合成系统与 5 个严苛的物理非线性真实数据集,微观对比与系统测试完备。
  • 写作质量: ⭐⭐⭐⭐⭐ 数学推导严密自洽,附录提供了完整的 2x2 代数根证明及 GPU 自定义算子细节。
  • 价值: ⭐⭐⭐⭐☆ 为神经状态空间模型在严苛物理控制与经典系统辨识场景的落地提供了坚实可靠的硬约束基石。