高斯核——Python实现

import numpy as np
from scipy.spatial.distance import cdist
from scipy.optimize import minimize

RBF kernel

$$k(x,y)=exp(-\frac{||x-y|{|}^2}{2{\sigma }^2})$$

The function rbf_kernel computes the radial basis function (RBF) kernel between two vectors. This kernel is defined as:
$$k(x,y)=exp(-\gamma ||x-y|{|}^2)$$
where $x$ and $y$ are the input vectors. If $\gamma ={\sigma }^{-2}$the kernel is known as the Gaussian kernel of variance${\sigma }^{-2}$.

The euclidean distance between a pair of rowvector x and y is computed as:
$$dist(x,y)=sqrt(dot(x,x)-2\ast dot(x,y)+dot(y,y))$$
This formulation has two advantages over other ways of computing distances.

1. First, it is computationally efficient when dealing with sparse data.
2. Second, if one argument varies but the other remains unchanged, then dot(x, x) and/or dot(y, y) can be pre-computed.

There are four ways to obtain kernel matrixs.

• According to the theory.

def rbf_kernel_1(X,Y,gamma = 0.1):
    """
    param X: ndarray of shape (n_samples_X,n_features)
    param Y: ndarray of shape (n_samples_Y,n_features)
    param gamma: if None, default to 0.1, refer to 1 / n_features
    Return: kernel_matrix: ndarray of shape(n_samples_X,n_samples_Y)
    """
    #compute X_norm_squared,Y_norm_squared
    if X.ndim == 1: #dim of X, if dim_x = 1
        X_norm_squared = X **2
        Y_norm_squared = Y **2
    else: #dim >= 2
        X_norm_squared = (X **2).sum(axis = 1).reshape(X.shape[0],1)
        Y_norm_squared = (Y **2).sum(axis = 1).reshape(Y.shape[0],1)
        #X_norm_squared : array-like of shape (n_samples_X,)
        #Y_norm_squared : array-like of shape (n_samples_Y,)
    squared_Euclidean_distances = Y_norm_squared[:,] + X_norm_squared.T - 2 * np.dot(Y,X.T)
    return np.exp(-squared_Euclidean_distances * gamma)

• by using scipy distance cdist to get Euclidean distances

def rbf_kernel_2(X,Y,gamma = 0.1):
    if X.ndim==X.ndim and X.ndim==2: # both matrices
        return np.exp(-gamma * cdist(X,Y)**2)
    else: # both vectors or a vector and a matrix
        return np.exp(- gamma * ( dot(X,X.T) + dot(Y,Y.T)- 2*dot(X,Y)))

• The most clear and rough way.

def rbf_kernel_3(X,Y, gamma = 0.1):
    dist_matrix = np.sum(X**2, 1).reshape(-1, 1) + np.sum(Y**2, 1) - 2 * np.dot(X,Y.T)
    return np.exp(-gamma * dist_matrix)

• rbf_kernel reference from sklearn

from sklearn.metrics.pairwise import rbf_kernel
def rbf_kernel_4(X,Y,gamma = 0.1):
    return rbf_kernel(X,Y,gamma)

Having a test as follows.

if "__name__" == "main":
    n = 100
    dim = 1
    X = np.array(np.linspace(1,10,n)).reshape(n,dim)
    Y = np.array(np.linspace(1,10,n)).reshape(n,dim)
print(rbf_kernel_1(X,Y,gamma = 0.1))
print("----------------------------------------------------")
print(rbf_kernel_2(X,Y,gamma = 0.1))
print("----------------------------------------------------")
print(rbf_kernel_3(X,Y,gamma = 0.1))
print("----------------------------------------------------")
print(rbf_kernel_4(X,Y,gamma = 0.1))

[[1.00000000e+00 9.99173895e-01 9.96699673e-01 ... 4.19673698e-04
  3.57208797e-04 3.03539138e-04]
 [9.99173895e-01 1.00000000e+00 9.99173895e-01 ... 4.92247497e-04
  4.19673698e-04 3.57208797e-04]
 [9.96699673e-01 9.99173895e-01 1.00000000e+00 ... 5.76417873e-04
  4.92247497e-04 4.19673698e-04]
 ...
 [4.19673698e-04 4.92247497e-04 5.76417873e-04 ... 1.00000000e+00
  9.99173895e-01 9.96699673e-01]
 [3.57208797e-04 4.19673698e-04 4.92247497e-04 ... 9.99173895e-01
  1.00000000e+00 9.99173895e-01]
 [3.03539138e-04 3.57208797e-04 4.19673698e-04 ... 9.96699673e-01
  9.99173895e-01 1.00000000e+00]]
----------------------------------------------------
[[1.00000000e+00 9.99173895e-01 9.96699673e-01 ... 4.19673698e-04
  3.57208797e-04 3.03539138e-04]
 [9.99173895e-01 1.00000000e+00 9.99173895e-01 ... 4.92247497e-04
  4.19673698e-04 3.57208797e-04]
 [9.96699673e-01 9.99173895e-01 1.00000000e+00 ... 5.76417873e-04
  4.92247497e-04 4.19673698e-04]
 ...
 [4.19673698e-04 4.92247497e-04 5.76417873e-04 ... 1.00000000e+00
  9.99173895e-01 9.96699673e-01]
 [3.57208797e-04 4.19673698e-04 4.92247497e-04 ... 9.99173895e-01
  1.00000000e+00 9.99173895e-01]
 [3.03539138e-04 3.57208797e-04 4.19673698e-04 ... 9.96699673e-01
  9.99173895e-01 1.00000000e+00]]
----------------------------------------------------
[[1.00000000e+00 9.99173895e-01 9.96699673e-01 ... 4.19673698e-04
  3.57208797e-04 3.03539138e-04]
 [9.99173895e-01 1.00000000e+00 9.99173895e-01 ... 4.92247497e-04
  4.19673698e-04 3.57208797e-04]
 [9.96699673e-01 9.99173895e-01 1.00000000e+00 ... 5.76417873e-04
  4.92247497e-04 4.19673698e-04]
 ...
 [4.19673698e-04 4.92247497e-04 5.76417873e-04 ... 1.00000000e+00
  9.99173895e-01 9.96699673e-01]
 [3.57208797e-04 4.19673698e-04 4.92247497e-04 ... 9.99173895e-01
  1.00000000e+00 9.99173895e-01]
 [3.03539138e-04 3.57208797e-04 4.19673698e-04 ... 9.96699673e-01
  9.99173895e-01 1.00000000e+00]]
----------------------------------------------------
[[1.00000000e+00 9.99173895e-01 9.96699673e-01 ... 4.19673698e-04
  3.57208797e-04 3.03539138e-04]
 [9.99173895e-01 1.00000000e+00 9.99173895e-01 ... 4.92247497e-04
  4.19673698e-04 3.57208797e-04]
 [9.96699673e-01 9.99173895e-01 1.00000000e+00 ... 5.76417873e-04
  4.92247497e-04 4.19673698e-04]
 ...
 [4.19673698e-04 4.92247497e-04 5.76417873e-04 ... 1.00000000e+00
  9.99173895e-01 9.96699673e-01]
 [3.57208797e-04 4.19673698e-04 4.92247497e-04 ... 9.99173895e-01
  1.00000000e+00 9.99173895e-01]
 [3.03539138e-04 3.57208797e-04 4.19673698e-04 ... 9.96699673e-01
  9.99173895e-01 1.00000000e+00]]