python - 如何绘制 Python 3 维水平集？

Question

我在绘制脑海中的图像时遇到了一些麻烦。
我想用支持向量机可视化内核技巧。所以我做了一些由两个圆(一个内圆和一个外圆)组成的二维数据，这两个圆应该被一个超平面分开。显然这在二维中是不可能的 - 所以我将它们转换成 3D。令 n 为样本数。现在我有一个 (n,3)-array(3 列，n 行)X 数据点和一个 (n,1)-array y 带标签。使用 sklearn 我通过获得线性分类器

clf = svm.SVC(kernel='linear', C=1000)
clf.fit(X, y)

我已经通过以下方式将数据点绘制为散点图

plt.scatter(X[:, 0], X[:, 1], c=y, s=30, cmap=plt.cm.Paired)

现在我想将分离超平面绘制为曲面图。我的问题是缺少超平面的显式表示，因为决策函数仅通过 decision_function = 0 产生隐式超平面.因此，我需要绘制 4 维对象的级别集(级别 0)。
由于我不是 Python 专家，如果有人能帮助我，我将不胜感激!而且我知道这并不是使用 SVM 的真正“风格”，但我需要这张图片作为我论文的插图。
编辑:我目前的“代码”

import numpy as np
import matplotlib.pyplot as plt
from sklearn import svm
from sklearn.datasets import make_blobs, make_circles

from tikzplotlib import save as tikz_save

plt.close('all')

# we create 50 separable points
#X, y = make_blobs(n_samples=40, centers=2, random_state=6)
X,  y  = make_circles(n_samples=50, factor=0.5, random_state=4, noise=.05)
X2, y2 = make_circles(n_samples=50, factor=0.2, random_state=5, noise=.08)

X = np.append(X,X2, axis=0)
y = np.append(y,y2, axis=0)

# shifte X to [0,2]x[0,2]
X = np.array([[item[0] + 1, item[1] + 1] for item in X])
X[X<0] = 0.01

clf = svm.SVC(kernel='rbf', C=1000)
clf.fit(X, y)

plt.scatter(X[:, 0], X[:, 1], c=y, s=30, cmap=plt.cm.Paired)


# plot the decision function
ax = plt.gca()
xlim = ax.get_xlim()
ylim = ax.get_ylim()

# create grid to evaluate model
xx = np.linspace(xlim[0], xlim[1], 30)
yy = np.linspace(ylim[0], ylim[1], 30)
YY, XX = np.meshgrid(yy, xx)
xy = np.vstack([XX.ravel(), YY.ravel()]).T
Z = clf.decision_function(xy).reshape(XX.shape)

# plot decision boundary and margins
ax.contour(XX, YY, Z, colors='k', levels=[-1, 0, 1], alpha=0.5, linestyles=['--','-','--'])
# plot support vectors
ax.scatter(clf.support_vectors_[:, 0], clf.support_vectors_[:, 1], s=100,
           linewidth=1, facecolors='none', edgecolors='k')


################## KERNEL TRICK - 3D ##################

trans_X = np.array([[item[0]**2, item[1]**2, np.sqrt(2*item[0]*item[1])] for item in X])

fig = plt.figure()
ax = plt.axes(projection ="3d")
 
# creating scatter plot
ax.scatter3D(trans_X[:,0],trans_X[:,1],trans_X[:,2], c = y, cmap=plt.cm.Paired)

clf2 = svm.SVC(kernel='linear', C=1000)
clf2.fit(trans_X, y)


ax = plt.gca(projection='3d')
xlim = ax.get_xlim()
ylim = ax.get_ylim()
zlim = ax.get_zlim()

### from here i don't know what to do ###
xx = np.linspace(xlim[0], xlim[1], 3)
yy = np.linspace(ylim[0], ylim[1], 3)
zz = np.linspace(zlim[0], zlim[1], 3)
ZZ, YY, XX = np.meshgrid(zz, yy, xx)
xyz = np.vstack([XX.ravel(), YY.ravel(), ZZ.ravel()]).T
Z = clf2.decision_function(xyz).reshape(XX.shape)
#ax.contour(XX, YY, ZZ, Z, colors='k', levels=[-1, 0, 1], alpha=0.5, linestyles=['--','-','--'])

期望输出
我想得到类似 that 的东西.
一般来说，我想重建他们在 this article 中所做的事情，尤其是“非线性变换”。

Answer 1

您的部分问题已在 this question on linear-kernel SVM 中得到解决.这是一个部分答案，因为只有线性内核可以用这种方式表示，即多亏了使用线性内核时可以通过估计器访问的超平面坐标。
另一种解决方案是使用 marching_cubes 找到等值面。
此解决方案涉及安装 scikit-image工具包 ( https://scikit-image.org ) 允许从 3D 坐标的网格中找到给定值的等值面(这里，我考虑了 0，因为它代表到超平面的距离)。
在下面的代码(从您的代码中复制)中，我为任何内核实现了这个想法(在示例中，我使用了 RBF 内核)，输出显示在代码下方。请考虑我关于使用 matplotlib 进行 3D 绘图的脚注，这可能是您的情况的另一个问题。

import numpy as np
import matplotlib.pyplot as plt
from sklearn import svm
from skimage import measure
from sklearn.datasets import make_blobs, make_circles

from tikzplotlib import save as tikz_save
from mpl_toolkits.mplot3d.art3d import Poly3DCollection

plt.close('all')

# we create 50 separable points
#X, y = make_blobs(n_samples=40, centers=2, random_state=6)
X,  y  = make_circles(n_samples=50, factor=0.5, random_state=4, noise=.05)
X2, y2 = make_circles(n_samples=50, factor=0.2, random_state=5, noise=.08)

X = np.append(X,X2, axis=0)
y = np.append(y,y2, axis=0)

# shifte X to [0,2]x[0,2]
X = np.array([[item[0] + 1, item[1] + 1] for item in X])
X[X<0] = 0.01

clf = svm.SVC(kernel='rbf', C=1000)
clf.fit(X, y)

plt.scatter(X[:, 0], X[:, 1], c=y, s=30, cmap=plt.cm.Paired)


# plot the decision function
ax = plt.gca()
xlim = ax.get_xlim()
ylim = ax.get_ylim()

# create grid to evaluate model
xx = np.linspace(xlim[0], xlim[1], 30)
yy = np.linspace(ylim[0], ylim[1], 30)
YY, XX = np.meshgrid(yy, xx)
xy = np.vstack([XX.ravel(), YY.ravel()]).T
Z = clf.decision_function(xy).reshape(XX.shape)

# plot decision boundary and margins
ax.contour(XX, YY, Z, colors='k', levels=[-1, 0, 1], alpha=0.5, linestyles=['--','-','--'])
# plot support vectors
ax.scatter(clf.support_vectors_[:, 0], clf.support_vectors_[:, 1], s=100,
           linewidth=1, facecolors='none', edgecolors='k')


################## KERNEL TRICK - 3D ##################

trans_X = np.array([[item[0]**2, item[1]**2, np.sqrt(2*item[0]*item[1])] for item in X])

fig = plt.figure()
ax = plt.axes(projection ="3d")

# creating scatter plot
ax.scatter3D(trans_X[:,0],trans_X[:,1],trans_X[:,2], c = y, cmap=plt.cm.Paired)

clf2 = svm.SVC(kernel='rbf', C=1000)
clf2.fit(trans_X, y)


z = lambda x,y: (-clf2.intercept_[0]-clf2.coef_[0][0]*x-clf2.coef_[0][1]*y) / clf2.coef_[0][2]

ax = plt.gca(projection='3d')
xlim = ax.get_xlim()
ylim = ax.get_ylim()
zlim = ax.get_zlim()

### from here i don't know what to do ###
xx = np.linspace(xlim[0], xlim[1], 50)
yy = np.linspace(ylim[0], ylim[1], 50)
zz = np.linspace(zlim[0], zlim[1], 50)
XX ,YY, ZZ = np.meshgrid(xx, yy, zz)
xyz = np.vstack([XX.ravel(), YY.ravel(), ZZ.ravel()]).T
Z = clf2.decision_function(xyz).reshape(XX.shape)

# find isosurface with marching cubes
dx = xx[1] - xx[0]
dy = yy[1] - yy[0]
dz = zz[1] - zz[0]
verts, faces, _, _ = measure.marching_cubes_lewiner(Z, 0, spacing=(1, 1, 1), step_size=2)
verts *= np.array([dx, dy, dz])
verts -= np.array([xlim[0], ylim[0], zlim[0]])

# add as Poly3DCollection
mesh = Poly3DCollection(verts[faces])
mesh.set_facecolor('g')
mesh.set_edgecolor('none')
mesh.set_alpha(0.3)
ax.add_collection3d(mesh)
ax.view_init(20, -45)
plt.savefig('kerneltrick')

运行代码使用 Matplotlib 生成以下图像，其中绿色半透明表面表示非线性决策边界。

脚注:使用 matplotlib 进行 3D 绘图
请注意，在某些情况下，Matplotlib 3D 无法管理对象的“深度”，因为它可能与 zorder 冲突。这个对象的。这就是为什么有时超平面看起来是在点的“顶部”绘制的原因，即使它应该在“后面”。此问题是 the matplotlib 3d documentation 中讨论的已知错误。并在 this answer .
如果你想有更好的渲染效果，你可能需要使用 Mayavi ，按照 Matplotlib 开发人员的推荐，或任何其他 3D Python 绘图库。