如何在 python 中绘制多项式逻辑回归的决策边界?

how to plot the decision boundary of a polynomial logistic regression in python?

我查看了这个网站上的示例:https://scipython.com/blog/plotting-the-decision-boundary-of-a-logistic-regression-model/

我了解他们如何绘制线性特征向量的决策边界。但是,如果我申请

,我将如何绘制决策边界
from sklearn.preprocessing import PolynomialFeatures
...
poly = PolynomialFeatures(degree = 3, interaction_only=False, include_bias=False)
X_poly = poly.fit_transform(X)
# Fit the data to a logistic regression model.
clf = sklearn.linear_model.LogisticRegression()
clf.fit(X_poly, Y)

获得弯曲的决策边界? (我知道对于网站上的例子来说意义不大,但谈论它可能更容易)。

我试图通过叠加多项式图来绘制生成的多项式决策边界,但只得到如下奇怪的结果:

那么我怎样才能做一个弯曲的决策边界图呢?

编辑后的代码:

from sklearn.preprocessing import PolynomialFeatures
import numpy as np
import matplotlib.pyplot as plt
import sklearn.linear_model
plt.rc('text', usetex=True)
plt.figure(dpi=1200) 
pts = np.loadtxt(r'C:\Users\stefa\OneDrive\Desktop\linpts.txt')
X = pts[:,:2]
Y = pts[:,2].astype('int')
poly = PolynomialFeatures(degree = 2, interaction_only=False, include_bias=False)
X_poly = poly.fit_transform(X)
# Fit the data to a logistic regression model.
clf = sklearn.linear_model.LogisticRegression()
clf.fit(X_poly, Y)

# Retrieve the model parameters.
b = clf.intercept_[0]
w1, w2,w3,w4,w5 = clf.coef_.T


# In[]
def PolyCoefficients(x, coeffs):
    """ Returns a polynomial for ``x`` values for the ``coeffs`` provided.

    The coefficients must be in ascending order (``x**0`` to ``x**o``).
    """
    o = len(coeffs)
    print(f'# This is a polynomial of order {ord}.')
    y = 0
    for i in range(o):
        y += coeffs[i]*x**i
    return y

x = np.linspace(0, 9, 100)
coeffs = [b, w1, w2, w3, w4, w5]
plt.plot(x, PolyCoefficients(x, coeffs))
plt.show()


# In[]
# Calculate the intercept and gradient of the decision boundary.
c = -b/w2
m = -w1/w2

# Plot the data and the classification with the decision boundary.
xmin, xmax = -1, 2
ymin, ymax = -1, 2.5
xd = np.array([xmin, xmax])
yd = m*xd + c
#plt.plot(xd, yd, 'k', lw=1, ls='--')
plt.plot(x, PolyCoefficients(x, coeffs))
plt.fill_between(xd, yd, ymin, color='tab:blue', alpha=0.2)
plt.fill_between(xd, yd, ymax, color='tab:orange', alpha=0.2)

plt.scatter(*X[Y==0].T, s=8, alpha=0.5)
plt.scatter(*X[Y==1].T, s=8, alpha=0.5)
plt.xlim(xmin, xmax)
plt.ylim(ymin, ymax)
plt.ylabel(r'$x_2$')
plt.xlabel(r'$x_1$')

plt.show()

您的 PolyCoefficients 函数的输出是一个四阶多项式,由以下内容组成:

coeffs[0]*x^0 + coeffs[1]*x^1 + coeffs[2]*x^2 + coeffs[3]*x^3 + coeffs[i]*x^4

相反,您需要的是二阶多项式(由 sklearn.preprocessing.PolynomialFeatures 对象内的参数 degree = 2 指定),实际上是以下内容:

(coeffs[0]*x1^1) + (coeffs[1]*x2^1) + (coeffs[2]*x1^2 + 2*coeffs[3]*x1*x2 + coeffs[4]*x2^2)

只要您使用两个特征 x1 和 x2,此公式就有效,否则您将需要使用 x1x2、...、xN 和所有源自 (x1 + x2 + ... + xN)^2.

的术语

您可以找到更多详细信息和示例 here

让我生成一个演示数据。

from sklearn.preprocessing import PolynomialFeatures
import numpy as np
import matplotlib.pyplot as plt
import sklearn.linear_model

