 验证)
Proximal Gradient Method 实战：LASSO回归 Python 实现与收敛速度验证
在机器学习领域，正则化技术是防止模型过拟合的关键手段。其中LASSO（Least Absolute Shrinkage and Selection Operator）回归因其能够产生稀疏解的特性，在特征选择和高维数据分析中具有独特优势。本文将深入探讨如何利用近似点梯度法（Proximal Gradient Method）高效解决LASSO回归问题，并通过Python实现验证其O(1/k)的收敛速度。
1. LASSO回归与优化挑战
LASSO回归可以表示为以下优化问题：
import numpy as np def lasso_objective(X, y, beta, alpha): """计算LASSO目标函数值""" residual = y - X @ beta loss = 0.5 * np.sum(residual**2) penalty = alpha * np.sum(np.abs(beta)) return loss + penalty该问题的核心难点在于ℓ₁正则项的非光滑特性，这使得传统梯度下降法无法直接应用。我们需要处理以下两个关键部分：
- 光滑部分：最小二乘损失函数f(β)=½‖Xβ-y‖²₂
- 非光滑部分：ℓ₁正则项h(β)=α‖β‖₁
为什么需要特殊优化方法？
| 优化方法 | 适用条件 | 处理LASSO的能力 |
|---|---|---|
| 梯度下降 | 光滑函数 | 无法直接应用 |
| 次梯度法 | 非光滑函数 | 收敛速度慢(O(1/√k)) |
| 近似点梯度法 | 复合函数 | 高效处理(O(1/k)) |
2. 近似点梯度法原理剖析
近似点梯度法的核心思想是将复合优化问题分解处理：
- 对光滑部分执行梯度下降
- 对非光滑部分应用邻近算子
算法迭代格式为： β⁽ᵏ⁺¹⁾ = prox_{tₖh}(β⁽ᵏ⁾ - tₖ∇f(β⁽ᵏ⁾))
关键组件实现：
2.1 软阈值算子实现
ℓ₁正则项的邻近算子就是著名的软阈值函数：
def soft_threshold(x, threshold): """软阈值算子实现""" return np.sign(x) * np.maximum(np.abs(x) - threshold, 0)这个函数的数学表达式为： [prox_h(z)]_i = sign(z_i) · max{|z_i| - α, 0}
2.2 梯度计算
对于最小二乘损失，梯度计算为：
def compute_gradient(X, y, beta): """计算最小二乘梯度""" return X.T @ (X @ beta - y)3. 完整算法Python实现
下面我们实现带有固定步长和回溯线搜索两种策略的近似点梯度法：
class ProximalGradientLASSO: def __init__(self, alpha=1.0, max_iter=1000, tol=1e-6): self.alpha = alpha # 正则化系数 self.max_iter = max_iter # 最大迭代次数 self.tol = tol # 收敛阈值 def fit(self, X, y, step_size=None, line_search=True): """拟合LASSO模型""" n_samples, n_features = X.shape self.beta = np.zeros(n_features) # 初始化系数 self.loss_history = [] for k in range(self.max_iter): grad = compute_gradient(X, y, self.beta) # 步长选择策略 if line_search: step_size = self.backtracking_line_search(X, y, grad) elif step_size is None: step_size = 1.0 / np.linalg.norm(X, ord=2)**2 # Lipschitz常数倒数 # 梯度步 gradient_step = self.beta - step_size * grad # 邻近算子步 self.beta_new = soft_threshold(gradient_step, step_size * self.alpha) # 记录损失 current_loss = lasso_objective(X, y, self.beta_new, self.alpha) self.loss_history.append(current_loss) # 检查收敛 if np.linalg.norm(self.beta_new - self.beta) < self.tol: break self.beta = self.beta_new.copy() return self def backtracking_line_search(self, X, y, grad, beta=0.8): """回溯线搜索确定步长""" step_size = 1.0 current_loss = lasso_objective(X, y, self.beta, self.alpha) while True: gradient_step = self.beta - step_size * grad beta_new = soft_threshold(gradient_step, step_size * self.alpha) new_loss = lasso_objective(X, y, beta_new, self.alpha) # 检查充分下降条件 if new_loss <= current_loss + grad.dot(beta_new - self.beta) + \ (1/(2*step_size)) * np.linalg.norm(beta_new - self.beta)**2: break step_size *= beta return step_size4. 收敛速度实验验证
为了验证算法的O(1/k)收敛速度，我们设计以下实验：
import matplotlib.pyplot as plt from sklearn.datasets import make_regression # 生成实验数据 X, y = make_regression(n_samples=100, n_features=50, noise=0.1, random_state=42) X = (X - X.mean(axis=0)) / X.std(axis=0) # 标准化 # 运行算法 pg_lasso = ProximalGradientLASSO(alpha=0.1, max_iter=500) pg_lasso.fit(X, y, line_search=False) # 固定步长 pg_lasso_ls = ProximalGradientLASSO(alpha=0.1, max_iter=500) pg_lasso_ls.fit(X, y, line_search=True) # 线搜索 # 绘制收敛曲线 plt.figure(figsize=(10, 6)) plt.semilogy(pg_lasso.loss_history - pg_lasso.loss_history[-1], label='固定步长', linewidth=2) plt.semilogy(pg_lasso_ls.loss_history - pg_lasso_ls.loss_history[-1], label='回溯线搜索', linewidth=2) plt.plot(1.0 / (np.arange(500)+1), '--', label='O(1/k)参考线') plt.xlabel('迭代次数', fontsize=12) plt.ylabel('目标函数值差距(f(x)-f*)', fontsize=12) plt.title('近似点梯度法收敛速度验证', fontsize=14) plt.legend(fontsize=12) plt.grid(True, which='both', linestyle='--') plt.show()实验结果分析：
- 两种步长策略都展现出O(1/k)量级的收敛速度
- 回溯线搜索通常能获得更快的实际收敛速度
- 固定步长需要估计合适的Lipschitz常数
- 实际应用中，线搜索虽然单次迭代成本略高，但总体更高效
5. 工程实践中的关键技巧
在实际应用中，我们还需要注意以下几个关键点：
5.1 特征缩放的重要性
LASSO回归对特征尺度敏感，因此必须进行特征标准化：
from sklearn.preprocessing import StandardScaler scaler = StandardScaler() X_scaled = scaler.fit_transform(X)5.2 正则化路径计算
通过 warm start 技术高效计算不同α值的解路径：
alphas = np.logspace(-3, 1, 50) beta_path = [] beta_init = np.zeros(X.shape[1]) for alpha in alphas: pg_lasso = ProximalGradientLASSO(alpha=alpha, max_iter=1000) pg_lasso.fit(X, y, line_search=True) beta_path.append(pg_lasso.beta) beta_init = pg_lasso.beta.copy() # warm start5.3 与其他算法的对比
下表比较了几种常见LASSO求解器的特性：
| 算法 | 收敛速度 | 内存需求 | 适合问题规模 | 稀疏性利用 |
|---|---|---|---|---|
| 坐标下降 | 线性 | 低 | 大 | 优秀 |
| 近似点梯度 | O(1/k) | 中 | 中到大 | 良好 |
| ADMM | O(1/k) | 中高 | 中 | 良好 |
| L-BFGS | 超线性 | 高 | 小到中 | 一般 |
提示：当特征维度非常高（>10⁵）时，坐标下降法通常是更好的选择，因为它能更好地利用稀疏性。
6. 扩展与变体算法
6.1 加速近似点梯度法（FISTA）
通过引入动量项，可以获得O(1/k²)的收敛速度：
def fista(X, y, alpha, max_iter=1000): beta = np.zeros(X.shape[1]) beta_prev = beta.copy() t = 1 t_prev = 1 loss_history = [] for k in range(max_iter): # 计算梯度 grad = compute_gradient(X, y, beta) # 计算步长 step_size = 1.0 / np.linalg.norm(X, ord=2)**2 # 计算动量项 y_k = beta + ((t_prev - 1)/t) * (beta - beta_prev) # 梯度步+邻近算子 gradient_step = y_k - step_size * grad beta_new = soft_threshold(gradient_step, step_size * alpha) # 更新参数 t_prev = t t = (1 + np.sqrt(1 + 4 * t_prev**2)) / 2 beta_prev = beta.copy() beta = beta_new.copy() # 记录损失 loss_history.append(lasso_objective(X, y, beta, alpha)) return beta, loss_history6.2 随机近似点梯度法
对于大规模数据，可以使用随机版本提高计算效率：
def stochastic_proximal_gradient(X, y, alpha, max_iter=1000, batch_size=32): n_samples = X.shape[0] beta = np.zeros(X.shape[1]) loss_history = [] for k in range(max_iter): # 随机采样小批量 idx = np.random.choice(n_samples, batch_size, replace=False) X_batch, y_batch = X[idx], y[idx] # 计算随机梯度 grad = compute_gradient(X_batch, y_batch, beta) # 自适应步长 step_size = 1.0 / (k + 1) # 更新参数 gradient_step = beta - step_size * grad beta = soft_threshold(gradient_step, step_size * alpha) # 记录完整损失（非必需） if k % 10 == 0: loss_history.append(lasso_objective(X, y, beta, alpha)) return beta, loss_history7. 实际应用案例
让我们通过一个实际案例展示近似点梯度法在特征选择中的应用：
from sklearn.datasets import load_diabetes from sklearn.model_selection import train_test_split # 加载数据 data = load_diabetes() X, y = data.data, data.target X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42) # 训练模型 pg_lasso = ProximalGradientLASSO(alpha=0.5, max_iter=1000) pg_lasso.fit(X_train, y_train, line_search=True) # 评估结果 train_score = pg_lasso.score(X_train, y_train) # 需实现score方法 test_score = pg_lasso.score(X_test, y_test) selected_features = np.where(np.abs(pg_lasso.beta) > 1e-4)[0] print(f"选中的特征数量: {len(selected_features)}/{X.shape[1]}") print(f"训练集R²: {train_score:.3f}, 测试集R²: {test_score:.3f}") print("重要特征:", selected_features)在这个糖尿病数据集上的典型输出可能如下：
选中的特征数量: 4/10 训练集R²: 0.527, 测试集R²: 0.512 重要特征: [2 3 6 8]这展示了LASSO回归如何自动选择最相关的特征，同时保持模型的预测性能。
