SVM - 如何对核化的 gram 矩阵进行矢量化?
Posted
技术标签:
【中文标题】SVM - 如何对核化的 gram 矩阵进行矢量化?【英文标题】:SVM - How can I vectorize a kernalized gram matrix? 【发布时间】:2018-06-05 07:43:37 【问题描述】:我使用 cvxopt qp 求解器在 python 中实现了一个支持向量机,我需要计算两个向量的 gram 矩阵,每个元素都有一个核函数。我使用 for 循环正确实现了它,但这种策略计算量很大。我想对代码进行矢量化处理。
例子:
这是我写的:
K = np.array( [kernel(X[i], X[j],poly=poly_kernel)
for j in range(m)
for i in range(m)]).reshape((m, m))
如何在没有 for 循环的情况下对上述代码进行矢量化以更快地获得相同的结果?
核函数计算一个高斯核。
Here is a quick explanation of an svm with kernel trick.第二页说明了问题。
这是我完整的code 上下文。
编辑:这是一个快速代码 sn-p,它运行我需要以非矢量化形式矢量化的内容
from sklearn.datasets import make_gaussian_quantiles;
import numpy as np;
X,y = make_gaussian_quantiles(mean=None, cov=1.0, n_samples=100, n_features=2, n_classes=2, shuffle=True, random_state=5);
m = X.shape[0];
def kernel(a,b,d=20,poly=True,sigma=0.5):
if (poly):
return np.inner(a,b) ** d;
else:
return np.exp(-np.linalg.norm((a - b) ** 2)/sigma**2)
# Need to vectorize these loops
K = np.array([kernel(X[i], X[j],poly=False)
for j in range(m)
for i in range(m)]).reshape((m, m))
谢谢!
【问题讨论】:
请通过包含相关位(在本例中为内核函数)和一些模拟数据来帮助那些试图帮助您的人,以便他们可以直接运行您的示例。 更新问题谢谢 【参考方案1】:这是一个矢量化版本。非多边形分支有两种变体,一种是直接的,另一种是节省内存的,以防特征数量很大:
from sklearn.datasets import make_gaussian_quantiles;
import numpy as np;
X,y = make_gaussian_quantiles(mean=None, cov=1.0, n_samples=100, n_features=2, n_classes=2, shuffle=True, random_state=5);
Y,_ = make_gaussian_quantiles(mean=None, cov=1.0, n_samples=200, n_features=2, n_classes=2, shuffle=True, random_state=2);
m = X.shape[0];
n = Y.shape[0]
def kernel(a,b,d=20,poly=True,sigma=0.5):
if (poly):
return np.inner(a,b) ** d;
else:
return np.exp(-np.linalg.norm((a - b) ** 2)/sigma**2)
# Need to vectorize these loops
POLY = False
LOW_MEM = 0
K = np.array([kernel(X[i], Y[j], poly=POLY)
for i in range(m)
for j in range(n)]).reshape((m, n))
def kernel_v(X, Y=None, d=20, poly=True, sigma=0.5):
Z = X if Y is None else Y
if poly:
return np.einsum('ik,jk', X, Z)**d
elif X.shape[1] < LOW_MEM:
return np.exp(-np.sqrt(((X[:, None, :] - Z[None, :, :])**4).sum(axis=-1)) / sigma**2)
elif Y is None or Y is X:
X2 = X*X
H = np.einsum('ij,ij->i', X2, X2) + np.einsum('ik,jk', X2, 3*X2) - np.einsum('ik,jk', X2*X, 4*X)
return np.exp(-np.sqrt(np.maximum(0, H+H.T)) / sigma**2)
else:
X2, Y2 = X*X, Y*Y
E = np.einsum('ik,jk', X2, 6*Y2) - np.einsum('ik,jk', X2*X, 4*Y) - np.einsum('ik,jk', X, 4*Y2*Y)
E += np.add.outer(np.einsum('ij,ij->i', X2, X2), np.einsum('ij,ij->i', Y2, Y2))
return np.exp(-np.sqrt(np.maximum(0, E)) / sigma**2)
print(np.allclose(K, kernel_v(X, Y, poly=POLY)))
【讨论】:
谢谢!你能解释一下 einsums 在做什么吗?einsum('ik,jk', X, X)
将X
的元素与X
的元素相乘,输出将具有与第一个X
中的第一个类似的轴,因为规范字符串中的i
是唯一的并且是轴就像第二个“X”中的第一个,因为“j”是唯一的。所以关于这两个轴,它的行为就像一个外部产品。这两个因素的第二个轴没有展开,并且单个产品沿该轴求和,因为“k”出现在两个规格中。那是爱因斯坦约定,因此是“einsum”,因此对于第二个轴,它的行为类似于内积。 out_ij = sum_k X_ik X_jk
如果在这种情况下还存在np.einsum('ij,ij->i', X2, X2)
中的输出规范,这意味着i
没有被求和,也没有被外部扩展,所以out_i = sum_j X2_ij X2_ij
我希望你能弄清楚其他人自己。
感谢您的解释。为什么要乘以 3*X2 和 4*X?
只是一些诡计(a-b)^4=a^4 - 4a^3 b + 6 a^2 b^2 - 4 a b^3 + b^4
。因为我们这样做“外部”并且 a 和 b 的向量都是 X
我们可以计算前半部分,然后后半部分将只是转置。由于术语的数量是奇数,我们需要将中间的一分为二,因此3
。以上是关于SVM - 如何对核化的 gram 矩阵进行矢量化?的主要内容,如果未能解决你的问题,请参考以下文章
如何获得在更大矩阵上训练的 SVM 以对不同大小的矩阵进行分类