《非线性优化算法

您所在的位置:网站首页 线性优化方法的特点是 《非线性优化算法

《非线性优化算法

2024-07-10 23:50| 来源: 网络整理| 查看: 265

文章目录 概述对比推导应用曲线拟合 参考文献

概述 寻找参数向量,使得函数值最小的非线性优化算法Levenberg-Marquardt, LM算法,列文伯格-马夸尔特法LM算法成为**Gauss-Newton算法与最速下降法(梯度下降法,GD)**的结合,μ很小时,类似于高斯牛顿法;μ很大时,类似于LM法阻尼因子μ对于过参数化问题不敏感,能有效处理冗余参数问题,使代价函数陷入局部极小值的机会大大减小 对比 几种优化方法对比(转载) 方法介绍迭代公式最速下降法负梯度方向,收敛速度慢 x k + 1 = x k − α g k x_{k+1}=x_k-\alpha g_k xk+1​=xk​−αgk​Newton 法保留泰勒级数一阶和二阶项,二次收敛速度,但每步都计算Hessian矩阵,复杂 x k + 1 = x k − H k − 1 g k x_{k+1}=x_k-H^{-1}_{k}g_k xk+1​=xk​−Hk−1​gk​高斯牛顿法法目标函数的Jacobian 矩阵近似H矩阵,提高算法效率,但H矩阵不满秩则无法迭代 x k + 1 = x k − ( J T J ) − 1 J T r = x k − ( J T J ) − 1 g k x_{k+1}=x_k-(J^TJ)^{-1}J^Tr=x_k-(J^TJ)^{-1}g_k xk+1​=xk​−(JTJ)−1JTr=xk​−(JTJ)−1gk​LM法信赖域算法,解决H矩阵不满秩或非正定 x k + 1 = x k − ( J T J + μ I ) − 1 g k x_{k+1}=x_k-(J^TJ+μI)^{-1}g_k xk+1​=xk​−(JTJ+μI)−1gk​ 推导 目标函数优化,转化为最小化误差f(x)平方和

m i n F ( x 0 + △ x ) = ∣ ∣ f ( x 0 + △ x ∣ ∣ 2 minF(x_0+\triangle x) = ||f(x_0+\triangle x||^2 minF(x0​+△x)=∣∣f(x0​+△x∣∣2

一阶泰勒展开

m i n F ( x 0 + △ x ) = ∣ ∣ f ( x 0 ) + J T △ x ∣ ∣ 2 minF(x_0+\triangle x) = ||f(x_0)+J^T\triangle x||^2 minF(x0​+△x)=∣∣f(x0​)+JT△x∣∣2

雅可比矩阵J

J = [ ∂ y 1 x 1 . . . ∂ y 1 x n ⋯ ∂ y m x 1 . . . ∂ y m x n ] J=\left[\begin{matrix} \frac{∂y_1}{x_1} ... \frac{∂y_1}{x_n} \\ \cdots \\ \frac{∂y_m}{x_1} ... \frac{∂y_m}{x_n} \\ \end{matrix} \right] J=⎣⎡​x1​∂y1​​...xn​∂y1​​⋯x1​∂ym​​...xn​∂ym​​​⎦⎤​

求极值,函数求导,导数为0

∂ ∣ ∣ f ( x 0 ) + J T △ x ∣ ∣ 2 ∂ △ x = 0 \frac{∂||f(x_0)+J^T\triangle x||^2}{∂\triangle x} = 0 ∂△x∂∣∣f(x0​)+JT△x∣∣2​=0

求解得到

J J T △ x = − J f ( x 0 ) JJ^T\triangle x = -Jf(x_0) JJT△x=−Jf(x0​)

迭代公式为

△ x = − ( J T J ) − 1 J T f ( x 0 ) \triangle x = -(J^TJ)^{-1}J^Tf(x_0) △x=−(JTJ)−1JTf(x0​)

J^TJ不一定可逆,需要引入单位阵

H ≈ J T J + μ I H≈J^TJ+μI H≈JTJ+μI

更新迭代公式为

x k + 1 = x k − ( J T J + μ I ) − 1 g k x_{k+1}=x_k-(J^TJ+μI)^{-1}g_k xk+1​=xk​−(JTJ+μI)−1gk​

应用 实现可以参照论文里面的伪代码(完全一致)

在这里插入图片描述

曲线拟合 一组数据 # y = k1 * exp(k2 / (x + k3)) x = [1, 0.9, 0.8, 0.7, 0.6, 0.5, 0.4, 0.3, 0.2, 0.1, 0, 0.01, 0.02, 0.03, 0.04, 0.05, 0.08, 0.15] y = [2.98, 3.06, 3.17, 3.39, 3.71, 4.17, 4.98, 6.41, 9.09, 15.73, 57.38, 49.78, 42.42, 35.34, 29.87, 24.94, 18.71, 11.49] 曲线拟合 def func(x, para): # y = k1 * exp(k2 / (x + k3)) return para[0][0] * np.exp(para[1][0] / (x + para[2][0])) def jacobi(x, para): J = np.zeros((len(x), para.shape[0])) for i, item in enumerate(x): J[i, 0] = np.exp(para[1][0] / (item + para[2][0])) J[i, 1] = para[1][0] * np.exp(para[1][0] / (item + para[2][0])) / (item + para[2][0]) J[i, 2] = para[1][0] * np.exp(para[1][0] / (item + para[2][0])) * (-para[1][0] / (item + para[2][0])**2) return np.matrix(J) def fit_curve(x, y, iter_num = 10000, eps=1e-8): num_paras = 3 para_past = np.ones((num_paras,1)) # 参数初始化 y_gj = func(x, para_past) J = jacobi(x, para_past) # jacobi r_past = np.matrix(y - y_gj).T # 残差矩阵 print(J.shape, r_past.shape) g = J.T * r_past tao = 1e-3 # (1e-8, 1) u = tao * np.max(J.T * J) # 阻尼因子初始化值 v = 2 norm_inf = np.linalg.norm(J.T * r_past, ord = np.inf) stop = norm_inf < eps num = 0 while (not stop) and num 0: para_past = para_cur r_past = r_cur J = jacobi(x, para_past) g = J.T * r_past stop = (np.linalg.norm(g, ord=np.inf)


【本文地址】


今日新闻


推荐新闻


CopyRight 2018-2019 办公设备维修网 版权所有 豫ICP备15022753号-3