Implement Kernel PCA using the RBF (Radial Basis Function) kernel. Project data into a higher-dimensional feature space implicitly through the kernel trick, then perform PCA in that space.
import numpy as np
from typing import Dict
def rbf_kernel(X: np.ndarray, Y: np.ndarray, gamma: float) -> np.ndarray:
# ||x - y||^2 = ||x||^2 + ||y||^2 - 2 * x.y
X_sq = np.sum(X ** 2, axis=1).reshape(-1, 1)
Y_sq = np.sum(Y ** 2, axis=1).reshape(1, -1)
dist_sq = X_sq + Y_sq - 2 * X @ Y.T
return np.exp(-gamma * dist_sq)
def kernel_pca(
X: np.ndarray,
n_components: int = 2,
gamma: float = 1.0
) -> Dict:
n = X.shape[0]
# Compute kernel matrix
K = rbf_kernel(X, X, gamma)
# Center the kernel matrix
one_n = np.ones((n, n)) / n
K_centered = K - one_n @ K - K @ one_n + one_n @ K @ one_n
# Eigendecomposition
eigenvalues, eigenvectors = np.linalg.eigh(K_centered)
# Sort by descending eigenvalue
idx = np.argsort(eigenvalues)[::-1]
eigenvalues = eigenvalues[idx]
eigenvectors = eigenvectors[:, idx]
# Select top components and normalize
top_eigenvalues = eigenvalues[:n_components]
top_eigenvectors = eigenvectors[:, :n_components]
# Project: alpha_i / sqrt(lambda_i)
for i in range(n_components):
if top_eigenvalues[i] > 1e-10:
top_eigenvectors[:, i] /= np.sqrt(top_eigenvalues[i])
X_transformed = K_centered @ top_eigenvectors
return {
"transformed": X_transformed.tolist(),
"eigenvalues": top_eigenvalues.tolist(),
"n_components": n_components,
"gamma": gamma
}