numerical linear algebra

This work addresses the gap between algorithmic prototypes and efficient implementations in scientific research by proposing a lightweight approach to translate statistical and machine learning algorithms—such as kernel ridge regression and stochastic gradient descent matrix factorization—from mathematical formulations into readable, high-performance C++ code. Leveraging the Eigen template library for core linear algebra operations—including kernel matrix construction, regularized solvers, and vectorized updates—the implementation seamlessly integrates into the Python ecosystem via pybind11, enabling efficient interoperability with NumPy arrays. The project provides concise, reproducible code examples that encapsulate common computational patterns in research, significantly lowering the barrier for researchers to adopt C++ for high-performance development while balancing performance, readability, and usability.

This work addresses the challenge of efficiently computing leading eigenvectors in dynamic graphs, where frequent updates to the adjacency or Laplacian matrix render traditional eigendecomposition methods computationally prohibitive. To overcome this limitation, the authors propose a fast spectral embedding update framework based on Rayleigh-Ritz projection. By leveraging eigenvector perturbation analysis, the method constructs a low-dimensional approximate invariant subspace that preserves high approximation accuracy while substantially reducing computational and memory costs. Experimental results demonstrate that the proposed approach outperforms existing techniques in both the quality of leading eigenvector approximation and performance on downstream tasks—such as influential node identification and node clustering—offering a compelling balance between efficiency and accuracy.

This paper addresses the weak theoretical foundations of matrix decomposition in machine learning by systematically constructing a self-consistent, comprehensive, and modern-application-oriented pedagogical framework. Methodologically, it grounds the exposition in numerical linear algebra and matrix analysis, unifying classical decompositions—including LU, QR, SVD, and block triangular factorizations—while integrating numerical stability analysis and Hermitian/Hilbert space theory. Crucially, it bridges traditional numerical analysis with deep learning’s backpropagation setting, emphasizing differentiability and computational robustness of decompositions in algorithm design and model optimization. The primary contribution is a compact, dual-purpose (teaching and research) knowledge system that fills critical gaps in both theoretical coherence and machine-learning relevance present in existing literature, thereby providing rigorous mathematical foundations for high-dimensional data modeling and efficient training.

This paper investigates the ability of Invariant Coordinate Selection (ICS) to recover the Fisher Discriminant Subspace (FDS) in multi-class settings (≥3 classes). Addressing both rank-deficient and full-rank group center matrices, we systematically extend the ICS theoretical framework—marking the first rigorous generalization of ICS guarantees from binary to multi-class clustering. Using joint diagonalization of dual scatter matrices, multivariate statistical analysis, and numerical simulations, we prove that ICS exactly recovers the FDS under broad combinations of generalized scatter matrices; recovery fails only under specific, rare degeneracy conditions. Empirical evaluations confirm high robustness and practical efficacy. Our work establishes a solid theoretical foundation for high-dimensional multi-class discriminant analysis and provides a principled, computationally tractable tool for subspace identification.

Parameter estimation in high-dimensional structured generalized linear models suffers from low efficiency, particularly under realistic design matrices exhibiting anisotropy and strong correlations. Method: This paper introduces a novel spectral estimation framework based on Approximate Message Passing (AMP). Contribution/Results: We provide the first exact asymptotic characterization of spectral estimators under correlated Gaussian designs. We identify a universally optimal covariance-adaptive preprocessing strategy, partially resolving a long-standing conjecture on optimal spectral estimation for rotationally invariant models. Theoretically and empirically, our approach substantially reduces sample complexity and achieves provably statistically optimal estimation accuracy—outperforming existing heuristic methods on canonical designs from computational imaging and genomics.

This work addresses the tendency of Laplacian-constrained graphical models—such as the Laplacian-constrained Gaussian graphical model (LCGGM) and the Hüsler–Reiss model—to yield overly dense graphs during structure learning, which compromises interpretability and scalability. The paper introduces, for the first time, spectral graph sparsification as a post-processing step: without requiring additional hyperparameter tuning, it replaces the original Laplacian estimate with a spectrally approximated sparse Laplacian and refits the model. This approach effectively enhances both the sparsity and accuracy of the estimated graph structure. Integrating spectral graph theory, Laplacian-constrained Gaussian graphical models, extreme-value graphical models, and graph sparsification techniques, the method demonstrates superior performance on Erdős–Rényi and stochastic block model simulations and validates its practical utility on real-world data.

This study addresses a fundamental challenge in system identification: distinguishing spurious eigenvalues arising from limited data from those genuinely reflecting the underlying system dynamics. To this end, the paper introduces—for the first time—the probabilistic sampling pseudospectrum \( P(\lambda) \) and its computationally efficient estimator \( \hat{P}(\lambda) \). By leveraging resampling and statistical inference, this framework quantifies the uncertainty of eigenvalues across the complex plane. The proposed approach provides a general and rigorous statistical criterion for data-driven methods such as Dynamic Mode Decomposition and subspace identification, substantially enhancing the reliability of identifying true dynamical modes from noisy, finite-length observations.

This work addresses the challenge that real-world data in generalized linear models often violate the independent and identically distributed (i.i.d.) assumption. Under the relaxed assumption that the design matrix is orthogonally invariant—meaning its singular vectors are uniformly distributed while singular values remain arbitrary—the paper proposes an efficient parameter estimation method combining optimal spectral initialization with Approximate Message Passing (AMP). The proposed approach achieves the information-theoretically optimal sample complexity for weak recovery and attains the fundamental lower bound on estimation error, thereby extending beyond the classical i.i.d. Gaussian design setting. Rigorous theoretical analysis provides strong performance guarantees, and numerical experiments confirm both the algorithm’s effectiveness and the accuracy of the theoretical predictions on orthogonally invariant as well as more general correlated data.

This work addresses the lack of dedicated Python tools for Invariant Coordinate Selection (ICS) in the current machine learning ecosystem, which has hindered its adoption in tasks such as clustering and anomaly detection. We propose a dimensionality reduction method based on the joint diagonalization of two scatter matrices, effectively extracting features by optimizing directions corresponding to extremal generalized kurtosis. To enhance accessibility and reproducibility, we introduce ICSpyLab—the first open-source implementation of ICS under the MIT license—featuring support for multiple scatter matrix types, joint diagonalization algorithms, and component selection criteria. The library strictly adheres to the scikit-learn estimator interface and is accompanied by comprehensive documentation and reproducible examples, significantly improving the usability and extensibility of ICS methods.

This work proposes a unified, imaging-modality-agnostic geometric framework to characterize the information structure of imaging operators. By mapping the normalized singular spectrum onto a probability simplex and introducing the Fisher–Rao information metric, the authors construct a Riemannian geometry of constant positive curvature. This approach yields, for the first time, an intrinsic description at the operator level that does not rely on optimization, stochastic modeling, or modality-specific assumptions, thereby revealing the geometric invariance of spectral equivalence classes. The study derives closed-form expressions for information distance and geodesics, proves their invariance under unitary transformations and global scaling, and elucidates the mechanism by which nonlinear redistribution of spectral weights influences information distance.