逻辑斯蒂回归(logistic regression)是最常见的机器学习算法之一,它可以用于预测事件发生的概率,例如基于给定的标记数据集预测一封来信是否为垃圾邮件或者一个肿瘤是否恶性。由于其简单性,逻辑斯蒂回归常被用作基准模型与其他更复杂的模型进行比较。该模型的名称中包含“逻辑”一词,因为它使用逻辑函数(sigmoid)将输入特征的线性组合转换为概率。逻辑斯蒂回归的输出是一个介于0和1之间的连续值,因此它的名称中还包含“回归”一词,但是通常作为二元分类器使用,通过选择一个阈值(通常为0.5),将概率大于阈值的输入分类为正类,低于阈值的分类为负类。
在这篇文章中,我们将深入讨论逻辑斯蒂回归模型,从零开始在Python中实现它,然后展示它在Scikit-Learn中的实现。在介绍该算法之前,我们需要考虑为什么不用线性回归进行分类。
线性回归模型试图找到输入特征和目标变量之间的线性关系,基本模型表示为:
其中,y
是目标变量(通常是连续的),
线性回归可以通过以下数学公式解释:
我们可以使用上述公式预测y
值,通过调整m和c的值,找到一条最适合数据的直线,这条线可以用于预测新输入数据的y
值,下图将有助于更好地理解这个概念:
通过观察图中的散点图和拟合直线,我们可以直观地理解线性回归模型是如何工作的。
x值由蓝色点表示(输入数据),现在我们可以使用输入数据计算斜率和y坐标,以确保我们的预测线(红色线)覆盖大部分位置,然后根据红色线预测给定x值的任何y值。
关于线性回归需要记住的一件事是,它只适用于连续数据。如果我们想在分类方法中使用线性回归,我们需要对算法进行调整。首先我们必须选择一个阈值,如果我们的预测值小于阈值,则属于类1;否则,属于类2。
现在你是不是在想,“哦,这很简单,只需创建带有阈值的线性回归就可以构建分类了“。
但是这样做有个很大的瓶颈:如何选择合适的阈值,以及这种方法是否对异常值非常敏感
我用吴恩达老师的例子来阐述线性回归分类的缺点吧。
如下图,假设我们有肿瘤大小和恶性程度的信息。因为这是一个分类问题,我们可以看到所有的值都将介于0和1之间。通过拟合最佳的回归线(蓝色线)并假设阈值为0.5,线性回归有很好的分类效果。
image-20240415204156075我们可以在 x 轴上选择一个点,使得该点左侧的所有值都被视为负类,而右侧的所有值则被视为正类。
image-20240415204347186如果数据中包含异常值,最佳回归拟合线会发生很大的变化,若我们仍然设置阈值0.5,则会出现分类错误,如下图:
image-20240415204726321上图中的绿色线对应阈值0.5,可知有3个恶性肿瘤错误的分类为良性肿瘤。若我们使用黄色线作为阈值,则能够准确的进行分类。
因此线性回归分类很难确定某一阈值,且阈值对异常值非常敏感。
逻辑回归是一种处理二元分类问题的概率分类器。对于给定的样本(x,y),它输出样本属于正类的概率p:
如果这个概率高于某个阈值(通常选择为0.5),则将样本分类为1,否则将其分类为0。
模型如何估计概率p?
逻辑回归的基本假设是样本属于正类的对数几率(log-odds,也称为logit)是其特征的线性组合。
对数几率(logit)的含义是样本属于正类的概率与它属于负类的概率之比的对数:
我们假设对数的底是自然对数e,你也可以使用其他对数的底。
上述logit公式的图表:
image-20240415210940833根据图表可得:logit函数将概率值的范围(0,1)映射到实数(-∞, +∞)之间。
如之前所述,逻辑回归假设对数几率是特征的线性组合,即方程为:
这里是模型的参数(权重),参数通常称为截距(或偏置)。
对于的点(对数几率等于0)定义两个类之间的分界超平面,方程为:
权重向量与是超平面正交,超平面上方的每个样本()被分类为正样本,超平面下方的每个样本()被分类为负样本。权重向量和超平面的关系如下图:
image-20240415214153120这使得逻辑回归成为一种线性分类器,因为它假设类别之间的边界是一个线性平面。其他线性分类器包括感知器和支持向量机(SVM)。
通过对对数几率方程两边取指数,我们可以找到 p 与参数 w 之间的直接关系:
这里是sigmoid函数,也称为逻辑函数:
sigmoid函数用于将对数几率()转换为概率,它具有特征性的"S"曲线。
image-20240415215722151正如上图所见,该函数将实数在(-∞, +∞)之间映射到概率值在(0,1)之间。
sigmoid函数具有一些有用的数学性质,这些性质将在后面有用:
下图总结了逻辑回归从输入到最终预测的计算过程:
image-20240415220001574我们的目标是计算模型参数,使预测概率尽可能的接近真实标签y,这里我们需要定义损失函数衡量预测概率和真实标签的接近程度。损失函数是可微的,这样就可以使用梯度下降等方法进行模型优化。
逻辑回归使用的损失函数称为对数损失(或逻辑损失)。其定义如下:
如何得到这个函数的呢?这个函数的推导基于最大似然原理。更具体地说,我们可以证明在假设标签服从伯努利分布(即二元随机变量以概率 p 取值 1,以概率 1 − p 取值 0)的情况下,对数损失是负对数似然。
数学上可以表示:
其中,是在给定模型预测p的情况下预测标签y的概率,即给定模型,求数据的似然性。
证明对数损失是负对数似然:
给定伯努利分布数据标签属于正类概率p:
负类概率:
结合上述两式,有:
简单验证上式,若y=1,有,,即P(y|p)=p;若y=0,有,,即。
因此给定伯努力分布的概率p,根据上式,数据的对数似然为:
根据上式和对数损失,可得:
即负的对数损失等于数据的对数似然,最小化对数损失等价于最大化数据的对数似然。
下图显示了当 y = 1 时的对数损失:
image-20240415224256618对数损失仅在完美预测的情况下等于 0(即 p = 1 且 y = 1,或 p = 0 且 y = 0),并且随着预测越糟糕(即当 y = 1且 p->0 或 y = 0 且 p->1)而趋近于无穷大。
损失函数计算整个数据集上的平均损失:
方程可以用向量的形式表示:
其中,是包含所有训练样本的标签向量,而是包含模型对所有训练样本的预测概率向量。
这个损失函数是凸的,即它具有唯一的全局最小值。然后由于对数函数引入了非线性,因此没有闭合形式的解来找到最优的模型参数。因此我们需要使用迭代优化方法(比如梯度下降),以找到最小值。
梯度下降法是一种迭代的方法,用于寻找函数的最小值。我们沿着梯度的反方向迈出小步,以接近最小值。
image-20240416125650850为了使用梯度下降法找到损失函数的最小值,我们需要计算其相对于每个权重的偏导数。损失函数对权重的偏导:
为了让大家有个更好的理解,推导上式:
因此,梯度向量可以用矢量化的方式表示:
梯度下降法更新参数:
这里是学习率,用于调控梯度更新的幅度()。
请注意,无论何时使用梯度下降,都必须确保数据集已经标准化(否则梯度下降可能在不同方向上采取不同大小的步骤,这将使其不稳定)。
可视化梯度下降过程:
445391_XRCJt-5yNXDfzrVbEbh4DA正则化就像在赛车道上加入护栏一样,它防止我们的模型偏离轨道和过拟合,使其保持稳定。
我们简要讨论一下最小化包含正则化的损失函数,以匹配训练数据集的模型参数。L1(Lasso)和L2(Ridge)是最常见的两种正则化类型。与简单地最大化上述损失函数不同,正则化对系数的大小施加了限制,以避免过拟合。L1和L2使用不同的方法来定义系数的上限,L1正则话通过将系数设置为0来进行特征选择并减少多重共线性,L2正则化对模型系数进行极大的惩罚,但不会将任何系数设置为0。还有一个参数,即约束权重的参数 λ,以确保系数不会受到过度惩罚导致欠拟合。
L1正则化和L2正则化分别以绝对值和平方表示,导致它们对模型参数的有不同的限制,以及 λ 如何影响正则化和原始拟合项的权重,这是一个很有趣的话题。我们在这里不会详细展开,但值得花时间和精力去学习。下面的步骤演示了如何将原始损失函数转换为正则化损失函数。
L1正则化:
L2正则化:
其中n表示样本个数,m表示特征个数
我们将从头开始使用Python实现逻辑回归模型,包括损失函数、梯度计算、使用梯度下降优化模型,评估模型并绘制最终的决策边界。
为了演示,我们将使用Iris数据集(BSD许可)。原始数据集包含150个Iris花的样本,属于三个品种之一(setosa,versicolor和virginica)。我们通过使用前两种类型的花(setosa和versicolor)将其转换为二元分类问题。此外,我们仅使用每个花的前两个特征(sepal width和sepal length)。
我们首先导入所需的库,并设置随机种子以获得可重现的结果:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
np.random.seed(0)
加载数据集:
from sklearn.datasets import load_iris
iris = load_iris()
X = iris.data[:, :2] # 获取前两个特征
y = iris.target
# 获取前两种花
X = X[(y == 0) | (y == 1)]
y = y[(y == 0) | (y == 1)]
绘制数据:
def plot_data(X, y):
sns.scatterplot(x=X[:, 0], y=X[:, 1], hue=iris.target_names[y], style=iris.target_names[y],
palette=['r','b'], markers=('s','o'), edgecolor='k')
plt.xlabel(iris.feature_names[0])
plt.ylabel(iris.feature_names[1])
plt.legend()
plot_data(X, y)
image-20240416215740845
如上所见,该数据集是线性可分的,因此逻辑回归应该能够找到两个类别之间的决策边界。
接下来,我们需要在特征矩阵 X 中添加一列全为1的列,以表示偏置项():
# 添加偏置列
n = X.shape[0]
X_with_bias = np.hstack((np.ones((n, 1)), X))
数据集划分为训练集和测试集:
from sklearn.model_selection import train_test_split
X_train, X_test, y_train, y_test = train_test_split(X_with_bias, y, random_state=0)
我们现在准备实现逻辑回归模型。我们首先定义一个辅助函数来计算sigmoid函数:
def sigmoid(z):
""" 计算z的逻辑函数,z可以是标量或向量 """
z = np.array(z)
return 1 / (1 + np.exp(-z))
接下来,我们实现损失函数,该函数返回给定数据集 (X, y) 上具有参数的逻辑回归模型的损失,以及其相对于 的梯度。
def cost_function(X, y, w):
# 计算损失
p = sigmoid(X @ w)
J = -(1/n) * (y @ np.log(p) + (1-y) @ np.log(1-p))
# 计算损失对参数的梯度
grad = (1/n) * X.T @ (p - y)
return J, grad
上述代码使用了损失函数和梯度函数的矢量化形式,若有不懂可参照之前的推导。
为了验证这个函数的正确性,我们计算模型在某个随机权重向量上的损失和梯度:
w = np.random.rand(X_train.shape[1])
cost, grad = cost_function(X_train, y_train, w)
print('w:', w)
print('Cost at w:', cost)
print('Gradient at w:', grad)
输出为:
w: [0.5488135 0.71518937 0.60276338]
Cost at w: 2.314505839067951
Gradient at w: [0.36855061 1.86634895 1.27264487]
我们现在实现梯度下降,以便在给定的训练集上找到最小化代价函数的最优参数,该算法最多会在训练集上运行 max_iter
次迭代(默认为 5000),除非当前迭代的损失相较于上次迭代小于 tol
(默认为 0.0001),这种情况下训练会立即停止。
def optimize_model(X, y, alpha=0.01, max_iter=5000, tol=0.0001):
""" 使用梯度下降优化模型.
X, y: 训练数据集
alpha: 学习率
max_iter: 最大迭代次数
tol: 迭代应大于损失减少的最小值
"""
w = np.random.rand(X.shape[1])
cost, grad = cost_function(X, y, w)
for i in range(max_iter):
w = w - alpha * grad
new_cost, grad = cost_function(X, y, w)
if new_cost > cost - tol:
print(f'Converged after {i} iterations')
return w, new_cost
cost = new_cost
print('Maximum number of iterations reached')
return w, cost
通常,在这个阶段你需要标准化你的数据集,因为梯度下降在特征尺度不同的情况下效果不佳。在我们特定的数据集中,标准化是不必要的,因为这两个特征的范围相似。
现在我们调用这个函数来优化我们的模型:
opt_w, cost = optimize_model(X_train, y_train)
print('opt_w:', opt_w)
print('Cost at opt_w:', cost)
根据打印结果可知,模型再1413次迭代获得模型最优参数:
nverged after 1413 iterations
opt_w: [ 0.28014029 0.80541854 -1.48367938]
Cost at opt_w: 0.28389717767222566
还有其他优化器可以用于优化,它们通常比梯度下降更快,例如共轭梯度法(CG)和截断牛顿法(TNC)。有关如何使用这些优化器的更多详细信息,请参阅 scipy.optimize.minimize
。
我们可以用模型的最优参数进行预测。
首先,我们编写一个函数,该函数接收一个新的样本矩阵X并返回它们属于正类的概率:
def predict_prob(X, w):
""" 返回样本矩阵X属于正样本的概率向量
X: 特征矩阵
w: 模型参数
"""
p = sigmoid(X @ w)
return p
该函数通过计算的sigmoid预测样本属于正类的概率。
例如,让我们找出一个位于(6, 2)的样本属于versicolor类的概率:
predict_prob([[1, 6, 2]], opt_w)
输出:
array([0.89522808])
这个样本有89.52%的概率是versicolor花。这是合理的,因为这个样本位于versicolor花的区域内,远离两个类别之间的边界。
另一方面,一个位于(5.5, 3)的样本属于versicolor类的概率是:
predict_prob([[1, 5.5, 3]], opt_w)
array([0.56436688])
这次概率要低得多(只有56.44%),因为这个样本靠近两个类别之间的边界。
我们编写另一个函数,该函数返回预测的类标签而不是概率:
def predict(X, w):
""" 预测概率大于等于0.5,属于0类;
预测概率小于0.5,属于1类
"""
p = sigmoid(X @ w)
y_pred = (p >= 0.5).astype(int)
return y_pred
该函数只在属于正类的概率至少为0.5时预测为1,否则为0。
让我们使用上面的样本测试这个函数:
predict([[1, 6, 2], [1, 5.5, 3]], opt_w)
输出:
array([1, 1])
不出所料,这两个样本都被分类为1。
我们编写一个函数来计算模型在给定数据集上的准确率:
def evaluate_model(X, y, w):
y_pred = predict(X, w)
accuracy = np.mean(y == y_pred)
return accuracy
该函数首先找到模型在给定数据集 X 上的预测标签,然后将其与真实标签 y 进行比较。准确率等于预测准确的个数与预测总数之比。
准确率预测准确个数预测总数我们使用上述函数来计算我们的模型在训练集和测试集上的准确率:
train_accuracy = evaluate_model(X_train, y_train, opt_w)
print(f'Train accuracy: {train_accuracy * 100:.3f}%')
输出:
Train accuracy: 98.667%
测试数据集评估:
test_accuracy = evaluate_model(X_test, y_test, opt_w)
print(f'Test accuracy: {test_accuracy * 100:.3f}%')
输出:
Test accuracy: 100.000%
正如预期的那样,由于数据集是线性可分的,得分非常高。
除了准确率之外,还有其他重要的指标用于评估分类模型,例如精确度、召回率和F1分数。
最后数据集是二维的,我们可以绘制类之间的边界线。为此,我们首先需要找到这条线的方程。
边界线定义为预测概率等于0.5,对应于边界线,即:
根据上式得到x1和x2的关系:
边界线的斜率是,截距是。
def plot_decision_boundary(X, y, w):
""" 绘制类的边界线 """
plot_data(X, y)
line_x = np.array(plt.gca().get_xlim())
line_y = -1 / w[2] * (w[1] * line_x + w[0])
plt.plot(line_x, line_y, c='k', ls='--')
plot_decision_boundary(X, y, opt_w)
image-20240416225049447
如上图,模型只有一个样本被错误分类。对模型进行更多次数的训练(大约20万次)将找到一条完全分隔两个类的线。使用固定的步长,梯度下降的最优收敛速度非常慢。通过使用自适应学习率(例如使用更激进的步长来补偿迅速消失的梯度),这可以得到改善。
本节介绍如何使用Scikit-Learn中的LogisticRegression类构建最优模型,这个类使用了比梯度下降更高效的优化器,并且还提供了额外的选项,如正则化和提前停止。
该类的重要超参数包括:
当使用LogisticRegression类时,您不需要手动向设计矩阵X中添加一列全1,因为Scikit-Learn会自动完成。因此,在构建模型之前,我们将原始数据(不包含额外的全1列)拆分为训练集和测试集:
from sklearn.model_selection import train_test_split
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=0)
使用LogisticRegression的默认设置创建一个实例,并将其拟合到训练集:
from sklearn.linear_model import LogisticRegression
clf = LogisticRegression()
clf.fit(X_train, y_train)
接下来,我们评估模型在训练集和测试集上的表现:
train_accuracy = clf.score(X_train, y_train)
print(f'Train accuracy: {train_accuracy * 100:.3f}%')
test_accuracy = clf.score(X_test, y_test)
print(f'Test accuracy: {test_accuracy * 100:.3f}%')
输出:
Train accuracy: 100.000%
Test accuracy: 100.000%
这次我们在训练集和测试集上都获得了完美的分数。我们还可以查询 n_iter_ 属性来检查收敛所需的迭代次数:
print(clf.n_iter_)
输出:
[13]
只需要15次迭代就能收敛!显然,LogisticRegression使用的优化器(默认为L-BFGS)比我们的梯度下降实现更高效。
我们可以像之前一样绘制模型的决策边界。然而,这次模型的最优系数存储在模型的两个不同属性中:
opt_w = np.insert(clf.coef_, 0, [clf.intercept_])
plot_decision_boundary(X, y, opt_w)
image-20240416230110532
正如预期的那样,LogisticRegression的边界线完美地分隔了两个类。
以下是与其他分类模型相比,逻辑回归的优缺点总结:
优点:
缺点:
欢迎扫码关注: