Technical Summary: Random features for Grassmannian kernel approximation with bounded rank-one projections

Problem Statement

In many supervised and unsupervised learning tasks, data classes or clusters are naturally represented by low-dimensional subspaces in high-dimensional spaces, rather than individual Euclidean points. This necessitates operating on the Grassmannian manifold $\mathcal{G}(k, n)$, the set of all $k$-dimensional subspaces in $\mathbb{R}^n$. While kernel machines (e.g., SVMs) can effectively learn on such data using positive-definite kernels (like the projection or Binet–Cauchy kernels), they suffer from poor scalability. Standard kernel methods require constructing and storing an $N \times N$ Gram matrix, leading to $O(N^2)$ memory and computational costs, which is prohibitive for large datasets.

Existing attempts to approximate Grassmannian kernels via random features face two primary limitations:

1. Dense Random Projections: Embedding subspaces by projecting their orthogonal projectors onto dense random Gaussian matrices requires $O(n^2)$ storage and computation per feature, making it impractical for high dimensions.
2. Unbounded Rank-One Projections (ROP): While computationally cheaper ($O(kn)$ per feature), ROPs based on quadratic measurements $a^\top P b$ exhibit heavy-tailed distributions (sub-Weibull). This results in slow concentration rates, preventing uniform approximation guarantees over the entire Grassmannian manifold with a number of features scaling linearly with the intrinsic dimension.

Methodology

The authors propose a family of random feature maps that combine rank-one projections (ROP) with bounded non-linear transforms to achieve scalable, uniform kernel approximation.

1. Rank-One Projections (ROP)

For a subspace $P \in \mathcal{G}(k, n)$ represented by an orthonormal basis $U$ (where $P = UU^\top$), the method probes the subspace using $m$ independent rank-one quadratic measurements:
$$\psi_{\text{rop}}(U) = \left( a_i^\top P b_i \right)_{i=1}^m$$
where $a_i, b_i \sim \text{i.i.d. } \mathcal{N}(0, I_n)$. This reduces the per-feature computational cost to $O(kn)$ and storage to $O(kn)$ (after factoring $U^\top a_i$ and $U^\top b_i$).

2. Bounded Non-Linear Transforms

To tame the heavy-tailed distribution of the unbounded ROPs and ensure uniform convergence, the authors compose $\psi_{\text{rop}}$ with bounded non-linear functions $g(\cdot)$:
$$\psi_g(U) = g(\psi_{\text{rop}}(U))$$
Two specific transforms are analyzed:

• Binary Features: Using the sign function, $g(\lambda) = \text{sign}(\lambda)$. This yields binary random features $\psi_\pm \in \{\pm 1\}^m$, enabling compact 1-bit storage and fast bitwise operations.
• Periodic Features: Using the complex exponential, $g(\lambda) = \exp(i\omega \lambda)$. This yields periodic random features $\psi_\circ$, analogous to Random Fourier Features.

3. Structured Projections

To further reduce computational complexity, the authors replace the dense Gaussian probing vectors with structured random matrices based on compositions of fast Walsh-Hadamard transforms and random diagonal sign matrices (inspired by Choromanski et al., 2016). This reduces the computational complexity from $O(kmn)$ to $O(km \log n)$ while maintaining statistical properties.

Key Contributions

Theoretical Guarantees

The paper establishes that the empirical inner products of these bounded random features are unbiased estimators of specific, well-defined, rotation-invariant Grassmannian kernels. Crucially, the authors prove uniform approximation bounds over all subspaces of fixed dimension $k$ in $\mathbb{R}^n$.

