Python Forum

from scipy import exp
from scipy.linalg import eigh
from scipy.spatial.distance import pdist, squareform
from sklearn.datasets import make_moons
from sklearn.decomposition import PCA       
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd

def rbf_kernel_pca(X, gamma, n_components):
    """
    RBF kernel PCA implementation.
    Parameters
    ------------
    X: {NumPy ndarray}, shape = [n_samples, n_features]
    gamma: float
    Tuning parameter of the RBF kernel
    n_components: int
    Number of principal components to return
    Returns
    ------------
    X_pc: {NumPy ndarray}, shape = [n_samples, k_features]
    Projected dataset
    """

    sq_dists = pdist(X, 'sqeuclidean')
    # Convert pairwise distances into a square matrix.
    mat_sq_dists = squareform(sq_dists)
    # Compute the symmetric kernel matrix.
    K = exp(-gamma * mat_sq_dists)
    # Center the kernel matrix.
    N = K.shape[0]
    one_n = np.ones((N,N)) / N
    K = K - one_n.dot(K) - K.dot(one_n) + one_n.dot(K).dot(one_n)
    # Obtaining eigenpairs from the centered kernel matrix
    # numpy.eigh returns them in sorted order
    eigvals, eigvecs = eigh(K)
    # Collect the top k eigenvectors (projected samples)
    X_pc = np.column_stack((eigvecs[:, -i] for i in range(1, n_components + 1)))
    return X_pc

X, y = make_moons(n_samples=100, random_state=123)

from matplotlib.ticker import FormatStrFormatter
fig, ax = plt.subplots(nrows=3, ncols=1, figsize=(8,5))

X_kpca = rbf_kernel_pca(X, gamma=0.1, n_components=2)

ax[0].scatter(X_kpca[y==0, 0], X_kpca[y==0, 1], color='red', marker='^', alpha=0.5)
ax[0].scatter(X_kpca[y==1, 0], X_kpca[y==1, 1], color='blue', marker='o', alpha=0.5)
ax[0].set_xlim([-0.6, 0.6])
ax[0].set_ylim([-0.3, 0.3])
ax[0].set_xlabel('PC1')
ax[0].set_ylabel('PC2')
plt.title('gamma=0.1')

X_kpca = rbf_kernel_pca(X, gamma=15, n_components=2)

ax[1].scatter(X_kpca[y==0, 0], X_kpca[y==0, 1], color='red', marker='^', alpha=0.5)
ax[1].scatter(X_kpca[y==1, 0], X_kpca[y==1, 1], color='blue', marker='o', alpha=0.5)
ax[1].set_xlim([-0.4, 0.4])
ax[1].set_ylim([-0.3, 0.3])
ax[1].set_xlabel('PC1')
ax[1].set_ylabel('PC2')
plt.title('gamma=15')

X_kpca = rbf_kernel_pca(X, gamma=150, n_components=2)

ax[2].scatter(X_kpca[y==0, 0], X_kpca[y==0, 1], color='red', marker='^', alpha=0.5)
ax[2].scatter(X_kpca[y==1, 0], X_kpca[y==1, 1], color='blue', marker='o', alpha=0.5)
ax[2].set_xlim([-0.2, 0.2])
ax[2].set_ylim([-0.3, 0.3])
ax[2].set_xlabel('PC1')
ax[2].set_ylabel('PC2')
plt.title('gamma=150')

plt.tight_layout()
plt.show()