ToB企服应用市场:ToB评测及商务社交产业平台

标题: Python 实现 LM 算法(Levenberg-Marquardt) [打印本页]

作者: 来自云龙湖轮廓分明的月亮    时间: 3 天前
标题: Python 实现 LM 算法(Levenberg-Marquardt)
博客:Python 实现 LM 算法(Levenberg-Marquardt)

目录


1. 弁言

什么是 Levenberg-Marquardt (LM) 算法?

Levenberg-Marquardt(LM)算法是求解非线性最小二乘问题的著名方法。它结合了高斯牛顿法和梯度下降法的优点,广泛用于非线性曲线拟合、机器学习、最优参数估计等领域。
LM 算法的一个主要特点是可以大概在问题接近线性时,像高斯牛顿法一样快速收敛;而在问题严峻非线性时,它通过梯度下降法的机制制止了求解发散。它是一个调和局部曲率信息与梯度信息的肴杂算法。
LM 算法的应用场景


LM 算法的优点与局限性

优点:
局限性:

2. LM 算法的原理

LM 算法的基本思想

LM 算法是一种在高斯牛顿法和梯度下降法之间动态调整的算法。当误差函数的二次近似足够好时,算法运动类似高斯牛顿法;当近似不好时,算法则趋向梯度下降法。这通过一个阻尼因子 (\lambda) 来调节,该因子控制了每次迭代时的步长和方向。
LM 算法的数学推导

LM 算法的目标是最小化非线性最小二乘问题:
                                                                min                                  ⁡                                          θ                                                 ∑                                           i                                  =                                  1                                          n                                                             (                                               y                                     i                                              −                                  f                                  (                                               x                                     i                                              ,                                  θ                                  )                                  )                                          2                                            \min_{\theta} \sum_{i=1}^{n} \left( y_i - f(x_i, \theta) \right)^2                     θmin​i=1∑n​(yi​−f(xi​,θ))2
迭代的更新公式为:
                                         Δ                            θ                            =                            −                            (                                       J                               T                                      J                            +                            λ                            I                                       )                                           −                                  1                                                            J                               T                                      r                            (                            θ                            )                                  \Delta \theta = -(J^T J + \lambda I)^{-1} J^T r(\theta)                     Δθ=−(JTJ+λI)−1JTr(θ)
此中:

                                    λ                              \lambda                  λ 决定了算法的特性。当                                    λ                              \lambda                  λ 较大时,                                   λ                         I                              \lambda I                  λI 的影响增大,算法运动趋向于梯度下降法;当                                    λ                              \lambda                  λ 较小时,算法运动接近高斯牛顿法。
与高斯牛顿法和梯度下降法的比较



3. Python 实现 LM 算法

面向对象的设计思绪

