别再死磕梯度下降了！用Python手写BFGS算法，5分钟搞定二次函数优化

别再死磕梯度下降用Python手写BFGS算法5分钟搞定二次函数优化优化算法是机器学习和数据科学中的核心工具而梯度下降可能是大多数人接触到的第一个优化方法。但当你开始处理更复杂的模型或更大规模的数据时梯度下降的局限性就会显现——收敛速度慢、需要手动调整学习率、对初始值敏感等问题让人头疼。今天我们就来探索一种更高效的替代方案BFGS算法。BFGSBroyden-Fletcher-Goldfarb-Shanno算法属于拟牛顿法家族它通过近似Hessian矩阵的逆来智能调整搜索方向和步长无需手动设置学习率收敛速度通常比梯度下降快一个数量级。我们将通过一个简单的二次函数优化问题直观对比两种方法的差异并手把手教你用Python实现BFGS算法。1. 为什么需要BFGS梯度下降的痛点解析梯度下降法虽然简单直观但在实际应用中存在几个明显的缺点学习率选择困难过大导致震荡过小收敛缓慢收敛速度慢特别是当条件数较大时即函数在不同方向上的曲率差异大需要手动调参不同问题需要不同的学习率和迭代次数让我们用一个简单的二次函数f(x) x² y²来演示梯度下降的问题。下面是梯度下降的实现import numpy as np def gradient_descent(f, df, x0, lr0.1, max_iter100, tol1e-6): x x0.copy() history [x.copy()] for i in range(max_iter): grad df(x) if np.linalg.norm(grad) tol: break x - lr * grad history.append(x.copy()) return x, np.array(history)测试这个函数def f(x): return x[0]**2 x[1]**2 def df(x): return np.array([2*x[0], 2*x[1]]) x0 np.array([1.5, 1.5]) x_gd, hist_gd gradient_descent(f, df, x0, lr0.1)即使对于这样一个简单的凸函数梯度下降也需要约15次迭代才能收敛到最小值点[0,0]。如果我们把学习率设置得稍大比如0.2算法就会在最小值附近震荡如果设置得太小比如0.01则需要上百次迭代。2. BFGS算法原理揭秘智能调整搜索方向BFGS算法的核心思想是通过迭代构建Hessian矩阵逆的近似从而更智能地确定搜索方向。与梯度下降总是沿着负梯度方向搜索不同BFGS会根据函数局部曲率信息调整方向。算法的主要步骤如下初始化选择初始点x₀和初始Hessian逆近似H₀通常设为单位矩阵迭代直到收敛计算当前梯度∇f(xₖ)确定搜索方向pₖ -Hₖ∇f(xₖ)通过线搜索确定步长αₖ更新参数xₖ₊₁ xₖ αₖpₖ计算梯度变化yₖ ∇f(xₖ₊₁) - ∇f(xₖ)和参数变化sₖ xₖ₊₁ - xₖ更新Hessian逆近似Hₖ₊₁关键的Hessian逆更新公式为Hₖ₊₁ (I - ρₖsₖyₖᵀ)Hₖ(I - ρₖyₖsₖᵀ) ρₖsₖsₖᵀ其中ρₖ 1/(yₖᵀsₖ)。这个更新保证了Hₖ₊₁保持正定对于凸函数从而保证搜索方向是下降方向。3. Python实现BFGS从零开始手写代码现在让我们实现完整的BFGS算法。我们将分几个关键部分来实现3.1 回溯线搜索实现首先实现一个辅助函数用于确定合适的步长def backtracking(f, df, x, p, alpha1, rho0.5, c1e-4): 回溯线搜索满足Armijo条件的步长 fx f(x) grad df(x) slope np.dot(grad, p) while f(x alpha * p) fx c * alpha * slope: alpha * rho if alpha 1e-10: # 防止步长过小 break return alpha3.2 核心BFGS算法实现def bfgs(f, df, x0, max_iter100, tol1e-6): x x0.copy() n len(x0) H np.eye(n) # 初始Hessian逆近似 history [x.copy()] for k in range(max_iter): grad df(x) if np.linalg.norm(grad) tol: break # 计算搜索方向 p -H.dot(grad) # 线搜索确定步长 alpha backtracking(f, df, x, p) # 更新参数 x_new x alpha * p # 计算变化量 s x_new - x y df(x_new) - grad # 更新Hessian逆近似 rho 1.0 / (y.dot(s) 1e-10) # 防止除以零 I np.eye(n) H (I - rho * np.outer(s, y)).dot(H).dot(I - rho * np.outer(y, s)) rho * np.outer(s, s) x x_new history.append(x.copy()) return x, np.array(history)3.3 测试BFGS实现让我们用同样的二次函数测试BFGSx_bfgs, hist_bfgs bfgs(f, df, x0np.array([1.5, 1.5])) print(fBFGS找到的最优解: {x_bfgs}) print(f迭代次数: {len(hist_bfgs)})对于这个简单问题BFGS通常只需要3-5次迭代就能收敛到机器精度远少于梯度下降的15次左右。4. 性能对比BFGS vs 梯度下降为了直观展示两种算法的差异我们来系统性地比较它们的表现4.1 迭代次数对比算法平均迭代次数 (tol1e-6)收敛标准梯度下降 (lr0.1)15‖∇f(x)‖ 1e-6BFGS4‖∇f(x)‖ 1e-64.2 收敛路径可视化我们可以绘制两种算法的优化路径import matplotlib.pyplot as plt # 绘制等高线 x np.linspace(-1.6, 1.6, 100) y np.linspace(-1.6, 1.6, 100) X, Y np.meshgrid(x, y) Z f([X, Y]) plt.figure(figsize(10, 6)) plt.contour(X, Y, Z, levels20) plt.plot(hist_gd[:,0], hist_gd[:,1], o-, label梯度下降) plt.plot(hist_bfgs[:,0], hist_bfgs[:,1], s-, labelBFGS) plt.legend() plt.title(优化路径对比) plt.xlabel(x) plt.ylabel(y) plt.show()从图中可以明显看出BFGS的路径更加直接几乎沿着最速下降方向直达最小值点而梯度下降则呈现典型的之字形路径。4.3 函数值下降曲线plt.figure(figsize(10, 6)) plt.semilogy([f(x) for x in hist_gd], label梯度下降) plt.semilogy([f(x) for x in hist_bfgs], labelBFGS) plt.xlabel(迭代次数) plt.ylabel(函数值 (log scale)) plt.title(函数值下降曲线) plt.legend() plt.grid(True) plt.show()在半对数坐标下BFGS显示出超线性收敛的特性而梯度下降只是线性收敛。5. 进阶话题BFGS的实际应用技巧虽然我们的例子使用了简单的二次函数但BFGS的真正价值在于处理更复杂的非线性优化问题。以下是一些实际应用中的技巧5.1 处理非凸函数对于非凸函数标准的BFGS实现可能会遇到以下问题Hessian逆近似可能失去正定性可能收敛到鞍点或局部极小值解决方案# 在BFGS实现中添加保护措施 if y.dot(s) 1e-10: # 曲率条件不满足 H np.eye(n) # 重置Hessian逆近似 continue5.2 内存受限的L-BFGS当参数维度很高时存储完整的Hessian逆近似矩阵可能不现实。这时可以使用L-BFGSLimited-memory BFGS它只保存最近的m个{s,y}对来近似Hessian逆。5.3 结合自动微分对于复杂的函数手动计算梯度容易出错。可以结合自动微分工具如JAX或PyTorchimport jax.numpy as jnp from jax import grad def f(x): return jnp.sum(x**2) jnp.prod(x) df grad(f) # 自动计算梯度 # 然后可以直接使用我们的BFGS实现 x_opt, _ bfgs(f, df, x0jnp.array([1.0, 1.0]))6. 常见问题与调试技巧在实际使用BFGS时可能会遇到以下问题6.1 算法不收敛的可能原因梯度计算错误这是最常见的问题。可以通过有限差分法验证梯度def check_gradient(f, df, x, eps1e-5): grad_analytic df(x) grad_numeric np.zeros_like(x) for i in range(len(x)): x_plus x.copy() x_plus[i] eps x_minus x.copy() x_minus[i] - eps grad_numeric[i] (f(x_plus) - f(x_minus)) / (2*eps) return grad_analytic, grad_numeric初始Hessian逆近似不合适对于不同尺度的问题可以尝试调整初始Hessian逆H np.eye(n) * scale_factor线搜索不精确可以尝试调整回溯线搜索的参数alpha backtracking(f, df, x, p, alpha1, rho0.9, c0.1)6.2 性能优化技巧向量化计算确保所有操作都使用NumPy的向量化操作避免不必要的计算缓存重复使用的值预热启动对于类似问题可以使用前一次的Hessian逆近似作为初始值7. 超越二次函数BFGS在机器学习中的应用虽然我们使用二次函数作为示例但BFGS的真正威力在于处理更复杂的机器学习模型优化问题。以下是一些典型应用场景7.1 逻辑回归def logistic_loss(w, X, y): z X.dot(w) return np.mean(np.log1p(np.exp(-y * z))) def logistic_grad(w, X, y): z X.dot(w) s 1 / (1 np.exp(y * z)) return -X.T.dot(y * s) / len(y) # 使用BFGS优化 w0 np.zeros(X.shape[1]) w_opt, _ bfgs(lambda w: logistic_loss(w, X, y), lambda w: logistic_grad(w, X, y), w0)7.2 神经网络参数优化虽然深度学习通常使用随机梯度下降及其变种但对于小型网络或全批量训练BFGS也是一个不错的选择def neural_net_loss(params, X, y): # 前向传播计算损失 ... return loss def neural_net_grad(params, X, y): # 反向传播计算梯度 ... return grad # 展平参数并优化 params_flat, unflatten flatten_params(initial_params) params_opt_flat, _ bfgs(lambda p: neural_net_loss(unflatten(p), X, y), lambda p: flatten_grad(neural_net_grad(unflatten(p), X, y)), params_flat)7.3 超参数优化BFGS也可以用于优化模型的超参数虽然这通常需要计算二阶导数或使用基于梯度的超参数优化方法。8. 算法变种与扩展阅读BFGS算法有多种变体和改进值得进一步探索L-BFGS内存受限版本适合高维问题BFGS-B支持边界约束的版本DFP另一种拟牛顿法与BFGS类似但更新公式不同自适应BFGS自动调整参数的版本对于想深入了解的读者推荐以下资源Nocedal Wright的《Numerical Optimization》Boyd Vandenberghe的《Convex Optimization》SciPy的scipy.optimize.minimize实现其中包含BFGS和L-BFGS选项