🤖 machine learning

Lightweight Gaussian Process Inference in C++ on Metal and CUDA

This paper introduces LightGP, a lightweight, dependency-free C++17 library for Gaussian process regression that leverages optimized CPU, Metal, and CUDA backends to achieve significant speedups over existing Python-based frameworks like GPyTorch across a wide range of dataset sizes.

Original authors: Yu-Hsueh Fang

Published 2026-05-19

📖 5 min read🧠 Deep dive

CC BY 4.0

Original authors: Yu-Hsueh Fang

Original paper licensed under CC BY 4.0 (http://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Imagine you are trying to predict the weather for the next week. You have a massive amount of historical data, and you want to use a "smart calculator" called a Gaussian Process (GP) to make that prediction. This calculator is famous for being very accurate and for telling you how confident it is in its answer.

However, there's a catch: this calculator is incredibly heavy. Traditionally, to run it, you need to load a giant, bloated software suite (like a full kitchen with every appliance imaginable) just to use one blender. This is what current popular tools like GPyTorch do; they are built on top of huge deep-learning frameworks that take up gigabytes of space and slow things down with extra steps.

Enter LightGP.

The authors of this paper built a new, ultra-lightweight version of this calculator. Think of it as a Swiss Army knife instead of a full kitchen. It's written in a lean, efficient language (C++) and doesn't require any giant software suites to run. It can fit on a phone, a laptop, or a server without needing to download 2GB of extra baggage.

Here is how LightGP works, explained through everyday analogies:

1. The Four Different "Paths"

The paper says LightGP offers four different ways to solve the problem, depending on how big your data is:

The Exact Path (Cholesky): Like solving a math problem by writing out every single step perfectly. It's super accurate but gets very slow and memory-heavy if you have too many numbers (like trying to solve a puzzle with a million pieces).
The "Guess-and-Check" Path (Conjugate Gradients): Instead of writing out the whole puzzle, this method takes a smart shortcut. It asks, "If I try this, am I getting closer?" and adjusts. It's like navigating a maze by feeling the walls rather than drawing the whole map first. This saves a huge amount of memory.
The "Sample" Path (Sparse Variational): Imagine you have a million photos of cats, but you only need to look at 200 of them to learn what a cat looks like. This method picks a few "representative" data points to do the heavy lifting, ignoring the rest.
The "Fast-Forward" Path (SKI with FFT): This uses a special mathematical trick (like a magic lens) to blur the details just enough to see the big picture instantly, making it incredibly fast for huge datasets.

2. The Race: LightGP vs. The Giants

The authors put LightGP in a race against the giants (GPyTorch) on two different tracks: an Apple M4 laptop and an NVIDIA RTX 3060 graphics card.

On the Apple Laptop: LightGP was like a sprinter compared to a marathon runner. Because Apple's chips have a special "assistant" (called AMX) that handles math very fast, LightGP tapped into it directly. GPyTorch, however, had to run through a long line of Python managers and dispatchers before it could even ask the assistant for help.
- Result: LightGP was 2.6 to 8.7 times faster for standard tasks.
On the NVIDIA Graphics Card: LightGP was also faster for small-to-medium tasks (up to about 2,000 data points). However, for very large tasks, the "giant" (GPyTorch) started to catch up because it has built-in features to handle massive, continuous streams of data better.
- Result: LightGP was 2.3 to 6.7 times faster for smaller jobs, but GPyTorch won on the very largest jobs.

3. The "Magic" Tricks

The paper highlights two specific "magic tricks" LightGP uses:

The Invisible Matrix: Usually, to do the math, you have to write down a giant grid of numbers (the matrix) that takes up a lot of memory. LightGP has a trick where it calculates the result without ever writing down the grid. It's like calculating the total cost of a shopping trip without writing down every single item's price on a piece of paper first. This allows it to handle huge amounts of data on computers with limited memory.
The Speed Boost: On Apple computers, the authors found that using the "CPU" (the main brain) was actually faster than the "GPU" (the graphics card) for certain heavy math tasks because the CPU had a specialized assistant (AMX) that the GPU couldn't beat. LightGP automatically knows when to switch to the CPU to save time.

4. The Bottom Line

LightGP is a new tool that proves you don't need a massive, heavy software framework to do advanced statistical predictions.

It is fast: It runs 2 to 8 times faster than the popular alternatives on many devices.
It is light: It has no heavy dependencies and can be installed with a single command (pip install lightgp).
It is flexible: It works on both Apple and NVIDIA hardware.

The authors conclude that for many everyday data sizes (up to 100,000 points), this lightweight tool is the most efficient way to get accurate predictions, especially on Apple devices, without the bloat of the older, heavier tools.

Technical Summary: LightGP

Problem Statement

Gaussian Process (GP) inference is a principled Bayesian framework for regression with calibrated uncertainty, but its computational cost— $O(N^3)$ for exact inference and $O(N^2)$ memory—has historically limited its application to small datasets. While scalable approximations (e.g., inducing-point methods, conjugate gradients, structured kernel interpolation) have extended these boundaries to hundreds of thousands of points, dominant implementations like GPyTorch and GPflow rely on heavy deep-learning frameworks (PyTorch, TensorFlow). These frameworks impose significant dependency footprints (often exceeding 2 GB), non-trivial Python dispatch overhead per operation, and limited support for deployment in mobile or embedded environments. Furthermore, existing libraries often lack native, dependency-free C++ implementations that can leverage specific hardware accelerators like Apple Metal or NVIDIA CUDA without framework bloat.

Methodology

LightGP is a dependency-free C++17 library designed to address these limitations. It provides Python bindings but operates on a native C++ core, eliminating the need for PyTorch or TensorFlow. The library is structured into four layers:

Tensor and BLAS: A Tensor class storing row-major float32 data. It dispatches matrix multiplication to Apple Accelerate (using the AMX coprocessor on macOS) or OpenBLAS (on Linux). Cholesky factorization utilizes LAPACK spotrf.
Kernels: An abstract Kernel base class supports composition via operators (+, *, Scale). Concrete implementations include RBF, Matérn, Periodic, and Linear kernels.
GPU Backends:
- Metal (Apple Silicon): Implements tiled kernel matrix construction, a fused matrix-free kernel-vector product, blocked Cholesky with tiled threadgroup-memory GEMM, and triangular solves.
- CUDA (NVIDIA): Mirrors Metal functionality using cuBLAS, cuSOLVER, and custom kernels for matrix-free operations.
- Auto-Selection: A Backend::Auto heuristic selects CPU or GPU based on problem size and dimensionality.
Inference Models:
- GPExact: Supports three solvers: Exact Cholesky ( $O(N^3)$ ), Conjugate Gradients (CG) with stochastic Lanczos log-determinant (matrix-free, $O(N^2k)$ ), and Structured Kernel Interpolation (SKI) via FFT ( $O(N \log N)$ ).
- GPSparse: Implements Titsias' Variational Free Energy (VFE) with farthest-point inducing initialization and Adam optimization. It includes a "warm-refit" feature to skip re-initialization on repeated fits.

The library compiles as a single static library with zero external dependencies and is installable via pip install lightgp.

Key Contributions

Native C++ Implementation: A complete implementation of four GP inference paths (exact, CG, sparse VFE, SKI) with composable kernels, deployable without Python or framework dependencies.
Hardware-Optimized Kernels: Hand-tuned Metal compute shaders for Apple Silicon and CUDA kernels for NVIDIA GPUs. Notably, it introduces a fused matrix-free kernel-vector product that enables $O(N)$ memory inference for CG, a feature not natively supported by PyTorch's Metal backend.
FFT-Accelerated SKI: An implementation of SKI using Apple Accelerate vDSP and cuFFT to achieve $O(N \log N)$ matrix-vector multiplication at scale.
Cross-Platform Benchmarks: Systematic evaluation revealing that on Apple Silicon, the CPU-based AMX coprocessor (via Accelerate) significantly outperforms Metal for dense Cholesky factorization, a finding that informs backend selection heuristics.
Zero-Dependency Python Bindings: Full API exposure via pybind11 with NumPy interoperability, allowing installation without the heavy footprint of deep-learning frameworks.

Results

Benchmarks were conducted on an Apple M4 (Metal 3) and a Linux desktop with an Intel i7-12700 and NVIDIA RTX 3060 (CUDA 12.0).

Performance vs. GPyTorch:
- Apple M4 (CPU): LightGP is 2.6–8.7× faster than GPyTorch for exact GP inference and ~1.5× faster for sparse GP across all tested scales. The speedup is attributed to eliminating Python and PyTorch dispatch overhead.
- Linux (CUDA): LightGP is 2.3–6.7× faster than GPyTorch for exact GP up to $N=2,048$ . At $N=4,096$ , GPyTorch closes the gap due to its fused autograd pipeline amortizing overhead. For sparse GP, LightGP is faster at $N=10,000$ (1.7×), but GPyTorch wins at $N=50,000$ due to persistent device tensors and compiled gradient computation.
Matrix-Free Operations: On Apple M4, the fused Metal matrix-free kernel-vector product achieves a 32× speedup over the explicit CPU path at $N=20,000$ and is the only feasible option beyond $N=30,000$ on 16 GB memory. On Linux, CUDA achieves up to 259× speedup at $N=10,000$ .
Memory Scaling: Matrix-free CG, sparse VFE, and SKI maintain $O(N)$ memory usage. In contrast, exact Cholesky ( $O(N^2)$ ) exhausts 16 GB RAM at $N \approx 42,000$ –45,000 on Mac and 12 GB VRAM at $N \approx 35,000$ on Linux.
Accuracy: Predictions from CPU and GPU backends agree within $10^{-5}$ for mean and variance across standard datasets (Motorcycle, Mauna Loa, kin40k).

Significance

LightGP demonstrates that native C++ GP inference can match or exceed the performance of framework-based libraries while requiring zero external dependencies. By removing the dispatch overhead of Python and deep-learning frameworks, LightGP offers consistent performance gains, particularly on Apple Silicon where the AMX coprocessor is leveraged directly. The library enables GP inference at scales ( $N > 100,000$ ) previously impractical on consumer hardware through $O(N)$ memory matrix-free methods and FFT-accelerated SKI.

The paper highlights a critical hardware insight: on Apple Silicon, the CPU (via Accelerate/AMX) is often superior to the GPU (Metal) for dense linear algebra tasks like Cholesky factorization, a nuance that LightGP's heuristics automatically manage. While GPyTorch remains competitive or superior at very large scales ( $N \ge 50,000$ ) on CUDA due to its mature autograd and persistent tensor infrastructure, LightGP fills a niche for lightweight, high-performance, and deployable GP inference without the burden of a full deep-learning stack.