为了实现 LM 算法,首先须要构建一个代表非线性模子的类,用于计算误差和雅可比矩阵;其次,设计一个 LM 算法的类,用于执行参数优化。通过这种面向对象的方式,我们可以清楚地封装模子和优化算法的功能。
代码实现

  1. import numpy as np
  2. class NonlinearModel:
  3.     """表示非线性模型的类,包含残差和Jacobian矩阵的计算。"""
  4.     def __init__(self, func, jacobian):
  5.         """
  6.         :param func: 非线性模型函数
  7.         :param jacobian: Jacobian矩阵的计算函数
  8.         """
  9.         self.func = func
  10.         self.jacobian = jacobian
  11.    
  12.     def residuals(self, x_data, y_data, theta):
  13.         """计算残差向量"""
  14.         return y_data - self.func(x_data, theta)
  15.    
  16.     def jacobian_matrix(self, x_data, theta):
  17.         """计算给定参数下的Jacobian矩阵"""
  18.         return self.jacobian(x_data, theta)
  19. class LevenbergMarquardt:
  20.     """Levenberg-Marquardt算法的实现类。"""
  21.     def __init__(self, model, tolerance=1e-6, max_iters=100, lambda_init=0.01):
  22.         """
  23.         :param model: 待拟合的非线性模型对象
  24.         :param tolerance: 收敛阈值
  25.         :param max_iters: 最大迭代次数
  26.         :param lambda_init: 初始阻尼因子lambda
  27.         """
  28.         self.model = model
  29.         self.tolerance = tolerance
  30.         self.max_iters = max_iters
  31.         self.lambda_init = lambda_init
  32.    
  33.     def fit(self, x_data, y_data, initial_theta):
  34.         """使用Levenberg-Marquardt算法拟合模型参数"""
  35.         theta = initial_theta
  36.         lambda_factor = self.lambda_init
  37.         
  38.         for i in range(self.max_iters):
  39.             residuals = self.model.residuals(x_data, y_data, theta)
  40.             jacobian = self.model.jacobian_matrix(x_data, theta)
  41.             
  42.             # 计算Hessian近似
  43.             H = jacobian.T @ jacobian
  44.             
  45.             # 更新公式中增加lambda项
  46.             delta_theta = np.linalg.inv(H + lambda_factor * np.eye(H.shape[0])) @ jacobian.T @ residuals
  47.             
  48.             # 更新参数
  49.             theta_new = theta + delta_theta
  50.             
  51.             # 计算新的残差
  52.             residuals_new = self.model.residuals(x_data, y_data, theta_new)
  53.             
  54.             # 判断是否收敛
  55.             if np.linalg.norm(delta_theta) < self.tolerance:
  56.                 print(f"迭代收敛,共迭代 {i+1} 次")
  57.                 return theta_new
  58.             
  59.             # 动态调整lambda
  60.             if np.linalg.norm(residuals_new) < np.linalg.norm(residuals):
  61.                 lambda_factor /= 10  # 减少lambda
  62.                 theta = theta_new
  63.             else:
  64.                 lambda_factor *= 10  # 增加lambda
  65.         
  66.         print("达到最大迭代次数,未能完全收敛。")
  67.         return theta
  68. # 使用示例
  69. if __name__ == "__main__":
  70.     # 定义非线性模型 y = a * exp(b * x)
  71.     def func(x, theta):
  72.         return theta[0] * np.exp(theta[1] * x)
  73.    
  74.     # 定义Jacobian矩阵
  75.     def jacobian(x, theta):
  76.         J = np.zeros((len(x), len(theta)))
  77.         J[:, 0] = np.exp(theta[1] * x)
  78.         J[:, 1] = theta[0] * x * np.exp(theta[1] * x)
  79.         return J
  80.    
  81.     # 创建模型和LM算法实例
  82.     model = NonlinearModel(func, jacobian)
  83.     lm_solver = LevenbergMarquardt(model)
  84.    
  85.     # 生成数据
  86.     x_data = np.linspace(0, 1,
  87. 10)
  88.     y_data = 2 * np.exp(3 * x_data) + np.random.normal(0, 0.1, size=x_data.shape)
  89.    
  90.     # 初始参数猜测
  91.     initial_theta = np.array([1, 1])
  92.    
  93.     # 执行拟合
  94.     optimal_theta = lm_solver.fit(x_data, y_data, initial_theta)
  95.     print("最优参数:", optimal_theta)
复制代码

4. LM 算法应用实例:非线性曲线拟合

场景形貌

我们将使用 Levenberg-Marquardt 算法来拟合非线性模子                                    y                         =                         a                         ⋅                                   e                                       b                               ⋅                               x                                                 y = a \cdot e^{b \cdot x}                  y=a⋅eb⋅x。生成一组具有噪声的样本数据,并通过算法找到最优的                                    a                              a                  a 和                                    b                              b                  b 参数,使得模子曲线与数据点的误差最小。
结果分析与可视化

通过 matplotlib 可视化拟合结果:
  1. import matplotlib.pyplot as plt
  2. # 绘制拟合曲线与真实数据点
  3. y_fit = func(x_data, optimal_theta)
  4. plt.scatter(x_data, y_data, label="数据点")
  5. plt.plot(x_data, y_fit, label="拟合曲线", color='r')
  6. plt.xlabel("x")
  7. plt.ylabel("y")
  8. plt.legend()
  9. plt.show()
复制代码
结果显示,拟合的曲线与生成的数据点非常吻合,验证了 Levenberg-Marquardt 算法在非线性拟合问题中的有效性。

5. LM 算法的改进与扩展

LM 算法中的参数调节

LM 算法中的阻尼因子 (\lambda) 是一个至关重要的调节参数。初始值的选择以及动态调整策略会影响算法的收敛速度和结果。通常的做法是根据误差的变化动态增大或减小 (\lambda) 的值,从而在梯度下降和高斯牛顿法之间做出均衡。
LM 算法的改进与其他变种



6. 总结

在本篇博客中,我们详细介绍了 Levenberg-Marquardt 算法的基本思想、数学原理以及其 Python 的面向对象实现。通过非线性曲线拟合的示例,我们展示了该算法在实际应用中的有效性。LM 算法以其快速收敛和较好的鲁棒性,成为处置惩罚非线性最小二乘问题的重要工具。

免责声明:如果侵犯了您的权益,请联系站长,我们会及时删除侵权内容,谢谢合作!更多信息从访问主页:qidao123.com:ToB企服之家,中国第一个企服评测及商务社交产业平台。




欢迎光临 ToB企服应用市场:ToB评测及商务社交产业平台 (https://dis.qidao123.com/) Powered by Discuz! X3.4