带有 \(\ell_1\) 正则化的逻辑回归

在这个例子中,我们使用CVXPY来训练一个带有 \(\ell_1\) 正则化的逻辑回归分类器。我们给定数据 \((x_i,y_i)\)\(i=1,\ldots, m\)。其中 \(x_i \in {\bf R}^n\) 是特征向量,而 \(y_i \in \{0, 1\}\) 是相关的布尔类别。

我们的目标是构建一个线性分类器 \(\hat y = \mathbb{1}[\beta^T x > 0]\),当 \(\beta^T x\) 是正时为 \(1\),否则为 \(0\)。我们线性地对给定数据建模,假设

\[\log \frac{\mathrm{Pr} (Y=1 \mid X = x)}{\mathrm{Pr} (Y=0 \mid X = x)} = \beta^T x.\]

这意味着

\[\mathrm{Pr} (Y=1 \mid X = x) = \frac{\exp(\beta^T x)}{1 + \exp(\beta^T x)}, \quad \mathrm{Pr} (Y=0 \mid X = x) = \frac{1}{1 + \exp(\beta^T x)}.\]

我们通过最大化数据的对数似然函数加上正则化项 \(\lambda \|\beta\|_1\) 来拟合 \(\beta\),其中 \(\lambda > 0\)

import cvxpy as cp
import numpy as np
import matplotlib.pyplot as plt

以下代码在随机选择 \(x_i\) 并提供稀疏的 \(\beta_{\mathrm{true}} \in {\bf R}^n\) 后,生成具有 \(n=50\) 个特征的数据。然后,我们设置 \(y_i = \mathbb{1}[\beta_{\mathrm{true}}^T x_i + z_i > 0]\), 其中 \(z_i\) 是 i.i.d. 的正态随机变量。我们将数据分为包含 \(m=50\) 个示例的训练集和测试集。

np.random.seed(1)
n = 50
m = 50
def sigmoid(z):
  return 1/(1 + np.exp(-z))

beta_true = np.array([1, 0.5, -0.5] + [0]*(n - 3))
X = (np.random.random((m, n)) - 0.5)*10
Y = np.round(sigmoid(X @ beta_true + np.random.randn(m)*0.5))

X_test = (np.random.random((2*m, n)) - 0.5)*10
Y_test = np.round(sigmoid(X_test @ beta_true + np.random.randn(2*m)*0.5))

我们使用 CVXPY 来表示优化问题。

beta = cp.Variable(n)
lambd = cp.Parameter(nonneg=True)
log_likelihood = cp.sum(
    cp.multiply(Y, X @ beta) - cp.logistic(X @ beta)
)
problem = cp.Problem(cp.Maximize(log_likelihood/m - lambd * cp.norm(beta, 1)))

我们解决优化问题,范围为 \(\lambda\) ,计算一个权衡曲线。然后绘制权衡曲线上的训练误差和测试误差。一个合理的 \(\lambda\) 的选择是使测试误差最小化的值。

def error(scores, labels):
  scores[scores > 0] = 1
  scores[scores <= 0] = 0
  return np.sum(np.abs(scores - labels)) / float(np.size(labels))
trials = 100
train_error = np.zeros(trials)
test_error = np.zeros(trials)
lambda_vals = np.logspace(-2, 0, trials)
beta_vals = []
for i in range(trials):
    lambd.value = lambda_vals[i]
    problem.solve()
    train_error[i] = error( (X @ beta).value, Y)
    test_error[i] = error( (X_test @ beta).value, Y_test)
    beta_vals.append(beta.value)
%matplotlib inline
%config InlineBackend.figure_format = "svg"

plt.plot(lambda_vals, train_error, label="Train error")
plt.plot(lambda_vals, test_error, label="Test error")
plt.xscale("log")
plt.legend(loc="upper left")
plt.xlabel(r"$\lambda$", fontsize=16)
plt.show()
../../_images/logistic_regression_9_0.svg

我们还绘制了正则化路径,即 \(\beta_i\)\(\lambda\) 的关系。注意,相对于其他特征来说,一些特征在较大的 \(\lambda\) 下保持非零较长的时间,这表明这些特征是最重要的。

for i in range(n):
    plt.plot(lambda_vals, [wi for wi in beta_vals])
plt.xlabel(r"$\lambda$", fontsize=16)
plt.xscale("log")
../../_images/logistic_regression_11_0.svg

我们绘制了真实的 \(\beta\) 与重建的 \(\beta\),选择使其在测试集上的误差最小化。非零系数被准确地重建出来。重建的 \(\beta\) 中有几个值是非零的,但实际上应该是零。

idx = np.argmin(test_error)
plt.plot(beta_true, label=r"True $\beta$")
plt.plot(beta_vals[idx], label=r"Reconstructed $\beta$")
plt.xlabel(r"$i$", fontsize=16)
plt.ylabel(r"$\beta_i$", fontsize=16)
plt.legend(loc="upper right")
<matplotlib.legend.Legend at 0x108adedd8>
../../_images/logistic_regression_13_1.svg