Approximate Linear Solvers for Scalable Statistical Algorithms

Abstract

Many statistical methods for large, multivariate datasets require solving linear systems whose dimensions are determined by either the number of observations or the size of individual data vectors. This often becomes a bottleneck when scaling these methods to accommodate large and complex datasets. In particular, the direct solution of such systems can become computationally expensive and memory-intensive, making it impractical for real-world applications involving massive data sizes. In this thesis, we explore the application of Krylov subspace methods to address this computational challenge, focusing on a source separation problem in cosmology where the data size is prohibitively large for conventional methods. Krylov subspace methods offer a way to iteratively approximate the solution of linear systems, reducing the computational complexity while maintaining accuracy. Two distinct approaches adapted from existing literature are presented: one applies the conjugate gradient method directly to the Kronecker-structured problem, and the other reformulates the problem as a Sylvester matrix equation. Both approaches are analyzed theoretically and their robust performance is demonstrated on large-scale simulated and real-world datasets.

The key contribution of this thesis is the development and implementation of approximate linear solvers that enables the scalability of statistical algorithms to very large datasets. By integrating Krylov subspace techniques with existing statistical models, we demonstrate a significant reduction in both memory requirements and computation time without compromising the accuracy of the solution. This contribution is particularly valuable for the cosmology domain, where the analysis of large-scale observational data demands highly efficient methods for solving linear systems. Moreover, this thesis highlights the broader applicability of these approximate solvers in other fields of study, offering a scalable solution to the linear systems problem that arises in many high-dimensional statistical algorithms. The techniques developed here can be easily adapted to other statistical applications involving large datasets, such as machine learning and signal processing, where scalability is a critical factor.