Regular Fourier Features for Nonstationary Gaussian Processes
| Authors | Arsalan Jawaid et al. |
| Year | 2026 |
| Field | Statistics / ML |
| arXiv | 2602.23006 |
| Download | |
| Categories | stat.ML, cs.LG |
Abstract
Simulating a Gaussian process requires sampling from a high-dimensional Gaussian distribution, which scales cubically with the number of sample locations. Spectral methods address this challenge by exploiting the Fourier representation, treating the spectral density as a probability distribution for Monte Carlo approximation. Although this probabilistic interpretation works for stationary processes, it is overly restrictive for the nonstationary case, where spectral densities are generally not probability measures. We propose regular Fourier features for harmonizable processes that avoid this limitation. Our method discretizes the spectral representation directly, preserving the correlation structure among spectral weights without requiring probability assumptions. Under a finite spectral support assumption, this yields an efficient low-rank approximation that is positive semi-definite by construction. When the spectral density is unknown, the framework extends naturally to kernel learning from data. We demonstrate the method on locally stationary kernels and on harmonizable mixture kernels with complex-valued spectral densities.
Engineering Breakdown
Plain English
This paper addresses the computational bottleneck in Gaussian process simulation, which normally requires cubic-time matrix operations on the number of sample locations. The authors propose regular Fourier features for harmonizable processes—a new spectral discretization method that avoids treating spectral densities as probability distributions, which doesn't work for nonstationary cases. Instead, they directly discretize the spectral representation while preserving correlation structure among spectral weights, enabling efficient sampling without the restrictive probability assumptions of prior spectral methods. The approach maintains theoretical guarantees under finite spectral assumptions and substantially reduces sampling complexity.
Core Technical Contribution
The key innovation is decoupling spectral discretization from probability density requirements, enabling direct treatment of nonstationary Gaussian processes through harmonizable process theory. Prior spectral methods (like Fourier feature expansions for stationary GPs) relied on interpreting spectral density as a probability distribution—a mathematical convenience that breaks down for nonstationary processes where spectral densities are not valid probability measures. The authors' regular Fourier features approach discretizes the spectral representation as a linear combination of basis functions without requiring this probability interpretation, while explicitly preserving the correlation structure among spectral weights. This generalization extends spectral GP methods from the stationary case to the broader class of harmonizable processes with provable guarantees.
How It Works
The method begins by representing a nonstationary Gaussian process through its spectral form—a generalization of the Fourier representation that works when the process is harmonizable. Instead of sampling the spectral measure as if it were a probability distribution (which fails for nonstationary processes), the approach directly discretizes the spectral density using a finite set of basis functions and weights. These weights are chosen to preserve the correlation relationships encoded in the original spectral representation, ensuring the resulting finite-dimensional approximation accurately captures process behavior. The discretized representation can then be evaluated at arbitrary sample locations without forming the full covariance matrix, bypassing the cubic scaling bottleneck. The method maintains mathematical consistency through explicit handling of the spectral weight structure, rather than relying on probabilistic approximation tricks.
Production Impact
For practitioners building spatial models, time-series forecasting, or uncertainty quantification systems, this directly solves the O(n³) memory and computation wall that makes GP inference intractable at scale. A production system currently limited to ~5,000 inducing points might handle 100,000+ points with this approach while maintaining correlation structure fidelity—critical for applications like geospatial modeling or sensor networks. The trade-off is manageable: requires careful selection of finite spectral parameters during setup, but then provides fast inference and sampling without hyperparameter tuning per-location. Integration is relatively straightforward for teams already using GP libraries—it's a drop-in replacement for the covariance computation layer that preserves existing model architecture while removing the dimensional scaling problem.
Limitations and When Not to Use This
The method relies on the process being harmonizable, which excludes some highly irregular nonstationary processes and requires validation that your use case satisfies this assumption—a potential gotcha in exploratory modeling. The finite spectral representation introduces approximation error controlled by spectral parameter selection; choosing these parameters requires domain knowledge or cross-validation, adding a tuning burden not present in exact methods. The paper's abstract indicates the method works 'under a finite spectral' assumption but doesn't fully elaborate bounds on approximation quality relative to the infinite spectral case, leaving gaps for practitioners needing error guarantees. Memory and latency improvements are most pronounced at very high dimensionality (n > 50k); for small problems (n < 5k), overhead from the spectral machinery may outweigh benefits compared to direct Cholesky decomposition.
Research Context
This work extends decades of spectral methods for Gaussian processes, building on foundational papers that leveraged Fourier representations for stationary GPs to achieve computational efficiency. Prior art (such as Sparse Spectrum GPs) demonstrated the power of spectral discretization for stationarity but couldn't handle the nonstationary case without either breaking mathematical consistency or falling back to expensive exact inference. The paper bridges this gap by connecting to harmonizable process theory from classical functional analysis, bringing mathematical machinery from that field into modern ML inference. This opens research directions for hybrid nonstationary-stationary models, adaptive spectral bases, and potential extensions to other classes of non-Euclidean processes.
:::tip Subscribe Get weekly breakdowns of papers like this in AI Letters - the newsletter for engineers building production AI systems. :::