X = np.random.normal(size=(1000, 2))
Y = ((X[:,0] - X[:,1] + 0.4*X[:,0]*X[:,1] + 0.7*X[:,0]**2 - 0.8*X[:,1]**2 +
     np.random.normal(scale=0.1, size=(1000,))) >= 0).astype(int)

flg = (Y > 0)
plt.scatter(X[flg,0], X[flg,1], alpha=0.3, marker="o")
plt.scatter(X[~flg,0], X[~flg,1], alpha=0.3, marker="x")

除了随机性之外,数据看起来是这样的。

像您一样训练模型。

poly = PolynomialFeatures(degree = 2, interaction_only=False, include_bias=False)
X_poly = poly.fit_transform(X)
# Fit the data to a logistic regression model.
clf = sklearn.linear_model.LogisticRegression()
clf.fit(X_poly, Y)
print(poly.powers_)
#[[1 0]
# [0 1]
# [2 0]
# [1 1]
# [0 2]]

这告诉我们特征的排序为:x1, x2, x1^2, x1*x2, x2^2。 所以收集系数和截距,给它们起个直观的名字。

w1, w2, w11, w12, w22 = clf.coef_[0]
b = clf.intercept_[0]

根据定义,决策边界是一组 (x1, x2),使得概率在两个 class 之间是偶数。从数学上讲,它们是以下问题的解决方案:

b + w1*x1 + x2*x2 + x11*x1^2 + w12*x1*x2 + w22x2^2 = 0 

如果我们固定x1,那么这是一个x2的二次方程,我们可以解析求解。以下函数完成这项工作。

def boundary(x1):
  # returns x2 on the boundary for a given x1
  # we solve square equation
  # a x^2 + b x + c = 0 
  # --> x = (-b +- sqrt(b^2 - 4ac)) / 2a
  a_ = w22
  b_ = w2 + w12 * x1
  c_ = b + w1*x1 + w11*x1**2
  tmp = b_**2 - 4*a_*c_
  if tmp < 0:
    return None
  ans = [(-b_ + tmp**0.5) / (2*a_), (-b_ - tmp**0.5) / (2*a_)]
  ans.sort()  # smaller first
  return ans
# compute the boundaries
xs = np.linspace(X[:,0].min(), X[:,0].max(), num=100)
ys_1 = []
ys_2 = []
for x1 in xs:
  tmp = boundary(x1)
  if tmp is None:
    ys_1.append(None)
    ys_2.append(None)
  else:
    ys_1.append(tmp[0])  # smaller boundary
    ys_2.append(tmp[1])  # larger boundary

现在我们将边界作为数据,我们可以轻松地将它们可视化。

flg = (Y > 0)
plt.scatter(X[flg,0], X[flg,1], alpha=0.3, marker="o")
plt.scatter(X[~flg,0], X[~flg,1], alpha=0.3, marker="x")

plt.plot(xs, ys_1, c="green")
plt.plot(xs, ys_2, c="gray")

# if ys contains None, need to skip them
plt.fill_between(xs, ys_1, ys_2, color='tab:blue', alpha=0.2)
plt.fill_between(xs, min(ys_1), ys_1, color='tab:orange', alpha=0.2)
plt.fill_between(xs, ys_2, max(ys_2), color='tab:orange', alpha=0.2)

请注意,由于模型是二次模型,因此可以显式计算边界。更通用、更复杂的 classifier 需要不同的方法。

一种更简单、普遍适用的方法是创建包含各种变量组合的虚拟数据,让 classifier 进行预测,并使用预测的 class.[=24 给出的颜色进行绘图=]

xs = np.linspace(X[:,0].min(), X[:,0].max(), num=100)
ys = np.linspace(X[:,1].min(), X[:,1].max(), num=100)

newX = []
for x1 in xs:
  for x2 in ys:
    newX.append((x1, x2))
newX = np.array(newX)
p = clf.predict(poly.transform(newX))

flg = (Y > 0)
plt.scatter(X[flg,0], X[flg,1], alpha=0.3, marker="o")
plt.scatter(X[~flg,0], X[~flg,1], alpha=0.3, marker="x")

flg = (p > 0)
plt.scatter(newX[flg,0], newX[flg,1], alpha=0.02, c="tab:blue", marker="s", s=20)
plt.scatter(newX[~flg,0], newX[~flg,1], alpha=0.02, c="tab:orange", marker="s", s=20)