第6章:逻辑斯谛回归与最大熵模型(Logistic Regression and Maximum Entropy Models)
一、导读
1.1 本章学习目标
第6章介绍逻辑斯谛回归——虽然名字里有"回归",但它实际上是一种分类算法。逻辑斯谛回归是最经典的线性分类方法之一,在工业界有着广泛的应用。同时,本章还介绍最大熵模型——一种基于信息论原理的通用分类框架。
本章的核心学习目标包括:
- 理解逻辑斯谛回归模型:sigmoid函数如何将线性输出映射为概率
- 掌握对数几率的概念:为什么 log(P/(1-P)) = w·x 是一个好的模型?
- 理解学习策略:为什么用对数似然损失而不是平方损失?
- 掌握优化算法:梯度下降法、牛顿法、拟牛顿法(BFGS)
- 理解最大熵原理:为什么在满足约束的条件下要使熵最大?
- 实现逻辑斯谛回归:从零开始编码,对比不同优化算法
1.2 为什么本章如此重要
逻辑斯谛回归是工业界最常用的分类算法之一。它的优势在于:简单、高效、可解释、概率输出。在广告推荐、信用评分、疾病预测等领域,逻辑斯谛回归往往是基线模型。
逻辑斯谛回归也是理解神经网络的起点。单个逻辑斯谛回归单元就是神经网络中最基本的"神经元"。理解逻辑斯谛回归,就是理解深度学习的基石。
1.3 核心问题引导
在阅读本章时,请带着以下问题思考:
逻辑斯谛回归的sigmoid函数有什么特殊性质?为什么选择它?
对数几率(log-odds)的直觉含义是什么?
为什么逻辑斯谛回归用对数似然损失而不是平方损失?
梯度下降法和牛顿法在优化逻辑斯谛回归时有什么区别?
最大熵模型和逻辑斯谛回归有什么关系?
逻辑斯谛回归如何处理多分类问题?
二、核心概念详解
6.1 逻辑斯谛回归模型
6.1.1 Sigmoid函数
逻辑斯谛回归的核心是sigmoid函数:
σ(z) = 1 / (1 + e^(-z))sigmoid函数将任意实数映射到 (0, 1) 区间,可以解释为概率。
import numpy as np
def sigmoid(z):
return 1 / (1 + np.exp(-np.clip(z, -500, 500)))
# sigmoid函数的性质
z = np.linspace(-10, 10, 100)
print(f"sigmoid(0) = {sigmoid(0):.4f}")
print(f"sigmoid(∞) → {sigmoid(100):.4f}")
print(f"sigmoid(-∞) → {sigmoid(-100):.4f}")
print(f"sigmoid'(0) = {sigmoid(0) * (1 - sigmoid(0)):.4f}")
print("\nsigmoid的关键性质:")
print("1. 输出范围 (0, 1),可以解释为概率")
print("2. 单调递增")
print("3. 导数 σ'(z) = σ(z)(1-σ(z)),计算方便")
print("4. z=0时输出0.5,是决策边界")6.1.2 模型定义
二分类逻辑斯谛回归模型:
P(Y=1|x) = 1 / (1 + e^(-(w·x + b)))
P(Y=0|x) = 1 - P(Y=1|x) = 1 / (1 + e^(w·x + b))决策规则:如果 P(Y=1|x) ≥ 0.5,则预测为正类。
6.1.3 对数几率
逻辑斯谛回归的一个优美性质是:
log(P(Y=1|x) / P(Y=0|x)) = w·x + b即对数几率(log-odds)是输入的线性函数。这意味着:
- w·x 每增加1,几率(odds)变为原来的 e 倍
- w_i 表示特征 x_i 每增加1个单位,对数几率的变化量
import numpy as np
def log_odds_interpretation():
"""对数几率的直观解释"""
w = np.array([0.5, -0.3, 1.0])
b = -0.5
# 基准点
x0 = np.array([0, 0, 0])
z0 = w @ x0 + b
p0 = sigmoid(z0)
print(f"基准点 x={x0}: P(Y=1) = {p0:.4f}")
print()
# 每个特征增加1个单位
for i in range(3):
x1 = x0.copy()
x1[i] = 1
z1 = w @ x1 + b
p1 = sigmoid(z1)
odds_ratio = (p1 / (1 - p1)) / (p0 / (1 - p0))
print(f"特征{i}增加1: P(Y=1) = {p1:.4f}, 几率比 = {odds_ratio:.4f} = e^{w[i]:.1f}")
print(f"\nw[0]=0.5: 特征0增加1,对数几率增加0.5")
print(f"w[1]=-0.3: 特征1增加1,对数几率减少0.3")
print(f"w[2]=1.0: 特征2增加1,对数几率增加1.0(影响最大)")
log_odds_interpretation()6.2 学习策略:对数似然损失
6.2.1 为什么不用平方损失
逻辑斯谛回归不使用平方损失,原因有两个:
非凸性:平方损失在逻辑斯谛回归中不是凸函数,存在局部最优
不匹配:平方损失假设误差是高斯分布的,但分类问题的输出是伯努利分布
6.2.2 对数似然损失
逻辑斯谛回归使用对数似然损失(Log Loss / Cross-Entropy Loss):
L(w, b) = -(1/N) * Σ_i [y_i * log(p_i) + (1-y_i) * log(1-p_i)]其中 p_i = P(Y=1|x_i) = sigmoid(w·x_i + b)。
这个损失函数是凸函数,有唯一的全局最优解。
import numpy as np
def log_loss(y_true, y_pred):
"""对数似然损失"""
eps = 1e-15
y_pred = np.clip(y_pred, eps, 1 - eps)
return -np.mean(y_true * np.log(y_pred) + (1 - y_true) * np.log(1 - y_pred))
def squared_loss(y_true, y_pred):
"""平方损失"""
return np.mean((y_true - y_pred) ** 2)
# 对比两种损失函数的梯度行为
y_true = np.array([1, 0, 1, 0])
y_pred_wrong = np.array([0.1, 0.9, 0.1, 0.9]) # 错误预测
y_pred_right = np.array([0.9, 0.1, 0.9, 0.1]) # 正确预测
print(f"{'':>15} {'错误预测':>12} {'正确预测':>12}")
print("-" * 42)
print(f"{'对数损失':>15} {log_loss(y_true, y_pred_wrong):>12.4f} {log_loss(y_true, y_pred_right):>12.4f}")
print(f"{'平方损失':>15} {squared_loss(y_true, y_pred_wrong):>12.4f} {squared_loss(y_true, y_pred_right):>12.4f}")6.3 学习算法
6.3.1 梯度下降法
逻辑斯谛回归的梯度:
∂L/∂w = (1/N) * X^T · (p - y)
∂L/∂b = (1/N) * Σ(p_i - y_i)其中 p = sigmoid(X·w + b)。
import numpy as np
class LogisticRegressionGD:
def __init__(self, lr=0.01, n_iter=1000):
self.lr = lr
self.n_iter = n_iter
self.loss_history = []
def fit(self, X, y):
self.w = np.zeros(X.shape[1])
self.b = 0.0
for _ in range(self.n_iter):
z = X @ self.w + self.b
p = sigmoid(z)
# 梯度
dw = X.T @ (p - y) / len(y)
db = np.mean(p - y)
# 更新
self.w -= self.lr * dw
self.b -= self.lr * db
# 记录损失
loss = log_loss(y, p)
self.loss_history.append(loss)
return self
def predict_proba(self, X):
return sigmoid(X @ self.w + self.b)
def predict(self, X):
return (self.predict_proba(X) >= 0.5).astype(int)
# 测试
from sklearn.datasets import make_classification
X, y = make_classification(n_samples=500, n_features=10, random_state=42)
X = (X - X.mean(axis=0)) / X.std(axis=0)
model = LogisticRegressionGD(lr=0.1, n_iter=500)
model.fit(X, y)
print(f"梯度下降法最终损失: {model.loss_history[-1]:.6f}")
print(f"训练准确率: {np.mean(model.predict(X) == y):.4f}")6.3.2 牛顿法
牛顿法利用二阶导数(Hessian矩阵)信息,收敛速度更快:
w_{t+1} = w_t - H^{-1} · g其中 g 是梯度,H 是 Hessian 矩阵。
import numpy as np
class LogisticRegressionNewton:
def __init__(self, n_iter=20):
self.n_iter = n_iter
self.loss_history = []
def fit(self, X, y):
n, p = X.shape
self.w = np.zeros(p)
self.b = 0.0
# 合并偏置
X_aug = np.column_stack([X, np.ones(n)])
theta = np.zeros(p + 1)
for _ in range(self.n_iter):
z = X_aug @ theta
pred = sigmoid(z)
# 梯度
grad = X_aug.T @ (pred - y) / n
# Hessian矩阵
W = np.diag(pred * (1 - pred))
H = X_aug.T @ W @ X_aug / n
# 牛顿更新
theta -= np.linalg.solve(H + 1e-8 * np.eye(p + 1), grad)
loss = log_loss(y, pred)
self.loss_history.append(loss)
self.w = theta[:-1]
self.b = theta[-1]
return self
def predict_proba(self, X):
return sigmoid(X @ self.w + self.b)
def predict(self, X):
return (self.predict_proba(X) >= 0.5).astype(int)
# 对比收敛速度
model_gd = LogisticRegressionGD(lr=0.1, n_iter=500)
model_gd.fit(X, y)
model_newton = LogisticRegressionNewton(n_iter=20)
model_newton.fit(X, y)
print(f"梯度下降法: {len(model_gd.loss_history)}次迭代, 最终损失={model_gd.loss_history[-1]:.6f}")
print(f"牛顿法: {len(model_newton.loss_history)}次迭代, 最终损失={model_newton.loss_history[-1]:.6f}")
print(f"牛顿法收敛快得多!但每次迭代计算Hessian矩阵的代价更高")6.3.3 拟牛顿法(BFGS)
BFGS是一种拟牛顿法,它通过近似Hessian矩阵的逆来避免直接计算Hessian:
import numpy as np
from scipy.optimize import minimize
class LogisticRegressionBFGS:
def fit(self, X, y):
n, p = X.shape
def objective(theta):
w, b = theta[:-1], theta[-1]
z = X @ w + b
p_pred = sigmoid(z)
loss = -np.mean(y * np.log(p_pred + 1e-10) + (1 - y) * np.log(1 - p_pred + 1e-10))
return loss
def gradient(theta):
w, b = theta[:-1], theta[-1]
z = X @ w + b
p_pred = sigmoid(z)
grad_w = X.T @ (p_pred - y) / n
grad_b = np.mean(p_pred - y)
return np.concatenate([grad_w, [grad_b]])
theta0 = np.zeros(p + 1)
result = minimize(objective, theta0, jac=gradient, method='BFGS')
self.w = result.x[:-1]
self.b = result.x[-1]
return self
def predict_proba(self, X):
return sigmoid(X @ self.w + self.b)
def predict(self, X):
return (self.predict_proba(X) >= 0.5).astype(int)
model_bfgs = LogisticRegressionBFGS()
model_bfgs.fit(X, y)
print(f"BFGS准确率: {np.mean(model_bfgs.predict(X) == y):.4f}")6.4 最大熵模型
6.4.1 最大熵原理
最大熵原理:在满足所有已知约束的条件下,选择熵最大的概率分布。
直觉:如果我们只知道部分信息,那么最"公平"的假设是让未知的部分尽可能均匀分布。
6.4.2 最大熵模型
最大熵模型是最大熵原理在分类问题中的应用。给定特征函数 f_i(x, y),最大熵模型的形式为:
P(y|x) = (1/Z(x)) * exp(Σ_i w_i * f_i(x, y))其中 Z(x) 是归一化因子。
当特征函数是线性函数时,最大熵模型等价于逻辑斯谛回归。
import numpy as np
def max_entropy_demo():
"""最大熵原理的直观演示"""
# 假设一个骰子,我们知道它的平均值是3.5
# 最大熵分布是什么?
# 约束:Σ p_i = 1, Σ i*p_i = 3.5
# 最大熵分布:均匀分布 p_i = 1/6
# 验证:均匀分布的均值
p_uniform = np.ones(6) / 6
mean_uniform = np.sum(np.arange(1, 7) * p_uniform)
entropy_uniform = -np.sum(p_uniform * np.log2(p_uniform))
print(f"均匀分布:")
print(f" 均值 = {mean_uniform:.2f}")
print(f" 熵 = {entropy_uniform:.4f} bits")
# 对比一个非均匀分布
p_biased = np.array([0.1, 0.1, 0.1, 0.3, 0.2, 0.2])
mean_biased = np.sum(np.arange(1, 7) * p_biased)
entropy_biased = -np.sum(p_biased * np.log2(p_biased))
print(f"\n非均匀分布:")
print(f" 均值 = {mean_biased:.2f}")
print(f" 熵 = {entropy_biased:.4f} bits")
print(f"\n均匀分布的熵更大,这就是最大熵原理")
max_entropy_demo()6.5 多分类扩展:Softmax回归
逻辑斯谛回归可以自然地扩展到多分类问题,使用softmax函数:
P(Y=k|x) = exp(w_k·x) / Σ_j exp(w_j·x)import numpy as np
def softmax(z):
"""softmax函数"""
z_max = np.max(z, axis=-1, keepdims=True)
exp_z = np.exp(z - z_max)
return exp_z / np.sum(exp_z, axis=-1, keepdims=True)
class SoftmaxRegression:
def __init__(self, lr=0.01, n_iter=500):
self.lr = lr
self.n_iter = n_iter
def fit(self, X, y):
n, p = X.shape
self.classes = np.unique(y)
K = len(self.classes)
self.W = np.zeros((p, K))
self.b = np.zeros(K)
# one-hot编码
Y_onehot = np.zeros((n, K))
for i, c in enumerate(self.classes):
Y_onehot[y == c, i] = 1
for _ in range(self.n_iter):
z = X @ self.W + self.b
P = softmax(z)
# 梯度
error = P - Y_onehot
self.W -= self.lr * X.T @ error / n
self.b -= self.lr * np.mean(error, axis=0)
return self
def predict_proba(self, X):
return softmax(X @ self.W + self.b)
def predict(self, X):
P = self.predict_proba(X)
return self.classes[np.argmax(P, axis=1)]
# 测试
from sklearn.datasets import load_iris
iris = load_iris()
X = (iris.data - iris.data.mean(axis=0)) / iris.data.std(axis=0)
y = iris.target
sm = SoftmaxRegression(lr=0.1, n_iter=1000)
sm.fit(X, y)
print(f"Softmax回归准确率: {np.mean(sm.predict(X) == y):.4f}")三、实战案例
3.1 正则化逻辑斯谛回归
import numpy as np
class RegularizedLogisticRegression:
def __init__(self, lr=0.01, n_iter=1000, reg_type='l2', C=1.0):
self.lr = lr
self.n_iter = n_iter
self.reg_type = reg_type
self.C = C
def fit(self, X, y):
n, p = X.shape
self.w = np.zeros(p)
self.b = 0.0
for _ in range(self.n_iter):
z = X @ self.w + self.b
pred = sigmoid(z)
# 梯度
dw = X.T @ (pred - y) / n
db = np.mean(pred - y)
# 正则化
if self.reg_type == 'l2':
dw += self.w / self.C
elif self.reg_type == 'l1':
dw += np.sign(self.w) / self.C
self.w -= self.lr * dw
self.b -= self.lr * db
return self
def predict(self, X):
return (sigmoid(X @ self.w + self.b) >= 0.5).astype(int)
# 对比不同正则化
np.random.seed(42)
X, y = make_classification(n_samples=200, n_features=20, n_informative=5, random_state=42)
X = (X - X.mean(axis=0)) / X.std(axis=0)
for reg in ['none', 'l1', 'l2']:
if reg == 'none':
model = LogisticRegressionGD(lr=0.1, n_iter=1000)
else:
model = RegularizedLogisticRegression(lr=0.1, n_iter=1000, reg_type=reg, C=0.1)
model.fit(X, y)
acc = np.mean(model.predict(X) == y)
n_nonzero = np.sum(np.abs(model.w) > 0.01)
print(f"{reg:>5}: 准确率={acc:.4f}, 非零系数={n_nonzero}")3.2 决策边界可视化
import numpy as np
def visualize_decision_boundary():
"""可视化逻辑斯谛回归的决策边界"""
np.random.seed(42)
# 生成二分类数据
X_pos = np.random.randn(50, 2) + np.array([1.5, 1.5])
X_neg = np.random.randn(50, 2) + np.array([-1, -1])
X = np.vstack([X_pos, X_neg])
y = np.array([1] * 50 + [0] * 50)
# 训练
model = LogisticRegressionGD(lr=0.1, n_iter=500)
model.fit(X, y)
print(f"w = [{model.w[0]:.4f}, {model.w[1]:.4f}]")
print(f"b = {model.b:.4f}")
print(f"决策边界: {model.w[0]:.4f}*x1 + {model.w[1]:.4f}*x2 + {model.b:.4f} = 0")
print(f"即 x2 = {-model.w[0]/model.w[1]:.4f}*x1 + {-model.b/model.w[1]:.4f}")
print(f"准确率: {np.mean(model.predict(X) == y):.4f}")
visualize_decision_boundary()四、本章小结
第6章介绍了逻辑斯谛回归与最大熵模型:
逻辑斯谛回归模型:通过sigmoid函数将线性输出映射为概率,对数几率是输入的线性函数。
学习策略:使用对数似然损失(交叉熵损失),它是凸函数,有唯一全局最优。
优化算法:梯度下降法简单但收敛慢;牛顿法利用二阶信息收敛快;BFGS拟牛顿法在实际中最常用。
最大熵模型:在满足约束条件下使熵最大,当特征为线性时等价于逻辑斯谛回归。
多分类扩展:通过softmax函数自然地扩展到多分类问题。
逻辑斯谛回归的核心启示是:简单的线性模型 + 非线性映射 = 强大的分类器。
五、思考与练习
思考题
sigmoid函数有什么特殊性质使得它特别适合分类问题?还有其他类似的函数吗?
为什么逻辑斯谛回归不用平方损失?从凸性和概率分布的角度解释。
牛顿法为什么比梯度下降法收敛快?它的缺点是什么?
最大熵原理的直觉是什么?为什么"最不确定"的分布反而最好?
逻辑斯谛回归在什么情况下会欠拟合?如何解决?
编程练习
从零实现逻辑斯谛回归(梯度下降法),在二分类数据上测试。
实现牛顿法优化逻辑斯谛回归,对比与梯度下降法的收敛速度。
实现Softmax回归,在鸢尾花数据集(三分类)上测试。