Paper deep dive

Residual Koopman Spectral Profiling for Predicting and Preventing Transformer Training Instability

Bum Jun Kim, Shohei Taniguchi, Makoto Kawano, Yusuke Iwasawa, Yutaka Matsuo

Year: 2026Venue: arXiv preprintArea: cs.LGType: PreprintEmbeddings: 77

Intelligence

Status: succeeded | Model: Gemma-4-26B-A4B | Prompt: intel-v1 | Confidence: 93%

Last extracted: 7/20/2026, 9:33:03 AM

Summary

The paper introduces Residual Koopman Spectral Profiling (RKSP), a method to predict transformer training instability at initialization by analyzing the spectral properties of layer-wise residual dynamics using whitened Dynamic Mode Decomposition (DMD). The central metric, near-unit spectral mass, quantifies the risk of divergence. The authors also propose Koopman Spectral Shaping (KSS), a regularization technique that stabilizes training by reshaping these spectra, significantly reducing divergence rates and allowing for higher learning rates in normalization-free regimes.

Entities (12)

Transformer · model-architecture · 99%Residual Koopman Spectral Profiling · method · 98%Koopman Spectral Shaping · method · 97%Near-unit spectral mass · metric · 96%Dynamic Mode Decomposition · algorithm · 95%CIFAR-10 · dataset · 92%WikiText-103 · dataset · 92%GPT-2 · model · 90%LLaMA-2 · model · 90%KAN · architecture · 85%Mamba-style SSMs · architecture · 85%MoE · architecture · 85%

Relation Signals (8)

Residual Koopman Spectral Profiling → predicts → Transformer Training Instability

confidence 98% · RKSP predicts divergence at initialization... captures instability risk.

Koopman Spectral Shaping → prevents → Divergence

confidence 97% · when RKSP flags high risk, turning on KSS successfully prevents divergence.

Near-unit spectral mass → quantifies → Instability Risk

confidence 96% · the near-unit spectral mass, quantifies the fraction of modes concentrated near the unit circle, which captures instability risk.

Residual Koopman Spectral Profiling → uses → Dynamic Mode Decomposition

confidence 95% · RKSP extracts Koopman spectral features by applying whitened dynamic mode decomposition to layer-wise residual snapshots.

Koopman Spectral Shaping → enables → Higher Learning Rates

confidence 90% · enables learning rates that are 50% to 150% higher.

Residual Koopman Spectral Profiling → generalizesto → GPT-2

confidence 90% · These findings generalize to ... pretrained language models, including GPT-2

Residual Koopman Spectral Profiling → generalizesto → LLaMA-2

confidence 90% · These findings generalize to ... LLaMA-2 up to 7B

Residual Koopman Spectral Profiling → appliesto →

Cypher Suggestions (0)

No Cypher suggestions yet.

Abstract

Abstract:Training divergence in transformers wastes compute, yet practitioners discover instability only after expensive runs begin. They therefore need an expected probability of failure for a transformer before training starts. Our study of Residual Koopman Spectral Profiling (RKSP) provides such an estimate. From a single forward pass at initialization, RKSP extracts Koopman spectral features by applying whitened dynamic mode decomposition to layer-wise residual snapshots. Our central diagnostic, the near-unit spectral mass, quantifies the fraction of modes concentrated near the unit circle, which captures instability risk. For predicting divergence across extensive configurations, this estimator achieves an AUROC of 0.995, outperforming the best gradient baseline. We further make this diagnostic actionable through Koopman Spectral Shaping (KSS), which reshapes spectra during training. We empirically validate that our method works in practice: RKSP predicts divergence at initialization, and when RKSP flags high risk, turning on KSS successfully prevents divergence. In the challenging high learning rate regime without normalization layers, KSS reduces the divergence rate from 66.7% to 12.5% and enables learning rates that are 50% to 150% higher. These findings generalize to WikiText-103 language modeling, vision transformers on CIFAR-10, and pretrained language models, including GPT-2 and LLaMA-2 up to 7B, as well as emerging architectures such as MoE, Mamba-style SSMs, and KAN.

PDF

Open source PDF →Open local PDF →

Trouble viewing inline? Open PDF directly →

Full Text

75,672 characters extracted from source content.

Expand or collapse full text

Residual Koopman Spectral Profiling for Predicting and Preventing Transformer Training Instability Bum Jun Kim Graduate School of Engineering, The University of Tokyo, Japan Shohei Taniguchi Graduate School of Engineering, The University of Tokyo, Japan Makoto Kawano Graduate School of Engineering, The University of Tokyo, Japan Yusuke Iwasawa Graduate School of Engineering, The University of Tokyo, Japan Yutaka Matsuo Graduate School of Engineering, The University of Tokyo, Japan Abstract Training divergence in transformers wastes compute, yet practitioners discover instability only after expensive runs begin. They therefore need an expected probability of failure for a transformer before training starts. Our study of Residual Koopman Spectral Profiling (RKSP) provides such an estimate. From a single forward pass at initialization, RKSP extracts Koopman spectral features by applying whitened dynamic mode decomposition to layer-wise residual snapshots. Our central diagnostic, the near-unit spectral mass, quantifies the fraction of modes concentrated near the unit circle, which captures instability risk. For predicting divergence across extensive configurations, this estimator achieves an AUROC of 0.995, outperforming the best gradient baseline. We further make this diagnostic actionable through Koopman Spectral Shaping (KSS), which reshapes spectra during training. We empirically validate that our method works in practice: RKSP predicts divergence at initialization, and when RKSP flags high risk, turning on KSS successfully prevents divergence. In the challenging high learning rate regime without normalization layers, KSS reduces the divergence rate from 66.7% to 12.5% and enables learning rates that are 50% to 150% higher. These findings generalize to WikiText-103 language modeling, vision transformers on CIFAR-10, and pretrained language models, including GPT-2 and LLaMA-2 up to 7B, as well as emerging architectures such as MoE, Mamba-style SSMs, and KAN. 1 Introduction Transformers [vaswani2017attention] exhibit unpredictable training dynamics despite stabilization techniques such as layer normalization [ba2016layernorm] and careful initialization. In particular, exploding-gradient instability remains a known failure mode in deep networks [pascanu2013difficulty]. This unpredictability in training divergence wastes compute: practitioners often discover divergence only after expensive runs have begun. Addressing this problem requires a calibrated, initialization-time risk estimate, which enables early pruning of risky configurations. Such an estimate should output a probability that matches empirical frequency rather than a heuristic score. To provide such estimates, we propose Residual Koopman Spectral Profiling (RKSP), which views transformer layers as discrete-time dynamical systems. From a single forward pass, RKSP applies whitened dynamic mode decomposition (DMD) [tu2014dmd] to estimate local linear operators that approximate the residual stream evolution 0→1→⋯→Lh_0 _1→·s _L. The resulting layer-wise spectra provide a compact, predictive summary of stability. Intuitively, when a layer’s local linearization is close to normal, eigenvalues near the unit circle imply near-isometric propagation and therefore weak damping of signals and gradients. In high learning rate regimes, such weak damping can leave perturbations and optimization noise to be less attenuated, increasing divergence risk; conversely, strongly contractive spectra tend to be more stable but may cause rapid gradient decay. We summarize this trade-off with the near-unit mass M≈1M_≈ 1, the fraction of modes near the unit circle, together with a separate measure of non-normality; Section 4 provides a theoretical explanation for this behavior. The near-unit mass M≈1M_≈ 1 exhibits a strong instability signal; using the monotone risk score M≈1M_≈ 1 yields an Area Under the Receiver Operating Characteristic curve (AUROC) of 0.995 for divergence prediction across various normalization strategies, tasks, and architectures. To make this signal further actionable, we introduce Koopman Spectral Shaping (KSS), which reshapes spectra during training to reduce divergence. Our contributions are as follows. • We propose RKSP, which estimates layer-wise Koopman spectra at initialization via whitened DMD (Section 3.2). • We propose KSS, which reshapes spectra during training, preventing instability and permitting larger learning rates (Section 3.3). • We provide theoretical bounds linking near-unit mass, near-normality, and the trade-offs between instability and expressivity (Section 4). • We validate RKSP across various normalization strategies, datasets, and architectures, including language models, vision transformer (ViT) models, and emerging neural network families (Section 5 and Appendix G). 2 Background A transformer with L layers updates its residual stream according to the residual formulation [he2016resnet]: ℓ+1=ℓ+fℓ(ℓ;θℓ)=Fℓ(ℓ),h_ +1=h_ +f_ (h_ ; _ )=F_ (h_ ), (1) where fℓf_ comprises self-attention or multi-layer perceptron sub-layers. This formulation reveals that each layer transition defines a discrete-time dynamical system. Also, residual networks can be viewed as discretizations of continuous-time dynamical systems, motivating stability analysis from an ordinary differential equation perspective [haber2017stable, chen2018neuralode]. Indeed, the local linearization of this system characterizes its stability. 2.1 Koopman Operator Theory and DMD Primer Consider a discrete-time dynamical system t+1=F(t)x_t+1=F(x_t) on state space ⊆ℝdX ^d. The Koopman operator :ℱ→ℱK:F acts on observable functions g:→ℂg:X via composition with the dynamics [koopman1931hamiltonian, mezic2005spectral, mezic2013analysis, tu2014dmd]: (g)()≜g(F()).(Kg)(x) g(F(x)). (2) For a nonlinear F, the operator K is infinite-dimensional but remains linear regardless of the nonlinearity in F. The spectral properties of K—its eigenvalues λj\ _j\ and eigenfunctions ϕj\ _j\—encode the intrinsic timescales and geometric structure of the dynamics: ϕj _j =λjϕj, = _j _j, (3) therefore ϕj(Fn()) _j(F^n(x)) =λjnϕj(),for all n≥0. = _j^n _j(x),for all n≥ 0. The spectrum admits direct interpretation: modes with |λj|>1| _j|>1 grow exponentially, modes with |λj|<1| _j|<1 decay exponentially, and modes with |λj|=1| _j|=1 persist or oscillate. The argument arg⁡(λj) ( _j) gives the oscillation frequency of mode j [rowley2009spectral]. In a layer-wise, non-autonomous setting, each layer has its own distinct operator, so each ^ℓ A_ must be estimated separately. Thus, the modulus |λ(^ℓ)||λ( A_ )| indicates a local, per-layer expansion or contraction tendency. DMD approximates these layer-wise Koopman operators from finite data [schmid2010dmd, tu2014dmd, kutz2016dmdbook]. 2.2 Related Work Koopman methods in machine learning. Existing Koopman methods approximate operators via DMD and its extensions or learn linearizing transformations for dynamics prediction [lusch2018deeplin, williams2015edmd, korda2018linear]. These methods primarily target representation learning. Edge of chaos. The edge of chaos hypothesis links optimal trainability to critical initialization and signal propagation regimes [schoenholz2017deepinfo, poole2016chaos, pennington2017dynamical], where signals neither explode nor vanish. Neural tangent kernel. The neural tangent kernel characterizes gradient descent in the infinite-width limit and yields kernel-regression behavior [jacot2018ntk]. In particular, linearizing the network around its initialization makes the kernel essentially constant, so training reduces to regression with this fixed kernel. Mean-field theory and μ . Mean-field theory and Maximal Update Parameterization (μ ) enable hyperparameter transfer across scales through asymptotic analysis [schoenholz2017deepinfo, yang2022tensorprograms]. Normalization-free training. Normalization-free residual networks can be stabilized with Fixup initialization [zhang2019fixup]. ReZero trains deep residual networks and transformers without normalization by introducing a residual scaling parameter initialized to zero, so the network starts near an identity map [bachlechner2020rezero]. Normalizer-Free Networks replace normalization with scaled activations and adaptive gradient clipping, enabling stable large-scale training without normalization layers [brock2021nfnet]. We discuss Fixup and ReZero-style identity initialization through the RKSP lens in Appendix J. Positioning of this study. RKSP uses Koopman spectra as an initialization-time diagnostic: from a single forward pass, we estimate layer-wise operators and predict divergence risk. Unlike Koopman representation learning for forecasting [lusch2018deeplin, williams2015edmd, korda2018linear], our near-unit mass M≈1M_≈ 1 provides a measurable handle on critical signal propagation [schoenholz2017deepinfo, poole2016chaos, pennington2017dynamical] and captures instabilities beyond neural tangent kernel analyses [jacot2018ntk]; we make it actionable via KSS and validate it in no-normalization regimes [zhang2019fixup, brock2021nfnet]. 3 Method 3.1 Problem Setup Divergence Definition A training run is marked as diverged if, at any step, the loss exceeds 50.0 or the gradient norm exceeds 500.0. This criterion defines the binary label D∈0,1D∈\0,1\ used throughout our experiments. Divergence Prediction Task Given a model architecture, normalization strategy, optimizer choice, and dataset, we collect N residual-stream snapshots at initialization from a single forward pass and compute the spectral profile S. A probabilistic predictor then maps S to P(D=1∣)P(D=1 ), estimating divergence risk before training begins. 3.2 Residual Koopman Spectral Profiling Algorithm 1 summarizes our RKSP procedure. The algorithm computes DMD for each layer to obtain its spectral profile: ρℓ _ denotes the spectral radius of ^ℓ A_ , κℓ _ the eigenvector condition number of its eigenbasis, and ℓK_ the Kreiss constant (Appendix B). Algorithm 1 Residual Koopman Spectral Profiling 1:Model ℳM with L layers; dataset D; the number of samples N 2:Spectral profile =(M≈1ℓ,ρℓ,κℓ,ηnlℓ,ℓ)ℓ=0L−1S=\(M_≈ 1 , _ , _ , _nl ,K_ )\_ =0^L-1 3:Collect residuals: for each batch ∈x , store ℓ()ℓ=0L\h_ (x)\_ =0^L 4:for ℓ=0,…,L−1 =0,…,L-1 do 5: Form snapshot matrices ℓ,ℓ∈ℝd×NX_ ,Y_ ^d× N 6: Whitening: ~ℓ=^ℓ−1/2(ℓ−¯ℓ) X_ = _ ^-1/2(X_ - X_ ), ~ℓ=^ℓ−1/2(ℓ−¯ℓ) Y_ = _ ^-1/2(Y_ - Y_ ) 7: DMD: ^ℓ=~ℓ~ℓ† A_ = Y_ X_ 8: Eigendecomposition: ^ℓ=ℓℓℓ−1 A_ =V_ _ V_ ^-1 9: Compute M≈1ℓ,ρℓ,κℓM_≈ 1 , _ , _ , nonlinearity ηnlℓ _nl , Kreiss ℓK_ 10:end for 11:Aggregate mean, max, min, std across layers 12:return S 3.2.1 Snapshot Construction RKSP applies DMD to each layer transition ℓ→ℓ+1h_ _ +1 as described in Eq. 1. Let ii=1N\x_i\_i=1^N denote the N inputs used to form the snapshots. For layer ℓ , define the paired residual vectors i(ℓ)=ℓ(i)x_i^( )=h_ (x_i) and i(ℓ)=ℓ+1(i)y_i^( )=h_ +1(x_i). The snapshot matrices are ℓ _ =[1(ℓ),…,N(ℓ)], =[x_1^( ),…,x_N^( )], (4) ℓ _ =[1(ℓ),…,N(ℓ)]∈ℝd×N, =[y_1^( ),…,y_N^( )] ^d× N, so each column pair corresponds to the same sample. Unlike standard DMD, which uses time-shifted trajectories, RKSP pairs columns across different samples. Concretely, column i in X is the residual snapshot ℓ(i)h_ (x_i) and column i in Y is the corresponding next-layer snapshot ℓ+1(i)h_ +1(x_i) for the same sample ix_i; columns index independent samples, not time steps of a single trajectory. For each layer, DMD yields a local linear approximation ^ℓ A_ whose spectrum characterizes the dynamics at that depth. To quantify how well this linear approximation fits the data, we define the nonlinearity ratio in whitened coordinates: ηnl(ℓ)≔‖~ℓ−^ℓ~ℓ‖F‖~ℓ−~ℓ‖F+εnl, _nl( ) \| Y_ - A_ X_ \|_F \| Y_ - X_ \|_F+ _nl, (5) where εnl>0 _nl>0 is a small constant that prevents division by zero. This ratio ηnl _nl normalizes the fit error by the update magnitude. When the residual update ‖~ℓ−~ℓ‖F \| Y_ - X_ \|_F is tiny, ηnl _nl can be large even for small absolute errors. We therefore use ηnl _nl primarily as a DMD reliability flag rather than as a pure measure of nonlinearity. 3.2.2 Whitened DMD and Reliability Filtering DMD approximates the Koopman operator from data snapshots. Given paired snapshot matrices =[1,…,N]∈ℝd×NX=[x_1,…,x_N] ^d× N and =[1,…,N]∈ℝd×NY=[y_1,…,y_N] ^d× N, DMD solves for the optimal linear operator: ^DMD=argmin∈ℝd×d‖−‖F2=†, A_DMD= *argmin_A ^d× d \|Y-AX \|_F^2=YX , (6) where †X denotes the Moore-Penrose pseudoinverse. To ensure scale-invariance and numerical stability, we apply X-based zero-phase component analysis whitening [kessy2018whitening]: ~ X =^X−1/2(−¯⊤), = _X^-1/2(X- X1 ), (7) ~ Y =^X−1/2(−¯⊤), = _X^-1/2(Y- Y1 ), (8) ^X _X =1N−1(−¯⊤)(−¯⊤)⊤+ϵ = 1N-1(X- X1 )(X- X1 ) + (9) where ¯=1N∑i=1Ni X= 1N _i=1^Nx_i, ¯=1N∑i=1Ni Y= 1N _i=1^Ny_i, and ϵ>0ε>0 ensures invertibility. The same whitening matrix is applied to both X and Y, so the regression operates within a single, X-normalized coordinate system. This whitening step ensures cross-model comparability and yields coordinate-invariant spectral estimates. From the whitened data, we form the DMD operator ^=~~† A= Y X and compute its eigendecomposition ^=−1 A=V V^-1. To report the eigenvector condition number κ()=‖2‖−1‖2κ(V)= \|V \|_2 \|V^-1 \|_2, we first normalize each right eigenvector to have unit Euclidean norm. This normalization fixes the otherwise arbitrary scaling of V and makes κ()κ(V) reproducible. To identify spurious eigenvalues, we apply residual DMD (ResDMD) reliability filtering [colbrook2023resdmd]. For each eigenvalue λj _j with left eigenvector ju_j satisfying j∗^=λjj∗u_j^* A= _ju_j^*, we compute the per-mode residual: rj=‖j∗(~−λj~)‖2‖j∗~‖2+εr,r_j= \|u_j^*( Y- _j X) \|_2 \|u_j^* X \|_2+ _r, (10) where εr>0 _r>0 prevents division by zero. Eigenvalues with rj>τr_j>τ are flagged as unreliable and potentially spurious; we use a default threshold of τ=0.1τ=0.1 to filter them out. 3.2.3 Spectral Mass Now, we define our criterion for divergence prediction. Definition 1 (Spectral Mass Partition). For DMD eigenvalues Λ=λjj=1m =\ _j\_j=1^m of an operator A, with m the number of eigenvalues used, we define the following bins; they are disjoint provided δc≥ϵn _c≥ _n: M>1(^) M_>1( A) ≜1m∑j=1m[|λj|>1+ϵu], 1m _j=1^m1[| _j|>1+ _u], (11) M≈1(^) M_≈ 1( A) ≜1m∑j=1m[|λj|∈[1−ϵn,1+ϵu]], 1m _j=1^m1[| _j|∈[1- _n,1+ _u]], (12) M<1(^) M_<1( A) ≜1m∑j=1m[|λj|<1−δc], 1m _j=1^m1[| _j|<1- _c], (13) where the three quantities denote the expansive mass, near-unit mass, and contractive mass, respectively. When δc>ϵn _c> _n, these three bins are not exhaustive; the remaining intermediate mass is Mmid(^)≜1−M>1(^)−M≈1(^)−M<1(^)M_mid( A) 1-M_>1( A)-M_≈ 1( A)-M_<1( A). For example, m=dm=d for full DMD, m=rm=r for rank-r randomized DMD [erichson2019randdmd], or m equals the count remaining after reliability filtering. We use default thresholds ϵu=0.05 _u=0.05, ϵn=0.10 _n=0.10, and δc=0.20 _c=0.20. Figure 1 shows a representative eigenvalue scatter that motivates these bins. Figure 1: Scatter plot of DMD eigenvalues across layers in a pre-layer normalization transformer. The color gradient indicates layer depth; blue is early and red is late. Early layers cluster near the unit circle; late layers exhibit an increased spectral radius. Metric Interpretation. For divergence prediction, we use M≈1M_≈ 1 itself as the scalar score and compute the AUROC against the divergence labels. Note that the expansive mass M>1M_>1 tracks eigenvalues outside the unit circle but does not map one-to-one with empirical divergence. Four factors explain this gap between M>1M_>1 and observed divergence. First, whitening rescales local coordinates, so raw eigenvalue magnitudes differ from unwhitened values. Second, DMD provides only a local linear approximation of the true nonlinear dynamics. Third, non-normal transient growth can trigger instability even when few eigenvalues exceed 1 [trefethen2005spectra]. Additionally, our divergence labels use coarse thresholds on loss or gradient norm, so finite-horizon training within the evaluation window can remain stable despite a nonzero M>1M_>1. These four factors together explain cases like Pre-LN, which shows a nonzero M>1M_>1 but 0% divergence in Table 1. 3.3 Koopman Spectral Shaping While RKSP diagnoses instability, KSS prevents it. KSS adds a differentiable spectral regularizer to the training objective that steers eigenvalues away from the unstable region while reducing excessive near-unit mass to restore damping without causing over-contraction. The total objective becomes ℒtotal=ℒtask+α∑ℓ∈ℒKSSℓ/||L_total=L_task+α _ L_KSS /|S|, where S samples 50% of layers per update. Definition 2 (KSS Regularization Loss). For layer ℓ with randomized DMD eigenvalues λjℓj=1r\ _j \_j=1^r, the KSS loss is ℒKSSℓ _KSS =∑j=1rσ(T(|λjℓ|−τu))⋅softplus(|λjℓ|−τu)2⏟Unstable penalty = _j=1^rσ(T(| _j |- _u))·softplus(| _j |- _u)^2_Unstable penalty (14) +β⋅(mℓsoft−γ)2⏟Near-unit target + β·(m_ ^soft-γ)^2_Near-unit target where σ(⋅)σ(·) denotes the sigmoid function and mℓsoft=1r∑j=1rσ(T(|λjℓ|−τl))⋅σ(T(τu−|λjℓ|)).m_ ^soft= 1r _j=1^rσ(T(| _j |- _l))·σ(T( _u-| _j |)). (15) This term facilitates reducing excessive near-unit mass while preventing it from becoming too small by nudging mℓsoftm_ ^soft toward the target band γ. We use the default hyperparameters T=20T=20, τu=1.05 _u=1.05, τl=0.90 _l=0.90, and γ∈[0.3,0.5]γ∈[0.3,0.5]. Full hyperparameter settings and the practical training recipe appear in Appendix D. 4 Theoretical Analysis We interpret the near-unit mass M≈1M_≈ 1 as an instability score. Two mechanisms support this view. First, when the layer linearization is approximately normal, eigenvalues accurately reflect singular values, so concentration near |λ|≈1|λ|≈ 1 implies near-isometric propagation and weak damping. Second, weak damping allows perturbations and optimization noise to persist across depth; in aggressive optimization regimes, this behavior raises divergence risk. We formalize the energy-preservation statement below and highlight non-normality as a caveat. Theorem 1 (Near-Unit Energy Preservation under Near-Normality). Let ∈ℂd×dA ^d× d be normal with eigenvalues λjj=1d\ _j\_j=1^d. For a unit vector x drawn uniformly on the sphere, ‖22=1d∑j=1d|λj|2.E \|Ax \|_2^2= 1d _j=1^d| _j|^2. (16) If ρ()≤1+ϵuρ(A)≤ 1+ _u and M≈1()M_≈ 1(A) denotes the fraction of eigenvalues with |λj|∈[1−ϵn,1+ϵu]| _j|∈[1- _n,1+ _u], then (1−ϵn)2M≈1()≤‖22≤(1+ϵu)2.(1- _n)^2M_≈ 1(A) \|Ax \|_2^2≤(1+ _u)^2. (17) Hence larger M≈1M_≈ 1 implies more energy-preserving and less damped propagation; this corresponds to higher instability risk. More generally, if A is diagonalizable with =−1A=V V^-1, then the same conclusion holds up to factors κ()±2κ(V)^± 2: (1−ϵn)2κ()2M≈1()≤‖22≤κ()2(1+ϵu)2. (1- _n)^2κ(V)^2M_≈ 1(A) \|Ax \|_2^2≤κ(V)^2(1+ _u)^2. (18) Corollary 2 (Depth-wise Damping and Gradient Flow). Consider a depth-L linearization ℓ+1=ℓℓh_ +1=A_ h_ where each ℓA_ is a normal and ρ(ℓ)≤1+ϵuρ(A_ )≤ 1+ _u. Assume that each layer has an isotropic second moment, that is, [ℓℓ∗]=1d‖ℓ‖22E[h_ h_ ^*]= 1dE \|h_ \|_2^2I for ℓ=0,…,L−1 =0,…,L-1. For example, ℓ‖ℓ‖2 h_ \|h_ \|_2 is uniform on the sphere. Let qℓ≜1d∑j=1d|λjℓ|2q_ 1d _j=1^d| _j |^2 be the average energy gain of layer ℓ . Then ‖L‖22=‖0‖22∏ℓ=0L−1qℓ,E \|h_L \|_2^2=E \|h_0 \|_2^2 _ =0^L-1q_ , (19) and qℓ≥(1−ϵn)2M≈1(ℓ)q_ ≥(1- _n)^2M_≈ 1(A_ ). Thus, a larger M≈1M_≈ 1 reduces the exponential contraction of signals and gradients; in high learning rate regimes, this weaker damping elevates instability risk, although a very small M≈1M_≈ 1 can hurt expressivity. The proof appears in Appendix C. Non-normality caveat. When κ()≫1κ(V) 1, non-normal transient growth can occur even if ρ()≤1ρ(A)≤ 1, and eigenvalues near the unit circle may be perturbation-sensitive. We therefore track κ()κ(V) alongside M≈1M_≈ 1. Appendix B and Appendix C characterize transient growth via the Kreiss theorem, provide normality-related bounds, and formalize the trade-off between instability and expressivity. 5 Experiments 5.1 Experimental Settings Our experiments target Generative Pretrained Transformer (GPT-2)-style transformers [radford2019language, vaswani2017attention] with d∈128,256,512,768,1024d∈\128,256,512,768,1024\ and L∈4,6,8,12,16,24L∈\4,6,8,12,16,24\, spanning 1M to 350M parameters. We compare six normalization strategies: pre-layer normalization (Pre-LN), post-layer normalization (Post-LN), root mean square normalization (RMSNorm) [zhang2019rmsnorm], DeepNorm [wang2022deepnet], sub-layer normalization (SubLN) [xiong2020layernorm], and no normalization (No-Norm). The evaluation tasks include an associative-recall classification task (Appendix D.1); language modeling (LM) tasks including our synthetic LM with next-token prediction (Appendix D.2), WikiText-103 [merity2016wikitext], and OpenWebText-style LM [gokaslan2020openwebtext]; and ViT experiments [dosovitskiy2020vit] on the Canadian Institute for Advanced Research (CIFAR-10) dataset [krizhevsky2009cifar]. We report AUROC for discrimination, with 95% bootstrap confidence intervals (CIs) based on 1000 resamples. Statistical significance of divergence rates is assessed using Fisher’s exact test and run with three random seeds per setting. Full hyperparameters and hardware details appear in Appendix D. 5.2 Prediction at Initialization Main Results: Normalization Comparison. Table 1 reveals three key findings. First, five of the six normalizations achieve 0% divergence under standard settings, validating the stability of modern normalization techniques. Second, No-Norm diverges in 96.4% of runs and also has high near-unit mass, with M≈1=0.80M_≈ 1=0.80; under its relatively low non-normality, with κ()=1.11κ(V)=1.11, this is consistent with M≈1M_≈ 1 acting as an instability indicator. Non-normality and expansive mass still matter, but within comparable regimes, a larger M≈1M_≈ 1 aligns with greater instability, as discussed in Section 4. Third, spectral signatures differ systematically across normalizations: Pre-LN and RMSNorm are more contractive with lower M≈1M_≈ 1, Post-LN retains a higher near-unit structure, and DeepNorm shifts mass into the unstable bin but remains bounded through scaling. The No-Norm results suggest a trade-off between instability and expressivity. Only three of 84 No-Norm runs converged, but their accuracy reached 54.2% ± 24.4%—higher than other methods. Among stable methods, DeepNorm achieves the best balance, with 20.0% accuracy and 0% divergence. Table 1: Comprehensive normalization comparison across different setups. RKSP metrics reveal distinct spectral signatures explaining stability differences. Accuracy is computed from only 3 converged runs out of 84. Statistical significance: No-Norm divergence compared with others, p<10−50p<10^-50 via Fisher’s exact test. AUROC for divergence prediction using M≈1M_≈ 1: 0.995 [95% CI: 0.986 to 1.00]. Measured on the associative-recall task. black Type n Div.% (lower) M≈1M_≈ 1 M>1M_>1 ρ Acc.% (higher) Pre-LN 25 0.00.0 0.16±0.010.16± 0.01 0.54±0.010.54± 0.01 2.29±0.042.29± 0.04 7.5±8.77.5± 8.7 Post-LN 69 0.00.0 0.66±0.020.66± 0.02 0.31±0.010.31± 0.01 7.48±0.127.48± 0.12 1.4±5.81.4± 5.8 RMSNorm 12 0.00.0 0.16±0.010.16± 0.01 0.54±0.010.54± 0.01 2.28±0.062.28± 0.06 13.1±8.613.1± 8.6 DeepNorm 12 0.00.0 0.00±0.000.00± 0.00 1.00±0.001.00± 0.00 3.94±0.213.94± 0.21 20.0±9.420.0± 9.4 SubLN 12 0.00.0 0.13±0.020.13± 0.02 0.54±0.010.54± 0.01 2.23±0.032.23± 0.03 9.0±8.19.0± 8.1 No-Norm 84 96.496.4 0.80±0.020.80± 0.02 0.19±0.010.19± 0.01 1.11±0.011.11± 0.01 54.2±24.4∗54.2± 24.4^* AUROC Analysis for Divergence Prediction. Table 5.2 compares spectral predictors against gradient baselines. We use the monotone risk score M≈1M_≈ 1 for divergence prediction. This achieves an AUROC of 0.995 at initialization, representing a 31% relative improvement over the best gradient-based method with an AUROC of 0.758. The superiority is statistically significant: the 95% CI lower bound of M≈1M_≈ 1 of 0.986 exceeds the upper bounds of all gradient-based methods. Table 2: AUROC for divergence prediction with bootstrap 95% confidence intervals. The M≈1M_≈ 1 CI lower bound of 0.986 exceeds gradient baselines’ upper bounds. Measured on the associative-recall normalization sweep. black AUROC (higher) 95% CI Timing M≈1M_≈ 1 at initialization 0.9950.995 [0.986,1.000][0.986,1.000] Initialization M≈1×log10⁡(κ())M_≈ 1× _10(κ(V)) 0.8730.873 [0.841,0.905][0.841,0.905] Initialization Spectral radius ρ at init 0.8450.845 [0.808,0.882][0.808,0.882] Initialization Eigenvector condition κ()κ(V) at init 0.6870.687 [0.638,0.736][0.638,0.736] Initialization Gradient-Based Baselines Initial gradient norm 0.6210.621 [0.568,0.674][0.568,0.674] After 1 step Gradient norm at step 100 0.6850.685 [0.635,0.735][0.635,0.735] Step 100 Gradient variance over steps 1 to 100 0.7020.702 [0.654,0.750][0.654,0.750] Through step 100 Loss spike count over steps 1 to 500 0.7580.758 [0.712,0.804][0.712,0.804] Through step 500 5.3 Effect of KSS on Stability KSS Results. Table 3 compares gradient clipping with KSS in the No-Norm setting. These two approaches differ fundamentally in their mechanism. Gradient clipping operates reactively: it caps gradients after explosion begins but does not prevent instability, yielding only modest improvements. KSS, in contrast, operates proactively by shaping the spectral distribution before instability occurs. With α=0.15α=0.15, KSS reduces divergence from 66.7% to 12.5%. Figure 2 visualizes the dose-response relationship between the KSS weight, M≈1M_≈ 1, and training stability. Table 3: Gradient clipping versus KSS in a challenging No-Norm setting with learning rate that ranges from 0.005 to 0.01 and 24 trials each. Measured on the associative-recall task. black Div.% (lower) Acc.% (higher) M≈1M_≈ 1 Overhead No control 66.766.7 28.528.5 0.850.85 — Gradient clip 0.5 50.050.0 32.132.1 0.820.82 <1%<1\% Gradient clip 1.0 58.358.3 30.830.8 0.830.83 <1%<1\% KSS α=0.10α=0.10 25.025.0 42.342.3 0.680.68 9.5%9.5\% KSS α=0.15α=0.15 12.512.5 48.248.2 0.580.58 10.8%10.8\% KSS α=0.20α=0.20 8.38.3 46.546.5 0.520.52 11.2%11.2\% Figure 2: KSS regularization effectiveness. (Left) The divergence rate decreases with KSS weight α. (Right) A dual axis shows accuracy improvement and M≈1M_≈ 1 shifting downward toward the target band. KSS shapes spectral properties, improving both stability and performance. Measured on the associative-recall task. Extended Baseline Comparison. Table 5.3 provides expanded baseline and optimizer results, including spectral normalization and weight normalization baselines [miyato2018spectralnorm, salimans2016weightnorm]. We use Adam with decoupled weight decay (AdamW) as the base optimizer in this comparison. Sharpness-Aware Minimization (SAM) [foret2020sam] reduces divergence to 33.3% but incurs a 2× computational cost, whereas KSS achieves a 2.7× lower divergence with 9× less overhead. The Lion optimizer [chen2023lion] yields 45.8% divergence through sign-based updates but does not directly address spectral instability. Combining KSS with SAM achieves the lowest divergence at 8.3% but with higher overhead. Table 4: Extended baseline comparison including SAM, spectral normalization, and the Lion optimizer. Measured on the associative-recall task. black Div.% (lower) Acc.% (higher) Overhead No stabilization, AdamW 66.766.7 28.528.5 — Gradient clipping at 1.0 58.358.3 30.830.8 <1%<1\% Spectral normalization 41.741.7 35.235.2 5.2%5.2\% Weight normalization 54.254.2 32.832.8 3.8%3.8\% Gradient penalty with λ=0.1λ=0.1 45.845.8 33.133.1 6.4%6.4\% Advanced Optimizers SAM, ρ=0.05ρ=0.05 37.537.5 38.638.6 about 100% SAM, ρ=0.10ρ=0.10 33.333.3 40.240.2 about 100% Lion optimizer 45.845.8 36.536.5 about 15% KSS, α=0.15α=0.15 12.512.5 48.248.2 10.8%10.8\% KSS + SAM, α=0.10α=0.10, ρ=0.05ρ=0.05 8.38.3 45.845.8 about 112% KSS Enables Higher Learning Rates. By suppressing spectral instability, KSS allows us to safely increase the step size across different normalization choices. Table 5 shows a 50% to 150% increase in the maximum stable learning rate under the same divergence criterion. Table 5: Maximum stable learning rate (LR). KSS enables learning rates that are 50% to 150% higher. The stability criterion is <<20% divergence across trials. Measured on the associative-recall task. black Type Max LR w/o KSS Max LR w/ KSS Increase Pre-LN 0.0050.005 0.0080.008 +60%+60\% RMSNorm 0.0030.003 0.0050.005 +67%+67\% No-Norm 0.0020.002 0.0050.005 +150%+150\% Mechanistic Evidence: KSS versus Random Regularization. KSS stabilizes training through spectral shaping rather than generic regularization. Table LABEL:tab:causal_ablation addresses this point via ablation studies with matched computational overhead. Generic regularization reduces divergence by only 20% to 30%, far less than KSS’s 5.3× reduction. Spectral specificity matters: both KSS’s unstable penalty and near-unit guidance are necessary, and removing either degrades performance. Let ΔM≈1≜M≈1KSS−M≈1base M_≈ 1 M_≈ 1^KSS-M_≈ 1^base; negative values indicate decreased near-unit mass. The correlation between ΔM≈1 M_≈ 1 and divergence reduction has r=−0.87r=-0.87 and p<0.001p<0.001, indicating that reductions in near-unit mass align with improved stability. (Supplementary Material) Appendix Table of Contents Appendix A List of Notation Table 12: List of notations used in the main paper and supplementary material. Symbol Meaning Core sizes and indices L The number of layers. ℓ Layer index. d Hidden dimension. N The number of snapshots. r Randomized DMD rank, the number of eigenvalues used in KSS. Dynamics and operators ℓh_ Residual stream at layer ℓ . Fℓ(⋅)F_ (·) Residual mapping Fℓ(ℓ)=ℓ+fℓ(ℓ;θℓ)F_ (h_ )=h_ +f_ (h_ ; _ ). K Koopman operator. ^ℓ A_ DMD estimate of the layer-ℓ Koopman operator. ,V, Eigenvectors and eigenvalues of ^ℓ A_ . λj _j Eigenvalue. ρ()ρ(A) Spectral radius. κ()κ(V) Eigenvector condition number, a measure of non-normality. ()K(A) Kreiss constant. Snapshots and whitening ℓ,ℓX_ ,Y_ Layer-ℓ snapshot matrices. ~ℓ,~ℓ X_ , Y_ Whitened snapshots. ^X _X Regularized sample covariance used for whitening. (⋅)†(·) Moore–Penrose pseudoinverse. Spectral diagnostics M>1M_>1 Unstable spectral mass. M≈1M_≈ 1 Near-unit spectral mass. M<1M_<1 Contractive spectral mass. ηnl _nl Nonlinearity ratio, a fit-error measure. ϵu,ϵn,δc _u, _n, _c Thresholds defining spectral-mass bins. KSS regularization ℒKSSℓL_KSS KSS loss for layer ℓ . α KSS regularization weight. τu,τl _u, _l Upper and lower band thresholds in KSS. γ Target near-unit mass level. mℓsoftm_ ^soft Soft near-unit mass estimate. Probability and norms D Divergence indicator. P(D=1∣)P(D=1 ) Predicted divergence probability. [⋅]E[·] Expectation. Tr(⋅)Tr(·) Trace. ∥⋅∥2,∥⋅∥F \|· \|_2, \|· \|_F Operator and Frobenius norms. |⋅| |· | Absolute value. ℝ,ℂR,C Real and complex number fields. Appendix B Additional Theoretical Results The following results extend the theoretical analysis presented in Section 4. B.1 Supplement to Theorem 1 Proposition 3 (Bauer–Fike: Non-normality Caveat). Let ∈ℂd×dA ^d× d be diagonalizable with eigendecomposition =−1A=V V^-1. For any perturbation E with ‖2≤δ \|E \|_2≤δ and any eigenvalue λ~∈spec⁡(+) λ (A+E), we have mink⁡|λ~−λk|≤κ()⋅δ. _k| λ- _k|≤κ(V)·δ. (20) Proof. Let λ~∈spec⁡(+) λ (A+E) with eigenvector ≠x 0, that is, (+)=λ~(A+E)x= λx. Write =−1A=V V^-1 and set ≜−1≠y ^-1x 0. Left-multiplying by −1V^-1 gives (−λ~)=−1.( - λI)y=-V^-1EVy. Taking Euclidean norms, ‖(−λ~)‖2≤‖−1‖2‖2‖2‖2=κ()δ‖2. \|( - λI)y \|_2≤ \|V^-1 \|_2 \|E \|_2 \|V \|_2 \|y \|_2=κ(V)δ \|y \|_2. Because −λ~ - λI is diagonal with diagonal entries (λj−λ~)( _j- λ), ‖(−λ~)‖2≥mink⁡|λk−λ~|⋅‖2. \|( - λI)y \|_2≥ _k| _k- λ|· \|y \|_2. Canceling ‖2 \|y \|_2 yields Eq. 20. ∎ B.2 DMD Convergence Analysis We first establish finite-sample convergence guarantees for DMD estimation in the presence of nonlinearity. Assumption 1 (Data Distribution). Let (,)∈ℝd×ℝd(x,y) ^d×R^d be a random pair drawn from the joint distribution induced by a layer transition. We assume centered covariates: []=E[x]=0 and [⊤]=E[xx ]= with σmin()≥σ0>0 _ ( )≥ _0>0. We also assume sub-Gaussian tails: ‖ψ2≤K \|x \|_ _2≤ K and ‖ψ2≤K \|y \|_ _2≤ K for some K>0K>0. Define the cross-covariance yx≜[⊤]C_yx [yx ]. Let ϵ≥0ε≥ 0 be the whitening regularizer in Eq. 7, and set ϵ≜+ϵ _ε + . Then σmin(ϵ)≥σ0 _ ( _ε)≥ _0. Define ϵ _ε ≜ϵ−1/2yxϵ−1/2, _ε^-1/2C_yx _ε^-1/2, (21) ϵ _ε ≜ϵ−1/2ϵ−1/2, _ε^-1/2 _ε^-1/2, (22) and the regularized whitened population least-squares operator w,ϵLS≜ϵϵ−1.A_w,ε^LS _εG_ε^-1. (23) The nonlinearity ratio ηnl _nl is a normalized linear-fit error that we use as a practical diagnostic for linear-approximation reliability. Theorem 4 (Whitened DMD Finite-Sample Convergence). Under Assumption 1, let σϵ≜σmin(ϵ) _ε _ ( _ε). Then σϵ≥σ0 _ε≥ _0. Let δfail∈(0,1) _fail∈(0,1). Consider the events ‖^X−ϵ‖2≤12σϵ, \| _X- _ε \|_2≤ 12 _ε, (24) and ‖^−ϵ‖2≤12σmin(ϵ), \| G-G_ε \|_2≤ 12 _ (G_ε), (25) where ^≜^X−1/2^^X−1/2 G _X^-1/2 _X^-1/2. Then G is invertible, and ‖^−1‖2≤2‖ϵ−1‖2 \| G^-1 \|_2≤ 2 \|G_ε^-1 \|_2. Assuming ~~⊤ X X is invertible, for example, when N≥dN≥ d and rank(~)=drank( X)=d, there exist absolute constants C,c>0C,c>0 such that if N≥c(d+log⁡(2/δfail))N≥ c(d+ (2/ _fail)), then on the event Eq. 24 and Eq. 25: ‖^−w,ϵLS‖2≤C‖ϵ−1‖2K2σϵΔN+C‖ϵ−1‖22K2‖yx‖2σϵ2ΔN, \| A-A_w,ε^LS \|_2≤ C \|G_ε^-1 \|_2 K^2 _ε _N+C \|G_ε^-1 \|_2^2 K^2 \|C_yx \|_2 _ε^2 _N, (26) where ΔN≜d+log⁡(2/δfail)N+d+log⁡(2/δfail)N _N d+ (2/ _fail)N+ d+ (2/ _fail)N. Proof. We bound the estimation error relative to the whitened population least-squares operator w,ϵLSA_w,ε^LS. Assume ~~⊤ X X is invertible, so that ~†=~⊤(~~⊤)−1 X = X ( X X )^-1 and ^=^^−1 A= M G^-1 with ^≜^X−1/2^yx^X−1/2 M _X^-1/2 C_yx _X^-1/2 and ^≜^X−1/2^^X−1/2 G _X^-1/2 _X^-1/2, where ¯=1N∑i=1Ni x= 1N _i=1^Nx_i, ¯=1N∑i=1Ni y= 1N _i=1^Ny_i, ^=1N−1∑i=1N(i−¯)(i−¯)⊤ = 1N-1 _i=1^N(x_i- x)(x_i- x) , ^yx=1N−1∑i=1N(i−¯)(i−¯)⊤ C_yx= 1N-1 _i=1^N(y_i- y)(x_i- x) , and ^X=^+ϵ _X= + . First, we estimate the covariance. For N independently and identically distributed samples with ‖ψ2≤K \|x \|_ _2≤ K, standard covariance concentration for the centered sample covariance yields [tropp2012tail] ‖^X−ϵ‖2≤CK2ΔN, \| _X- _ε \|_2≤ CK^2 _N, (27) with probability ≥1−δfail/2≥ 1- _fail/2 for N≳d+log⁡(1/δfail)N d+ (1/ _fail). Second, we bound the whitening perturbation. Standard perturbation theory for matrix square roots gives [higham2008functions]: ‖^X−1/2−ϵ−1/2‖2≤2σϵ3/2‖^X−ϵ‖2, \| _X^-1/2- _ε^-1/2 \|_2≤ 2 _ε^3/2 \| _X- _ε \|_2, (28) for ‖^X−ϵ‖2≤12σϵ \| _X- _ε \|_2≤ 12 _ε. Now, we combine the cross-covariance and covariance estimation errors. Let ^yx=1N−1∑i=1N(i−¯)(i−¯)⊤ C_yx= 1N-1 _i=1^N(y_i- y)(x_i- x) . A similar sub-exponential matrix concentration bound gives [tropp2012tail] ‖^yx−yx‖2≤CK2ΔN \| C_yx-C_yx \|_2≤ CK^2 _N with probability ≥1−δfail/2≥ 1- _fail/2. Decomposing ^−w,ϵLS=^^−1−ϵϵ−1=(^−ϵ)ϵ−1+^(^−1−ϵ−1), A-A_w,ε^LS= M G^-1-M_εG_ε^-1=( M-M_ε)G_ε^-1+ M ( G^-1-G_ε^-1 ), and using a standard matrix inverse perturbation bound, ‖^−1−ϵ−1‖2≤‖^−1‖2‖^−ϵ‖2‖ϵ−1‖2≤2‖ϵ−1‖22‖^−ϵ‖2 \| G^-1-G_ε^-1 \|_2≤ \| G^-1 \|_2 \| G-G_ε \|_2 \|G_ε^-1 \|_2≤ 2 \|G_ε^-1 \|_2^2 \| G-G_ε \|_2 on Eq. 25 yields a cross-covariance term scaling as ‖ϵ−1‖2σϵ−1 \|G_ε^-1 \|_2 _ε^-1 and a whitening term plus a covariance term scaling as ‖ϵ−1‖22‖yx‖2σϵ−2 \|G_ε^-1 \|_2^2 \|C_yx \|_2 _ε^-2, giving Eq. 26. Combining the bounds yields Eq. 26 under the event Eq. 24 and Eq. 25. ∎ Remark 1 (On ϵ−1G_ε^-1 for ϵ=+ϵ _ε= + ). Because and ϵ _ε commute, ϵG_ε has eigenvalues λi/(λi+ϵ) _i/( _i+ε), hence ‖ϵ−1‖2=σmin()+ϵσmin(). \|G_ε^-1 \|_2= _ ( )+ε _ ( ). If σmin() _ ( ) is treated as a fixed constant bounded away from 0, the factors of ‖ϵ−1‖2 \|G_ε^-1 \|_2 can be absorbed into the constant C. Remark 2 (Sample Complexity). Theorem 4 suggests that, under the stability events Eq. 24–Eq. 25, N=O~(d(‖ϵ−1‖2K2σϵ+‖ϵ−1‖22K2‖yx‖2σϵ2)2ε−2)N= O (d ( \|G_ε^-1 \|_2 K^2 _ε+ \|G_ε^-1 \|_2^2 K^2 \|C_yx \|_2 _ε^2 )^2 ^-2 ) samples are sufficient for ε -accurate DMD estimation, up to logarithmic factors. For typical transformers with d=256d=256 to 768768, N≈2048N≈ 2048 provides reliable estimates. Remark 3 (Modeling Mismatch and ηnl _nl). Theorem 4 is an estimation bound for the whitened population least-squares operator w,ϵLSA_w,ε^LS. When the layer transition is nonlinear, w,ϵLSA_w,ε^LS can be a poor proxy for other targets, such as a local Jacobian or a richer Koopman approximation, even if it is well-estimated. We use the empirical nonlinearity ratio ηnl _nl defined in Eq. 5 as a practical diagnostic for when linear DMD features are less reliable. B.3 Non-Normality and Transient Growth Spectral radius bounds alone are insufficient for analyzing non-normal matrices. The Kreiss matrix theorem provides tight bounds on transient behavior [kreiss1962stability, trefethen2005spectra]. Theorem 5 (Kreiss Constant Characterization). The Kreiss constant of ∈ℂd×dA ^d× d is: ()≜sup|z|>1(|z|−1)‖(z−)−1‖2K(A) _|z|>1(|z|-1) \|(zI-A)^-1 \|_2 (29) Assume A is power-bounded, that is, supn≥0‖n‖2<∞ _n≥ 0 \|A^n \|_2<∞; equivalently, ()<∞K(A)<∞. Then the Kreiss matrix theorem states: ()≤supn≥0‖n‖2≤e⋅d⋅()K(A)≤ _n≥ 0 \|A^n \|_2≤ e· d·K(A) (30) For a diagonalizable =−1A=V V^-1: ()≤κ()⋅sup|z|>1maxj⁡|z|−1|z−λj|K(A)≤κ(V)· _|z|>1 _j |z|-1|z- _j| (31) Proof. We first prove the two inequalities in Eq. 30 and then Eq. 31. Lower bound: ()≤supn≥0‖n‖2K(A)≤ _n≥ 0 \|A^n \|_2. Let M≜supn≥0‖n‖2<∞M _n≥ 0 \|A^n \|_2<∞. Because A is power-bounded, supn≥0‖n‖2<∞ _n≥ 0 \|A^n \|_2<∞, so the Neumann series converges in operator norm for any |z|>1|z|>1. For any |z|>1|z|>1, the Neumann series gives the operator-norm expansion (z−)−1=z−1∑n=0∞nz−n,(zI-A)^-1=z^-1 _n=0^∞A^nz^-n, hence ‖(z−)−1‖2≤1|z|∑n=0∞‖n‖2|z|n≤M|z|∑n=0∞|z|−n=M|z|−1. \|(zI-A)^-1 \|_2≤ 1|z| _n=0^∞ \|A^n \|_2|z|^n≤ M|z| _n=0^∞|z|^-n= M|z|-1. Multiplying by (|z|−1)(|z|-1) and taking the supremum over |z|>1|z|>1 yields ()≤MK(A)≤ M. Upper bound: supn≥0‖n‖2≤ed() _n≥ 0 \|A^n \|_2≤ edK(A). This is the finite-dimensional Kreiss matrix theorem: the resolvent bound sup|z|>1(|z|−1)‖(z−)−1‖2<∞ _|z|>1(|z|-1) \|(zI-A)^-1 \|_2<∞ is equivalent to power-boundedness, and quantitatively implies supn≥0‖n‖2≤Cd() _n≥ 0 \|A^n \|_2≤ C_dK(A) for an explicit dimension-dependent constant CdC_d; one standard choice is Cd=edC_d=ed. Diagonalizable case. If =−1A=V V^-1, then for any z∉spec⁡()z (A), (z−)−1=(z−)−1−1.(zI-A)^-1=V(zI- )^-1V^-1. Taking norms gives ‖(z−)−1‖2≤‖2‖(z−)−1‖2‖−1‖2=κ()‖(z−)−1‖2. \|(zI-A)^-1 \|_2≤ \|V \|_2 \|(zI- )^-1 \|_2 \|V^-1 \|_2=κ(V) \|(zI- )^-1 \|_2. Because (z−)−1(zI- )^-1 is diagonal with diagonal entries (z−λj)−1(z- _j)^-1, its spectral norm is maxj⁡|z−λj|−1 _j|z- _j|^-1, hence (|z|−1)‖(z−)−1‖2≤κ()⋅maxj⁡|z|−1|z−λj|.(|z|-1) \|(zI-A)^-1 \|_2≤κ(V)· _j |z|-1|z- _j|. Taking the supremum over |z|>1|z|>1 yields Eq. 31. ∎ The interpretation is as follows: a high ()K(A) indicates hidden instability. Even when ρ()≤1ρ(A)≤ 1, non-orthogonal eigenvectors produce a transient growth ‖n‖2≫1 \|A^n \|_2 1 for intermediate n. Appendix C PROOF OF THEOREM 1 Proof. For x uniform on the unit sphere, rotational invariance implies [∗]=1dE[xx^*]= 1dI. Therefore ‖22=[∗∗]=Tr(∗[∗])=1dTr(∗)=1d‖F2.E \|Ax \|_2^2=E[x^*A^*Ax]=Tr (A^*AE[xx^*] )= 1dTr(A^*A)= 1d \|A \|_F^2. For a normal A, ‖F2=∑j=1d|λj|2 \|A \|_F^2= _j=1^d| _j|^2, giving ‖22=1d∑j=1d|λj|2.E \|Ax \|_2^2= 1d _j=1^d| _j|^2. If ρ()≤1+ϵuρ(A)≤ 1+ _u, then |λj|2≤(1+ϵu)2| _j|^2≤(1+ _u)^2 for all j, giving the upper bound. For the lower bound, at least a fraction M≈1()M_≈ 1(A) of the eigenvalues satisfy |λj|≥1−ϵn| _j|≥ 1- _n, so ‖22≥(1−ϵn)2M≈1().E \|Ax \|_2^2≥(1- _n)^2M_≈ 1(A). For the diagonalizable extension =−1A=V V^-1, note that for an isotropic x we always have ‖22=1d‖F2E \|Ax \|_2^2= 1d \|A \|_F^2. Moreover, 1κ()‖F≤‖F≤κ()‖F 1κ(V) \| \|_F≤ \|A \|_F≤κ(V) \| \|_F, so ‖22E \|Ax \|_2^2 is within factors κ()±2κ(V)^± 2 of 1d∑j|λj|2 1d _j| _j|^2. Combining this with ρ()≤1+ϵuρ(A)≤ 1+ _u and the definition of M≈1()M_≈ 1(A) yields the stated bound. ∎ Proof of Corollary 2. Proof. Assume [ℓℓ∗]=1d‖ℓ‖22E[h_ h_ ^*]= 1dE \|h_ \|_2^2I for each ℓ=0,…,L−1 =0,…,L-1. Let qℓ≜1d∑j=1d|λjℓ|2q_ 1d _j=1^d| _j |^2. Then ‖ℓ+1‖22=Tr(ℓ∗ℓ[ℓℓ∗])=qℓ‖ℓ‖22,E \|h_ +1 \|_2^2=Tr (A_ ^*A_ E[h_ h_ ^*] )=q_ E \|h_ \|_2^2, and recursion yields ‖L‖22=‖0‖22∏ℓ=0L−1qℓE \|h_L \|_2^2=E \|h_0 \|_2^2 _ =0^L-1q_ . By Theorem 1, qℓ≥(1−ϵn)2M≈1(ℓ)q_ ≥(1- _n)^2M_≈ 1(A_ ), proving the stated bound. ∎ Appendix D Experimental Details We ran experiments on 4× NVIDIA A100-SXM4-40GB graphics processing unit (GPU) devices. We ran experiments across six normalization strategies, including Pre-LN, Post-LN, RMSNorm, DeepNorm, SubLN, and No-Norm, with additional architecture-specific studies on MoE [shazeer2017moe], Mamba [gu2023mamba], and KAN [liu2024kan]. We use the following hyperparameters. Models use d∈128,256,512,768,1024d∈\128,256,512,768,1024\, nheads∈4,8,16n_heads∈\4,8,16\, and L∈4,6,8,12,16,24L∈\4,6,8,12,16,24\. For training, the mini-batch size is 16 to 32 with 5 to 20 epochs and 100 to 1000 warmup steps. We use the AdamW and sweep learning rates within standard ranges. We use seeds 42,123,456\42,123,456\ for reproducibility. Our default recipe estimates randomized DMD eigenvalues with rank r=32r=32 using N=2048N=2048 snapshots with ϵ=10−5ε=10^-5. KSS is applied every 10 to 20 steps, sampling 50% of layers per update to reduce overhead. We sweep the regularization weight α∈0.01,0.05,0.10,0.15,0.20α∈\0.01,0.05,0.10,0.15,0.20\; the resulting overhead ranges from 8% to 12% in practice. D.1 Associative-Recall Task The associative-recall task refers to a synthetic key-value retrieval classification task commonly used to probe associative memory and recall in long-context sequence models [fu2023hungry, arora2024zoology]. For each sample, we generate npairsn_pairs key-value pairs (ki,vi)i=1npairs\(k_i,v_i)\_i=1^n_pairs and construct =[k1,v1,…,knpairs,vnpairs,p1,…,pm,q],x=[k_1,v_1,…,k_n_pairs,v_n_pairs,p_1,…,p_m,q], where q is a query key chosen from ki\k_i\ and the label is the corresponding matched value y=vjy=v_j. The model predicts only this final target with cross-entropy on the last-position logits, implemented as F.cross_entropy(logits[:, -1, :], y). D.2 Synthetic LM Task This synthetic token-level language modeling setup is inspired by prior work using controlled synthetic sequences to analyze recall and long-range behavior in efficient sequence models [fu2023hungry, arora2024zoology]. Each sample is a length-T token sequence over a vocabulary of size V; unless stated otherwise, we use T=256T=256 and V=10,000V=10,000. A sequence is generated by concatenating randomly sampled segments until reaching length T, then truncating. Segments come in three types: repetition segments repeat a short pattern of length 2 to 5 for 2 to 4 repeats, sequential segments are contiguous integer runs of length 5 to 15, and random segments are independently and identically distributed tokens of length 5 to 15. All tokens are sampled from 10,…,V−1\10,…,V-1\ so that a small identifier range remains available for special tokens. Sequences are generated directly at the token level by these rules. The learning objective is standard next-token prediction. Given tokens (x1,…,xT)(x_1,…,x_T), the model predicts xt+1x_t+1 from the prefix (x1,…,xt)(x_1,…,x_t) and is trained with token-level cross-entropy over t=1,…,T−1t=1,…,T-1. We report validation token accuracy and perplexity exp⁡(mean cross-entropy) (mean cross-entropy) on a held-out synthetic validation split. For large-scale runs, we use 20K training sequences and 2K validation sequences per trial, regenerated deterministically from the run seed. D.3 Pretrained LM Fixed-Prompt Protocol We use a forward-only profiling protocol with a deterministic text set and no fine-tuning updates, in the same spirit as prompt-based evaluation and activation-probing analyses of pretrained transformers [brown2020language, elhage2021mathematical]. The fixed short-prompt setting is implemented as explicit prompt lists with 32 total prompts per run, either 4 prompts repeated 8 times or 8 prompts repeated 4 times, and both variants use a token length cap of 64. For each model, we run a single batched forward pass with hidden-state outputs enabled and collect residual-stream activations from the embedding and transformer layer outputs. For each layer transition (ℓ,ℓ+1)( , +1), we flatten token positions, subsample up to N∈1024,2048N∈\1024,2048\ token states, and apply whitened DMD. Spectral partitions use the same thresholds as the analysis code: unstable when |λ|>1.05|λ|>1.05, near-unit when 0.90≤|λ|≤1.050.90≤|λ|≤ 1.05, and over-damped when |λ|<0.80|λ|<0.80. We report early, middle, and late summaries by splitting layers into depth thirds and averaging each metric within the corresponding group. Appendix E Computational Cost Whitened DMD requires O(d2N+d3)O(d^2N+d^3) operations. Randomized singular value decomposition reduces this cost to O(dNr+r3)O(dNr+r^3) for rank-r approximation. With typical values d=768d=768, N=2048N=2048, and r=32r=32, full RKSP analysis completes in 2.5 to 3.5 seconds per layer on a single GPU. Appendix F Extended Baseline and Optimizer Comparisons F.1 Extended Optimizer Baselines We extend the baseline comparisons to include μ and the Layer-wise Adaptive Moments for Batch training (LAMB) optimizer [you2019lamb]. μ . The μ enables hyperparameter transfer across different model widths by appropriately scaling learning rates appropriately. Table 13 summarizes the μ comparisons and the combined μ + KSS setting. We observe moderate stability gains from μ , with a 38% divergence reduction through initialization scaling. The combined μ + KSS approach achieves the lowest divergence, 8.3%, and the highest accuracy, 51.5%. The μ setting enables transfer, while KSS provides stability, and the distinct mechanisms suggest orthogonal benefits. Table 13: Comparison of μ , standard parameterization, and KSS. Results use the No-Norm setting with learning rate ranges from 0.005 to 0.01 across 24 trials. Measured on the associative-recall task. black Div.% (lower) Acc.% (higher) M≈1M_≈ 1 Transfer Overhead Standard Init, AdamW 66.766.7 28.528.5 0.850.85 ✗ — μ , width transfer 41.741.7 35.835.8 0.720.72 ✓ <1%<1\% μ + higher learning rate 50.050.0 38.238.2 0.780.78 ✓ <1%<1\% KSS, α=0.15α=0.15 12.512.5 48.248.2 0.580.58 ✗ 10.8%10.8\% μ + KSS 8.38.3 51.551.5 0.520.52 ✓ 11.2%11.2\% LAMB Optimizer. LAMB normalizes updates per layer, which changes spectral dynamics. Table 14 reports the optimizer comparison, including LAMB and Lion. LAMB outperforms AdamW: its layer-wise normalization provides implicit stability with a 31% divergence reduction. LAMB enables higher learning rates, reaching 3× to 4× AdamW, or up to 7× with learning rate warmup. KSS still provides orthogonal benefits, and the LAMB + KSS combination achieves the best results at 8.3% divergence and 52.8% accuracy. From a spectral perspective, LAMB suppresses expansive mass and reduces near-unit mass, from 0.85 with AdamW to 0.75, suggesting that its stability gains come from damping unstable modes and increasing damping; within comparable unstable-mass regimes, a lower M≈1M_≈ 1 aligns with greater stability. Table 14: Optimizer comparison among AdamW, LAMB, and Lion. Results use the No-Norm setting, 24 trials each. Measured on the associative-recall task. black Div.% (lower) Acc.% (higher) M≈1M_≈ 1 Best LR Overhead AdamW 66.766.7 28.528.5 0.850.85 0.0030.003 — AdamW + grad clip 58.358.3 30.830.8 0.830.83 0.0050.005 <1%<1\% LAMB 45.845.8 36.236.2 0.750.75 0.010.01 about 5% LAMB + warmup 37.537.5 40.840.8 0.680.68 0.020.02 about 5% Lion 45.845.8 36.536.5 0.780.78 0.0010.001 about 15% KSS, α=0.15α=0.15 12.512.5 48.248.2 0.580.58 0.0080.008 10.8%10.8\% LAMB + KSS 8.38.3 52.852.8 0.520.52 0.0150.015 about 16% Appendix G Large-Scale Pretrained Model Analysis GPT-2 Analysis. Table 15 reports layer-group spectral statistics for GPT-2 [radford2019language], and Figure 4 reveals a universal pattern. The normalized linear-fit error ηnl _nl increases with depth, rising from [0.48,0.52][0.48,0.52] in the early layers to [0.68,0.71][0.68,0.71] in late layers. Simultaneously, the near-unit mass decreases from M≈1∈[0.68,0.72]M_≈ 1∈[0.68,0.72] to M≈1∈[0.58,0.60]M_≈ 1∈[0.58,0.60]. This depth-wise trend has a clear implication: early layers are more linearly approximable, making DMD features more reliable, whereas late layers are less so. Table 15: GPT-2 layer-wise spectral analysis. Start Linear, End Nonlinear pattern, shorthand for increasing ηnl _nl. Spectral statistics are computed from residual-stream activations on a fixed set of short prompt sentences. Measured on a fixed short-prompt set. black Layers ρ κ()κ(V) ηnl _nl M≈1M_≈ 1 GPT-2 124M Early 0 to 3 1.151.15 8.38.3 0.520.52 0.680.68 GPT-2 124M Middle 4 to 7 1.231.23 10.110.1 0.610.61 0.620.62 GPT-2 124M Late 8 to 11 1.311.31 11.211.2 0.710.71 0.580.58 GPT-2 355M Early 0 to 7 1.121.12 7.57.5 0.480.48 0.720.72 GPT-2 355M Middle 8 to 15 1.191.19 9.29.2 0.580.58 0.650.65 GPT-2 355M Late 16 to 23 1.271.27 10.810.8 0.680.68 0.600.60 Figure 4: Start Linear, End Nonlinear pattern. Layer-wise normalized linear-fit error ηnl _nl across four pretrained models. All models exhibit a monotonically increasing ηnl _nl with depth, suggesting a consistent linear-approximation signature across models. Computed from residual-stream activations on a fixed set of short prompt sentences. Appendix H Calibration Beyond discrimination measured by AUROC, we assess calibration quality. Figure 5 shows that the risk score M≈1M_≈ 1 achieves an Expected Calibration Error (ECE) of 0.283, indicating moderate calibration. The reliability diagram reveals deviations between predicted probabilities and observed frequencies, while the distribution plots show clear separation between converged runs with lower M≈1M_≈ 1, corresponding to lower risk, and diverged runs with higher M≈1M_≈ 1, corresponding to higher risk. This calibration quality matters for deployment: practitioners can interpret RKSP’s probability estimates for early termination decisions while accounting for the moderate calibration. Figure 5: Calibration reliability diagram. (Left) Predicted divergence probability versus observed frequency, with an ECE of 0.283. (Right) Distribution of predictions separated by actual outcome. RKSP provides moderately calibrated probability estimates. Based on associative-recall runs, calibration compares predictions to divergence outcomes from that task. Appendix I Novel Architecture Case Studies To demonstrate RKSP’s value beyond standard transformers, we analyze three emerging architectures: MoE [shazeer2017moe], SSMs including Mamba [gu2023mamba], and KAN [liu2024kan]. MoE Transformers Table 16 presents a comparison of MoE routing and stability [shazeer2017moe]. MoE routing induces a higher normalized linear-fit error: ηnl _nl increases 15% to 20% compared to dense transformers due to discrete routing decisions. The choice of top-k affects spectral stability—higher k shifts spectral mass and changes M≈1M_≈ 1 alongside non-normality and unstable modes. The divergence reductions are consistent with suppressing unstable modes, and within comparable non-normality regimes, larger M≈1M_≈ 1 aligns with greater instability. KSS stabilizes MoE effectively, yielding a 4× divergence reduction with 6% accuracy improvement, thereby validating RKSP and KSS for novel architectures. Table 16: MoE transformer with RKSP analysis. Routing instability revealed via spectral signatures. Results use d=256d=256, L=6L=6, and 24 trials. Load balancing loss λ=0.01λ=0.01. Measured on the synthetic LM task with random-token next-token prediction. black Top-k M≈1M_≈ 1 ρ ηnl _nl Div. of 24 Acc.% MoE-Small, 8 experts k=1k=1 0.720.72 2.852.85 0.680.68 6 of 24 32.432.4 MoE-Small, 8 experts k=2k=2 0.580.58 2.122.12 0.550.55 2 of 24 41.741.7 MoE-Small, 8 experts k=4k=4 0.450.45 1.781.78 0.480.48 1 of 24 38.238.2 MoE-Medium, 16 experts k=2k=2 0.650.65 2.452.45 0.610.61 4 of 24 38.938.9 MoE-Medium, 16 experts k=2k=2 + KSS 0.480.48 1.921.92 0.580.58 1 of 24 44.544.5 State Space Models: Mamba Table 17 compares SSM and transformer spectral properties for Mamba [gu2023mamba]. The theoretical explanation is straightforward: SSMs are designed with stable discrete-time dynamics via highly structured polynomial projection operator initialization. RKSP reveals this design choice explicitly in the spectral signature: Mamba exhibited M<1≈0.85≫M≈1≈0.12M_<1≈ 0.85 M_≈ 1≈ 0.12. This separation indicates strongly contractive dynamics with short memory and weak near-isometry; stability here comes from suppressed unstable modes in a highly contractive regime. In transformer regimes that are closer to near-normal, larger M≈1M_≈ 1 corresponds to weaker damping, longer-range signal retention, and higher instability risk. Table 17: Mamba with RKSP analysis. Inherently stable spectral structure. Results use L=6L=6, 24 trials. Measured on the synthetic LM task with random-token next-token prediction. black M≈1M_≈ 1 M<1M_<1 ρ Div. of 24 Acc.% Mamba-Small, d=256d=256 0.120.12 0.850.85 0.950.95 0 of 24 48.248.2 Mamba-Medium, d=512d=512 0.150.15 0.820.82 0.970.97 0 of 24 52.652.6 Transformer, Pre-LN, comparable 0.420.42 0.380.38 1.851.85 1 of 24 45.845.8 KAN Table 18 reports KAN spectral diagnostics and KSS outcomes [liu2024kan]. KAN shows high normalized linear-fit error: B-spline basis functions produce ηnl≈[0.78,0.82] _nl≈[0.78,0.82], higher than the typical transformer layers with [0.4,0.7][0.4,0.7]. Despite this high ηnl _nl, RKSP remains informative—ResDMD filtering enables spectral analysis for 68% to 78% of modes. KSS benefits KAN with a 3× to 4× divergence reduction, suggesting that spectral shaping is architecture-agnostic. Table 18: KAN transformer with RKSP analysis. The B-spline nonlinearity challenges linear approximation. The B-spline order is B. We use d=256d=256, L=6L=6, and 24 trials. Measured on the synthetic LM task with random-token next-token prediction. black ηnl _nl M≈1M_≈ 1 ρ Div. of 24 DMD reliability KAN-Transformer, B=4B=4 0.780.78 0.520.52 2.152.15 3 of 24 Marginal, 68% KAN-Transformer, B=8B=8 0.820.82 0.480.48 2.352.35 4 of 24 Low, 52% KAN-Transformer + KSS 0.750.75 0.420.42 1.921.92 1 of 24 Improved, 78% Cross-Architecture Summary Table LABEL:tab:arch_compare summarizes cross-architecture metrics, while Figure 6 provides a normalized radar-chart view of the same comparison. Figure 6: Cross-architecture spectral radar chart. Comparison of five architectures across five normalized metrics. Mamba exhibits strong contraction with low M≈1M_≈ 1 and short memory; stability is maintained via suppressed unstable modes, while in near-normal transformer regimes, higher M≈1M_≈ 1 aligns with more unstable, near-isometric propagation. The No-Norm transformer shows high memory capacity but poor stability. KAN exhibits high ηnl _nl. Metrics are derived from Table LABEL:tab:arch_compare. Appendix J Practical Notes RKSP and KSS are most valuable in three scenarios. First, when mechanistic understanding matters, RKSP explains why Pre-LN outperforms Post-LN through spectral signatures. Second, when pushing training limits, KSS enables learning rates that are 50% to 150% higher for faster convergence. Third, when deploying novel architectures, RKSP verifies stability before expensive training runs. Edge cases benefit most from these diagnostics—situations where standard normalization fails or where training operates near stability boundaries. Fixup and ReZero-style identity initialization. A common stabilization trick in deep residual networks and transformers is to initialize the final projection of each residual branch to zero, for example the attention and MLP output weights, so that the network starts close to an identity map [zhang2019fixup, bachlechner2020rezero]. In our notation, this yields a residual-off regime with a vanishing layer update ℓ+1−ℓ≈h_ +1-h_ 0, so the snapshot pairs satisfy ℓ≈ℓY_ _ and DMD returns ^ℓ≈ A_ . Consequently, M≈1ℓM_≈ 1 can be close to 11 across layers even though training is often stable under Fixup and ReZero at initialization. Taken alone, a near-identity spectrum might seem to imply maximal instability risk. However, our instability mechanism assumes two conditions: weak damping with large M≈1M_≈ 1 under near-normality, and non-degenerate layer-wise dynamics with appreciable updates so that perturbations and optimization noise are repeatedly injected and propagated across depth. Fixup and ReZero violate the second condition at initialization. When ‖~ℓ−~ℓ‖F \| Y_ - X_ \|_F is near zero, there is essentially no layer-wise update to analyze, and the resulting DMD spectrum is not informative about the noisy training-time regime we target. Practically, this degeneracy is detectable from the same quantities RKSP already computes. When ‖~ℓ−~ℓ‖F≈0 \| Y_ - X_ \|_F≈ 0, the normalization in the nonlinearity ratio Eq. 5 becomes ill-conditioned, so ηnl(ℓ) _nl( ) should be interpreted as a DMD reliability flag rather than as a meaningful nonlinearity estimate. For Fixup and ReZero, RKSP becomes informative after a small amount of training, once the zero-initialized residual projections move away from zero and layer-wise updates become observable; at that point, RKSP can again capture whether the residual stream exhibits excessive near-isometric propagation (large M≈1M_≈ 1) that correlates with high-learning-rate divergence. Practical deployment is straightforward. We recommend using RKSP in four scenarios: first, as a fast filter during architecture search; second, before expensive hyperparameter grid search; third, for periodic spectral monitoring during training; and fourth, for debugging checkpoints before divergence. Figure 7 provides an actionable decision process. Run RKSP at initialization. Compute M≈1M_≈ 1 and κ()κ(V) for the decision. Meets safe-region criteria M≈1<0.3M_≈ 1<0.3 and modest κ()κ(V) Proceed without KSS. Under this criterion, training is more stable; optionally, monitor with periodic RKSP snapshots. Decide on KSS by context. If M≈1>0.5M_≈ 1>0.5 or κ()κ(V) is large, treat as high risk. Use KSS for aggressive learning rates above 2×2× standard, for no-normalization settings, or for novel architectures with RKSP monitoring. Skip KSS for standard Pre-LN and RMSNorm at conservative learning rates. yesno Figure 7: Decision flowchart for when to use RKSP and KSS in practice.

Instant research discovery