一、牛顿法
1、牛顿法在一维搜索中的应用
在一维搜索中我们所要解决的问题是如何找函数f(x)的最小值。
牛顿法的核心思想是用二次函数拟合函数f(x)的某一邻域区间,用二次函数的极小值点作为下一次的迭代点。通过多次迭代使得二次函数的极小值逼近函数f(x)的极小值
g ( x ) = f ( x ( k ) ) + f ′ ( x ( k ) ) ( x − x ( k ) ) + 1 2 f ′ ′ ( x ( k ) ) ( x − x ( k ) ) 2 g ( x ) ≈ f ( x ) , 求 f ( x ) 最小值 ≈ 求 g ( x ) 最小值 g ′ ( x ) = f ′ ( x ( k ) ) + f ′ ′ ( x ( k ) ) ( x − x ( k ) ) 令 g ′ ( x ) = 0 , x = x ( k ) − f ′ ( x ( k ) ) f ′ ′ ( x ( k ) ) 只有在 f ′ ′ ( x ) > 0 时成立, f ′ ( x ) = 0 只能保证该点为极值点, f ′ ′ ( x ) > 0 保证该点为极小值点 \begin{aligned} &g(x) = f(x^{(k)})+f'(x^{(k)})(x-x^{(k)})+\frac{1}{2}f''(x^{(k)})(x-x^{(k)})^2 \\ &g(x)\approx f(x),求f(x)最小值 \approx 求g(x)最小值 \\ &g'(x)=f'(x^{(k)})+f''(x^{(k)})(x-x^{(k)})\\ &令g'(x)=0,x = x^{(k)}-\frac{f'(x^{(k)})}{f''(x^{(k)})}\\ &只有在f''(x)>0时成立,f'(x)=0只能保证该点为极值点,f''(x)>0保证该点为极小值点 \end{aligned} g(x)=f(x(k))+f′(x(k))(x−x(k))+21f′′(x(k))(x−x(k))2g(x)≈f(x),求f(x)最小值≈求g(x)最小值g′(x)=f′(x(k))+f′′(x(k))(x−x(k))令g′(x)=0,x=x(k)−f′′(x(k))f′(x(k))只有在f′′(x)>0时成立,f′(x)=0只能保证该点为极值点,f′′(x)>0保证该点为极小值点
2、牛顿法在多维函数中的应用
多维的情况与一维类似,如果是二维函数,拟合的是一个二次曲面,用二次曲面的最低点作为下一次的迭代点。
g ( X ) = f ( X ( k ) ) + ( X − X ( k ) ) ∇ f ′ ( X ( k ) ) + 1 2 ( X − X ( k ) ) T ∇ 2 f ( X ( k ) ) ( X − X ( k ) ) g ( X ) ≈ f ( X ) , 求 f ( X ) 最小值 ≈ 求 g ( X ) 最小值 令 ∇ g ( X ) = ∇ f ( X ( k ) ) + ∇ f ( X ( k ) ) ( X − X ( k ) ) = 0 如果 ∇ 2 f ( X ) > 0 ( 正定矩阵 ) , X = X ( k ) − [ ∇ f ( x ( k ) ) ] − 1 ∇ f ( x ( k ) ) \begin{aligned} &g(X) = f(X^{(k)})+(X-X^{(k)})\nabla f'(X^{(k)})+\frac{1}{2}(X-X^{(k)})^T \nabla ^2f(X^{(k)})(X-X^{(k)}) \\ &g(X)\approx f(X),求f(X)最小值 \approx 求g(X)最小值 \\ &令\nabla g(X)=\nabla f(X^{(k)})+\nabla f(X^{(k)})(X-X^{(k)}) = 0\\ &如果\nabla^2f(X)>0(正定矩阵),X = X^{(k)}-[\nabla f(x^{(k)})]^{-1}\nabla f(x^{(k)})\\ \end{aligned} g(X)=f(X(k))+(X−X(k))∇f′(X(k))+21(X−X(k))T∇2f(X(k))(X−X(k))g(X)≈f(X),求f(X)最小值≈求g(X)最小值令∇g(X)=∇f(X(k))+∇f(X(k))(X−X(k))=0如果∇2f(X)>0(正定矩阵),X=X(k)−[∇f(x(k))]−1∇f(x(k))
3、Levenberg-Marquardt修正
上述方法只有在 H e s s Hess Hess矩阵正定是成立,如果 H e s s Hess Hess矩阵不是正定的要怎么办?
H s e e Hsee Hsee矩阵是实对称矩阵( ∂ f 2 ( X ) ∂ x j ∂ x i = ∂ f 2 ( X ) ∂ x i ∂ x j \frac{\partial f^2(X)}{\partial x_j \partial x_i} = \frac{\partial f^2(X)}{\partial x_i \partial x_j} ∂xj∂xi∂f2(X)=∂xi∂xj∂f2(X)),而实对称矩阵一定可以三角化
∇ 2 f ( X ( k ) ) = U T Λ U = [ λ 1 0 ⋯ 0 0 λ 2 ⋯ 0 ⋮ ⋮ ⋱ ⋮ 0 0 ⋯ λ n ] , U T U = I 如果 ∇ 2 f ( X ( k ) ) 非正定,说明 λ 1 ∼ λ n 中有若干个特征值小于 0 \begin{aligned} \end{aligned}\begin{aligned} &\nabla^2f(X^{(k)})=U^T\Lambda U= \begin{bmatrix} \lambda_1 & 0 & \cdots & 0\\ 0 & \lambda_2 & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0&\cdots&\lambda_n \end{bmatrix} ,U^TU=I \\ &如果\nabla^2f(X^{(k)})非正定,说明\lambda_1 \sim \lambda_n中有若干个特征值小于0 \end{aligned} ∇2f(X(k))=UTΛU= λ10⋮00λ2⋮0⋯⋯⋱⋯00⋮λn ,UTU=I如果∇2f(X(k))非正定,说明λ1∼λn中有若干个特征值小于0
用最小的特征值 λ m i n ( λ m i n < 0 ) \lambda {min}(\lambda {min}<0) λmin(λmin<0)对 H e s s Hess Hess矩阵进行修正
∇ 2 f ( X ( k ) ) = U T Λ U + ( ε − λ m i n ) I = U T [ λ 1 0 ⋯ 0 0 λ 2 ⋯ 0 ⋮ ⋮ ⋱ ⋮ 0 0 ⋯ λ n ] U + [ ε − λ m i n 0 ⋯ 0 0 ε − λ m i n ⋯ 0 ⋮ ⋮ ⋱ ⋮ 0 0 ⋯ ε − λ m i n ] U T U = U T [ λ 1 + ε − λ m i n 0 ⋯ 0 0 λ 2 + ε − λ m i n ⋯ 0 ⋮ ⋮ ⋱ ⋮ 0 0 ⋯ λ n + ε − λ m i n ] U = U T [ Λ + ( ε − λ m i n ) I ] U \begin{aligned} \end{aligned}\begin{aligned} &\nabla^2f(X^{(k)})=U^T\Lambda U + (\varepsilon - \lambda_{min})I\\ &=U^T \begin{bmatrix} \lambda_1 & 0 & \cdots & 0\\ 0 & \lambda_2 & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0&\cdots&\lambda_n \end{bmatrix}U+ \begin{bmatrix} \varepsilon -\lambda_{min} & 0 & \cdots & 0\\ 0 & \varepsilon -\lambda_{min} & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0&\cdots& \varepsilon - \lambda_{min} \end{bmatrix} U^TU \\ & = U^T \begin{bmatrix} \lambda_1 + \varepsilon -\lambda_{min}& 0 & \cdots & 0\\ 0 & \lambda_2 + \varepsilon -\lambda_{min} & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0&\cdots&\lambda_n+\varepsilon -\lambda_{min} \end{bmatrix}U \\ &=U^T[\Lambda+(\varepsilon-\lambda_{min})I]U \end{aligned} ∇2f(X(k))=UTΛU+(ε−λmin)I=UT λ10⋮00λ2⋮0⋯⋯⋱⋯00⋮λn U+ ε−λmin0⋮00ε−λmin⋮0⋯⋯⋱⋯00⋮ε−λmin UTU=UT λ1+ε−λmin0⋮00λ2+ε−λmin⋮0⋯⋯⋱⋯00⋮λn+ε−λmin U=UT[Λ+(ε−λmin)I]U
为了不和原始的 H e s s Hess Hess矩阵偏差太大, ε \varepsilon ε越小越好
Levenberg-Marquardt修正后,即保证了特征值都是正数,也保留原 H e s s Hess Hess矩阵尽可能多的信息
在工程中,为了减少算法计算的复杂度,不会计算特征值特征向量,而是根据经验值手动设置一个 μ k \mu_k μk,同时还会引入一个步长因子 α \alpha α
X = X ( k ) − α ( k ) [ ∇ 2 f ( X ( k ) ) + μ k I ] − 1 ∇ f ( X ( k ) ) X=X^{(k)}-\alpha^{(k)}[\nabla^2f(X^{(k)})+\mu_kI]^{-1}\nabla f(X^{(k)}) X=X(k)−α(k)[∇2f(X(k))+μkI]−1∇f(X(k))
通过手动调节 μ k \mu_k μk的值,使得 [ ∇ 2 f ( X ( k ) ) + μ k I ] > 0 [\nabla^2f(X^{(k)})+\mu_kI] >0 [∇2f(X(k))+μkI]>0
( μ k → 0 \mu_k \rightarrow 0 μk→0:趋近原牛顿法;$\mu_k \rightarrow \infty $:趋近步长很小的梯度下降法)
二、高斯-牛顿法
1、应用范围
高斯-牛顿法用于解决什么问题?
有一个函数 y = A s i n ( α t + β ) y=\color{red}A\color{black} sin(\color{red}\alpha\color{black} t+\color{red}\beta\color{black}) y=Asin(αt+β) ,其中 A 、 α 、 β \color{red} A、\alpha、\beta A、α、β未知,已知一些输入输出数据 [ t 1 , y 1 ] , [ t 2 , y 2 ] , ⋯ , [ t n , y n ] [t_1,y_1],[t_2,y_2],\cdots,[t_n,y_n] [t1,y1],[t2,y2],⋯,[tn,yn]
高斯-牛顿法想要解决的问题是如何根据已知数据,估计未知参数
这是一个非线性最小二乘问题, min A ^ , α ^ , β ^ ∑ i = 1 n ( A ^ s i n ( α ^ t i + β ^ ) − y i ) 2 \min_{\hat{A},\hat{\alpha},\hat{\beta}} \sum_{i=1}^{n}(\hat{A}sin(\hat{\alpha} t_i+\hat {\beta})-y_i)^2 minA^,α^,β^∑i=1n(A^sin(α^ti+β^)−yi)2
2、高斯-牛顿法原理
考虑更加一般的情况:
min ∑ i = 1 m ( r i ( X ) ) 2 令 r = [ r 1 , r 2 , ⋯ , r m ] T , 则目标函数为 f ( X ) = r ( X ) T r ( X ) , 为了使用牛顿法求解,需要计算梯度和 H e s s 矩阵 梯度 ∇ f ( X ) 的第 j 个元素为 : ( ∇ f ( X ) ) j = ∂ f ∂ x j ( X ) = 2 ∑ i = 1 m r i ( X ) ∂ r i ∂ x i ( X ) r 的 J a c o b i 矩阵为: J ( X ) = [ ∂ r 1 ∂ x 1 ( X ) ∂ r 1 ∂ x 2 ( X ) ⋯ ∂ r 1 ∂ x n ( X ) ∂ r 2 ∂ x 1 ( X ) ∂ r 2 ∂ x 2 ( X ) ⋯ ∂ r 2 ∂ x n ( X ) ⋮ ⋮ ⋱ ⋮ ∂ r m ∂ x 1 ( X ) ∂ r m ∂ x 2 ( X ) ⋯ ∂ r m ∂ x n ( X ) ] 因此,函数 f 的梯度可表示为: ∇ f ( X ) = 2 J ( X ) T r ( X ) 函数 f 的 H e s s 矩阵的第 ( k , j ) 个元素为: ∂ 2 f ∂ x k ∂ x j ( X ) = ∂ ∂ x k ( ∂ f ∂ x j ( X ) ) = ∂ ∂ x k ( 2 ∑ i = 1 m r i ( X ) ∂ r i ∂ x i ( X ) ) = 2 ∑ i = 1 m ( ∂ r i ∂ x k ( X ) ∂ r i ∂ x j ( X ) + r i ( X ) ∂ 2 r i ∂ x k ∂ x j ( X ) ) 令 S ( X ) 表示一个矩阵其中 ( k , j ) 的元素为 ∑ i = 1 m r i ( X ) ∂ 2 r i ∂ x k ∂ x j ( X ) f ( x ) 的 H e s s 矩阵可以表示为: ∇ 2 f = 2 ( J ( X ) T J ( X ) + S ( X ) ) 迭代公式为: X = X ( k ) − ( J ( X ) T J ( X ) + S ( X ) ) − 1 J ( X ) T r ( X ) 由于 S ( X ) 包含函数 r 的二阶导,数值较小可以忽略,所以迭代公式可变为: X = X ( k ) − ( J ( X ) T J ( X ) ) − 1 J ( X ) T r ( X ) \begin{aligned} &\min \sum_{i=1}^{m}(r_i(X))^2 \\ &令r=[r_1,r_2,\cdots,r_m]^T,则目标函数为f(X)=r(X)^Tr(X),为了使用牛顿法求解,需要计算梯度和Hess矩阵\\ &梯度\nabla f(X)的第j个元素为:(\nabla f(X))_j = \frac{\partial f}{\partial x_j}(X) = 2\sum_{i=1}^{m}r_i(X)\frac{\partial r_i}{\partial x_i}(X)\\ &r的Jacobi矩阵为:J(X)= \begin{bmatrix} \frac{\partial r_1}{\partial x_1}(X) & \frac{\partial r_1}{\partial x_2}(X) & \cdots &\frac{\partial r_1}{\partial x_n}(X) \\ \frac{\partial r_2}{\partial x_1}(X) & \frac{\partial r_2}{\partial x_2}(X) & \cdots &\frac{\partial r_2}{\partial x_n}(X) \\ \vdots &\vdots &\ddots &\vdots\\ \frac{\partial r_m}{\partial x_1}(X) & \frac{\partial r_m}{\partial x_2}(X) & \cdots &\frac{\partial r_m}{\partial x_n}(X) \\ \end{bmatrix}\\ &因此,函数f的梯度可表示为:\nabla f(X) = 2J(X)^Tr(X) \\ \\ &函数f的Hess矩阵的第(k,j)个元素为:\\ & \frac{\partial^2f}{\partial x_k \partial x_j}(X) =\frac{\partial}{\partial x_k}\left ( \frac{\partial f}{\partial x_j}(X) \right ) = \frac{\partial}{\partial x_k}\left ( 2\sum_{i=1}^{m}r_i(X)\frac{\partial r_i}{\partial x_i}(X)\right ) = 2\sum_{i=1}^{m}\left( \frac{\partial r_i}{\partial x_k}(X)\frac{\partial r_i}{\partial x_j}(X) +\color{blue} r_i(X)\frac{\partial^2r_i}{\partial x_k \partial x_j}(X) \color{black} \right) \\ &令\color{blue}S(X)\color{black}表示一个矩阵其中(k,j)的元素为\color{blue}\sum_{i=1}^{m} r_i(X)\frac{\partial^2r_i}{\partial x_k \partial x_j}(X) \color{black} \\ &f(x)的Hess矩阵可以表示为:\nabla^2f=2\left(J(X)^TJ(X)+\color{blue} S(X)\color{black} \right)\\ &迭代公式为:X = X^{(k)}-\left(J(X)^TJ(X)+\color{blue} S(X)\color{black} \right)^{-1}J(X)^Tr(X)\\ &由于S(X)包含函数r的二阶导,数值较小可以忽略,所以迭代公式可变为:\\ &X = X^{(k)}-\left(J(X)^TJ(X)\right)^{-1}J(X)^Tr(X)\\ \end{aligned} mini=1∑m(ri(X))2令r=[r1,r2,⋯,rm]T,则目标函数为f(X)=r(X)Tr(X),为了使用牛顿法求解,需要计算梯度和Hess矩阵梯度∇f(X)的第j个元素为:(∇f(X))j=∂xj∂f(X)=2i=1∑mri(X)∂xi∂ri(X)r的Jacobi矩阵为:J(X)= ∂x1∂r1(X)∂x1∂r2(X)⋮∂x1∂rm(X)∂x2∂r1(X)∂x2∂r2(X)⋮∂x2∂rm(X)⋯⋯⋱⋯∂xn∂r1(X)∂xn∂r2(X)⋮∂xn∂rm(X) 因此,函数f的梯度可表示为:∇f(X)=2J(X)Tr(X)函数f的Hess矩阵的第(k,j)个元素为:∂xk∂xj∂2f(X)=∂xk∂(∂xj∂f(X))=∂xk∂(2i=1∑mri(X)∂xi∂ri(X))=2i=1∑m(∂xk∂ri(X)∂xj∂ri(X)+ri(X)∂xk∂xj∂2ri(X))令S(X)表示一个矩阵其中(k,j)的元素为i=1∑mri(X)∂xk∂xj∂2ri(X)f(x)的Hess矩阵可以表示为:∇2f=2(J(X)TJ(X)+S(X))迭代公式为:X=X(k)−(J(X)TJ(X)+S(X))−1J(X)Tr(X)由于S(X)包含函数r的二阶导,数值较小可以忽略,所以迭代公式可变为:X=X(k)−(J(X)TJ(X))−1J(X)Tr(X)
3、Levenberg-Marquardt修正
高斯-牛顿法在使用的过程中也会出现 J ( X ) T J ( X ) J(X)^TJ(X) J(X)TJ(X)不是正定矩阵的问题,所以也可以使用Levenberg-Marquardt修正就加以解决:
X = X ( k ) − ( J ( X ) T J ( X ) + μ k I ) − 1 J ( X ) T r ( X ) X = X^{(k)}-\left(J(X)^TJ(X) + \mu_k I\right)^{-1}J(X)^Tr(X) X=X(k)−(J(X)TJ(X)+μkI)−1J(X)Tr(X)