Cómo calcular un núcleo gaussiano efectivamente en numpy [cerrado]

Tengo una matriz numpy con m columnas yn filas, siendo las columnas dimensiones y los puntos de datos de las filas.

Ahora necesito calcular los valores del núcleo para cada combinación de puntos de datos.

Para una lineal del núcleo $K(\mathbf{x}_i,\mathbf{x}_j) = \langle \mathbf{x}_i,\mathbf{x}_j \rangle$ puedo simplemente hacerdot(X,X.T)

$K(\mathbf{x}_i,\mathbf{x}_j) = \exp{-\frac{\|\mathbf{x}_i-\mathbf{x}_j\|_2^2}{s^2}}$

Creo que el problema principal es obtener las distancias por pares de manera eficiente. Una vez que tienes eso, el resto es elemento sabio.

Para hacer esto, probablemente quieras usar scipy. La función scipy.spatial.distance.pdisthace lo que necesita y scipy.spatial.distance.squareformposiblemente le facilitará la vida.

Entonces, si quieres la matriz del núcleo, haz

from scipy.spatial.distance import pdist, squareform
  # this is an NxD matrix, where N is number of items and D its dimensionalites
X = loaddata() 
pairwise_dists = squareform(pdist(X, 'euclidean'))
K = scip.exp(-pairwise_dists ** 2 / s ** 2)

La documentación se puede encontrar aquí . 

Como una pequeña adición a la respuesta de bayerj, la pdistfunción de scipy puede calcular directamente las normas euclidianas cuadráticas llamándolas como pdist(X, 'sqeuclidean'). El código completo se puede escribir de manera más eficiente como

from scipy.spatial.distance import pdist, squareform
  # this is an NxD matrix, where N is number of items and D its dimensionalites
X = loaddata() 
pairwise_sq_dists = squareform(pdist(X, 'sqeuclidean'))
K = scip.exp(-pairwise_sq_dists / s**2)

También puedes escribir forma cuadrada a mano:

import numpy as np
def vectorized_RBF_kernel(X, sigma):
    # % This is equivalent to computing the kernel on every pair of examples
    X2 = np.sum(np.multiply(X, X), 1) # sum colums of the matrix
    K0 = X2 + X2.T - 2 * X * X.T
    K = np.power(np.exp(-1.0 / sigma**2), K0)
    return K

PD pero esto funciona un 30% más lento

def my_kernel(X,Y):
    K = np.zeros((X.shape[0],Y.shape[0]))
    for i,x in enumerate(X):
        for j,y in enumerate(Y):
            K[i,j] = np.exp(-1*np.linalg.norm(x-y)**2)
    return K

clf=SVR(kernel=my_kernel)

que es igual a

clf=SVR(kernel="rbf",gamma=1)

Puede calcular efectivamente el RBF a partir de la nota del código anterior que indica que el valor gamma es 1, ya que es una constante, la que solicitó también es la misma constante.

Creo que esto ayudará:

def GaussianKernel(v1, v2, sigma):
    return exp(-norm(v1-v2, 2)**2/(2.*sigma**2))