Kernel basis transforms

Exact, Nystrom, and random-Fourier feature maps

from _api_doc_utils import *

1 Where it fits

Group: Transforms

KernelBasis stores a training design and transforms new rows into kernel similarities against that basis. For a Gaussian kernel, the transformed feature for training row \(j\) is

\[ \phi_j(x) = \exp\{-\|x-x_j\|^2/h\}. \]

The resulting feature matrix can be fed into any downstream regression estimator.

2 Kernel definitions

The shared kernels are

\[ k_{\mathrm{Gaussian}}(x,z) = \exp\{-\|x-z\|_2^2/h\}, \qquad k_{\mathrm{linear}}(x,z)=x'z, \]

and

\[ k_{\mathrm{poly}}(x,z)=(x'z+c_0)^d. \]

Here the bandwidth \(h\) is the direct denominator of squared distance. It is not a Gaussian standard deviation and there is no extra factor of two.

KernelBasis stores all \(n\) training rows and returns the exact cross-kernel matrix

\[ \Phi(x)=K(x,X_{\mathrm{train}})\in\mathbb R^{n_{\mathrm{new}}\times n}. \]

It is a feature map only; a downstream linear estimator determines the eventual loss and regularization.

3 Approximate bases

NystromBasis chooses \(m\) landmark rows, forms \(K_{mm}\), and returns

\[ \Phi_{\mathrm{Nys}}(x) = K(x,X_m)(K_{mm}+\rho I)^{-1/2}. \]

The inverse square root uses a symmetric eigendecomposition and floors eigenvalues at \(10^{-14}\). With no seed, the first \(m\) rows are landmarks; with a seed, a reproducibly shuffled subset is used. Then \(\Phi_{\mathrm{Nys}}\Phi_{\mathrm{Nys}}'\) approximates the full kernel matrix.

RandomFourierFeatures is available only for the Gaussian kernel above. It draws

\[ \omega_j\sim N(0,2I/h), \qquad b_j\sim\operatorname{Uniform}(0,2\pi), \]

and returns

\[ \phi_j(x)=\sqrt{\frac2m}\cos(x'\omega_j+b_j). \]

Its feature inner products approximate the Gaussian kernel in expectation. An omitted seed still uses a fixed internal seed, so repeated fits are reproducible rather than independently randomized.

4 Implementation walkthrough

These three classes share kernel conventions but have deliberately different implementation boundaries.

4.1 Exact `KernelBasis`

fit() asks Linfa’s kernel builder to materialize the square training Gram matrix and stores both that dense matrix and an owned copy of the training rows. This is the only delegated numerical primitive on this page.
transform() is package-owned. It checks the incoming feature count and uses nested row loops to evaluate the selected Linfa kernel method between every new row and every stored training row. Thus training and out-of-sample transformations use the same kernel definition but not the same matrix-building path.
fit_transform() returns the already materialized training Gram matrix. summary() reads its diagonal rather than recomputing self-similarities.

4.2 `NystromBasis`

The landmark indices begin as 0..n. With a seed they are shuffled once; without a seed they remain in data order. The vector is truncated to \(m\) and then sorted, so summary() reports landmarks in original row order even when the subset was random.
Native cross-kernel loops form \(K_{mm}\). The implementation adds ridge to its diagonal, runs a symmetric eigendecomposition, floors each resulting eigenvalue at \(10^{-14}\), and reconstructs \(U\operatorname{diag}(\lambda_j^{-1/2})U'\). The floor is applied after the requested ridge shift.
A transform forms \(K(x,X_m)\) with the same nested kernel loops and right-multiplies by the stored inverse square root. It retains landmarks and the dense \(m\times m\) factor, but not the original full training sample.

4.3 `RandomFourierFeatures`

fit() seeds a Rust standard RNG from the supplied seed or the fixed hexadecimal constant 0xF0471E5. It draws standard-normal pairs with a Box-Muller transform, fills the \(p\times m\) frequency matrix in row-major logical order, and scales every draw by \(\sqrt{2/h}\).
It draws a separate length-\(m\) phase vector uniformly on \([0,2\pi)\). Neither frequencies nor phases depend on the observed values of \(X\); fitting uses the input only to validate a nonempty sample and record \(p\).
transform() computes one dense matrix product \(X\Omega\), adds the stored phase columnwise, applies cosine, and multiplies by \(\sqrt{2/m}\). The transform is reproducible and stateless after those random arrays have been stored.

The exact basis makes sample size the feature dimension. Nystrom makes selected observations the basis and pays an \(m^3\) eigendecomposition once. Random Fourier features make random spectral draws the basis and reduce transformation to matrix multiplication plus cosine. None centers or standardizes \(X\), so distance scale remains the caller’s responsibility.

5 Inference and performance

None of these transforms estimates a response model or supplies statistical inference. KernelBasis fitting and storage are \(O(n^2p)\) time and \(O(n^2+np)\) memory because it materializes the full Gram matrix. Transforming \(r\) new rows costs \(O(rnp)\) and returns \(n\) features.

NystromBasis costs \(O(m^2p+m^3)\) to fit and approximately \(O(rmp+rm^2)\) to transform, with \(O(mp+m^2)\) fitted storage. RandomFourierFeatures stores \(O(pm)\) random weights and transforms in \(O(rpm)\). Kernel and random-feature outputs are dense; downstream ridge or OLS can still become expensive when \(m\) or the exact training-basis width is large.

6 Python API

Constructors: cm.KernelBasis, cm.NystromBasis, and cm.RandomFourierFeatures

All three classes expose fit, transform, fit_transform, and summary. KernelBasis uses the full training sample as its basis. NystromBasis(n_components, ...) selects landmark rows and takes a kernel ridge for its inverse square root. RandomFourierFeatures(n_components, bandwidth, seed) draws a fixed-width Gaussian-kernel approximation.

print(inspect.signature(cm.KernelBasis))

(kernel='gaussian', bandwidth=0.5, coef0=1.0, degree=2.0)

cls = cm.KernelBasis
display(HTML(html_table(["Public method"], public_methods(cls))))

Public method
`fit(self, /, x)`
`fit_transform(self, /, x)`
`summary(self, /)`
`transform(self, /, x)`

for cls in (cm.NystromBasis, cm.RandomFourierFeatures):
    print(cls.__name__, inspect.signature(cls))
    display(HTML(html_table(["Public method"], public_methods(cls))))

NystromBasis (n_components, kernel='gaussian', bandwidth=0.5, coef0=1.0, degree=2.0, ridge=1e-10, seed=None)

Public method
`fit(self, /, x)`
`fit_transform(self, /, x)`
`summary(self, /)`
`transform(self, /, x)`

RandomFourierFeatures (n_components, bandwidth=0.5, seed=None)

Public method
`fit(self, /, x)`
`fit_transform(self, /, x)`
`summary(self, /)`
`transform(self, /, x)`

7 Minimal example

rng = np.random.default_rng(23)
x = rng.normal(size=(80, 2))
y = np.sin(x[:, 0]) + rng.normal(scale=0.1, size=80)
basis = cm.KernelBasis(kernel='gaussian', bandwidth=0.8)
z = basis.fit_transform(x)
reg = cm.Ridge(penalty=0.1)
reg.fit(z, y)
print(basis.summary())
print(reg.predict(basis.transform(x[:3])))

nystrom = cm.NystromBasis(16, kernel='gaussian', bandwidth=0.8, seed=23)
random_features = cm.RandomFourierFeatures(64, bandwidth=0.8, seed=23)
print(nystrom.fit_transform(x).shape)
print(random_features.fit_transform(x).shape)

{'kernel': 'gaussian', 'n_train': 80, 'n_features': 2, 'bandwidth': 0.8, 'coef0': 1.0, 'degree': 2.0, 'diagonal': array([1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1.,
       1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1.,
       1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1.,
       1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1.,
       1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1.])}
[ 0.53255368 -0.02676873  0.27494531]
(80, 16)
(80, 64)

8 `summary()` contract

The table below is generated by fitting the live class in this repository and then inspecting summary(). Shapes are shown because most values are plain NumPy arrays or scalars.

rng = np.random.default_rng(123)
x = rng.normal(size=(50, 2))
model = cm.KernelBasis(kernel='gaussian', bandwidth=0.8)
model.fit(x)
summary = model.summary()
display(HTML(html_table(["summary() key", "shape"], summary_shape_rows(summary))))

summary() key	shape
`kernel`	`()`
`n_train`	`()`
`n_features`	`()`
`bandwidth`	`()`
`coef0`	`()`
`degree`	`()`
`diagonal`	`(50,)`

nystrom = cm.NystromBasis(10, kernel='gaussian', bandwidth=0.8, seed=123)
nystrom.fit(x)
display(HTML(html_table(
    ["NystromBasis summary key", "shape"],
    summary_shape_rows(nystrom.summary()),
)))

random_features = cm.RandomFourierFeatures(32, bandwidth=0.8, seed=123)
random_features.fit(x)
display(HTML(html_table(
    ["RandomFourierFeatures summary key", "shape"],
    summary_shape_rows(random_features.summary()),
)))

NystromBasis summary key	shape
`kernel`	`()`
`n_components`	`()`
`n_features`	`()`
`bandwidth`	`()`
`coef0`	`()`
`degree`	`()`
`ridge`	`()`
`landmark_indices`	`(10,)`

RandomFourierFeatures summary key	shape
`kernel`	`()`
`n_components`	`()`
`n_features`	`()`
`bandwidth`	`()`
`weights`	`(2, 32)`
`bias`	`(32,)`