• Sample Complexity: They show that $m = O(kn)$ random features (up to logarithmic factors) suffice to control the approximation error uniformly. This scales linearly with the ambient dimension $n$ and the subspace dimension $k$, a significant improvement over the $O((kn)^2)$ requirements for unbounded ROPs.
• Binary Kernel ($\kappa_\pm$): The expected kernel for binary features is shown to be positive-definite and rotationally invariant. While no closed-form expression exists for general $k$, it is characterized semi-explicitly as:
  $$\kappa_\pm(U, V) = 1 - \frac{2}{\pi} \mathbb{E}[\angle(P g, Q g)]$$
  where $g \sim \mathcal{N}(0, I_n)$ and $\angle(\cdot, \cdot)$ is the angle between vectors. For $k=1$, a closed form $(1 - 2\theta/\pi)^2$ is derived.
• Periodic Kernel ($\kappa_\circ$): The expected kernel for periodic features admits a closed-form expression depending on the principal angles $\{\theta_j\}_{j=1}^k$:
  $$\kappa_\circ(U, V) = \prod_{j=1}^k (1 + \omega^2 \sin^2 \theta_j)^{-1}$$
  The authors analyze the behavior of $\kappa_\circ$ with respect to the frequency parameter $\omega$:
  • Large $\omega$: $\kappa_\circ$ behaves like an inverse Binet–Cauchy kernel (sensitive to small principal angles).
  • Small $\omega$: $\kappa_\circ$ approximates a projection Radial Basis Function (RBF) kernel, $\exp(-\frac{1}{2}\omega^2 \|P-Q\|_F^2)$.

Computational Complexity

The proposed methods offer substantial efficiency gains:

• Memory: The storage for unstructured bounded ROPs is $O(mn)$ for the probing vectors. When substituting $m = O(kn)$ to ensure uniform approximation, the overall memory complexity becomes $O(kn^2)$ (specifically derived as $O(\delta^{-2}(kn^2))$ in the paper). This is distinct from dense projections which require $O(mn^2)$ storage (scaling as $O(k^3 n^3)$ for uniform approximation).
• Computation: The computational complexity for unstructured bounded ROPs is $O(kmn)$. With $m = O(kn)$, this scales as $O(k^2 n^2)$. This is significantly lower than the $O(k^2 n^3)$ complexity associated with dense projection methods (which require $O(\delta^{-2} k^3 n^3)$).
• Storage: Binary features require only $O(m)$ bits per subspace.

Experimental Results

The authors validate their approach using synthetic experiments and the ETH-80 dataset (visual object images from different viewpoints).

1. Kernel Approximation: Experiments confirm that bounded ROPs converge to their target kernels with significantly lower computational and memory costs than dense projections.
2. Structured Embeddings: Empirical analysis of the structured Hadamard-diagonal sketches shows that using $S=3$ blocks is sufficient to approximate a Gaussian distribution, validating the use of fast transforms.
3. Classification Performance:
   • On superclass classification (8-way) and object classification (80-way) tasks, the proposed random features (both unstructured and structured) achieve classification accuracies comparable to exact kernel methods, though with some trade-offs depending on the specific task and kernel.
   • Specifically, structured embeddings with an oversampling ratio $\rho = m/(nk)$ as low as 0.05 (5% of the original basis size) yield high accuracy while reducing computation time by orders of magnitude compared to exact kernel evaluation. For instance, in the 80-way classification task, structured methods reduced total runtime from ~52.9 seconds (exact projection kernel) to ~1.9 seconds. However, this speedup came with a notable drop in accuracy from 98.75% (exact) to 83.00% (structured), highlighting the trade-off between computational efficiency and classification performance in this specific setting.

Significance and Claims

The paper claims to provide a practical and scalable alternative to classical Grassmannian kernels. By combining rank-one projections with bounded non-linearities, the authors demonstrate that it is possible to:

1. Recover Grassmannian geometry with high accuracy.
2. Achieve uniform approximation guarantees over the manifold with a number of features scaling linearly with the intrinsic dimension.
3. Drastically reduce memory and computational requirements, making kernel-based learning on subspace data feasible for large-scale, high-dimensional applications.

The work bridges the gap between theoretical kernel methods and practical scalability, offering a framework where the trade-off between approximation error and computational cost can be explicitly controlled via the number of random features and the choice of non-linear transform.