fmincon函数浅析（原创）

clarkyeah

fmincon命令浅析，有错漏之处望大家不吝赐教
命令格式：
[x,fval,exitflag,output,lambda,grad,hessian] = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon,options)

如matlab帮助文档中所述，fmincon命令使用的算法对于大规模优化问题和中等问题是有所区分的：
Large-Scale Optimization
The large-scale algorithm is a subspace trust region method and is based on the interior-reflective Newton method described in [1] and [2]. Each iteration involves the approximate solution of a large linear system using the method of preconditioned conjugate gradients (PCG)..
Medium-Scale Optimization
fmincon uses a sequential quadratic programming (SQP) method. In this method, the function solves a quadratic programming (QP) subproblem at each iteration. An estimate of the Hessian of the Lagrangian is updated at each iteration using the BFGS formula. A line search is performed The QP subproblem is solved using an active set strategy.
这里试图回答三个问题：
1.
什么Large-Scale Optimization，什么是Medium-Scale Optimization？
2.
fimcon提供的subspace trust region和sequential quadratic programming方法原理？
3.
BFGS公式和线性搜索是什么？

问题1
所谓大规模问题指的是出现在工程，化学等领域中有大量优化变量的问题。由于自变量的维数很高，这样的问题是被分解成多个低维子问题来求解的。Medium－Scale优化问题实际上是matlab自己提出和大规模问题对应的一个概念，就是通常一般的优化算法，如牛顿法，最速下降法之类的处理优化变量不是很多的问题。
问题2
对于大规模问题，fmincon采用了subspace trust region优化算法。这种算法是把目标函数在点x的邻域泰勒展开（x可以认为是人为提供的初始猜测），这个展开的邻域就是所谓的trust region，泰勒展开进行到二阶项为止：
Q(x) = 1/2*<x,Hx> + <f,x>






（1）
这时目标函数在某一个局部的特性就可以“看出来了”。在这样的一个邻域里，我们求一个新的点x1，使得目标函数值减小，这个问题相比于原来的问题要简单。然而实际上对于存在非常大规模优化变量的问题，直接对这个子问题的求解仍然是不可忍受的。
同时我们注意到，由于泰勒展开要进行的第二项，这就要求我们能够提供一阶导计算的函数。如果我们不能提供一阶导表达式，二阶导（Hessian矩阵）matlab是无法计算的。所以我们使用fmincon命令而不给一阶导表达式，fmincon会放弃使用大规模算法。
如前所述，原问题转化后的直接求解仍然是无法忍受的，通过进一步近似subspace trust region将这个问题局限在trust region的二维子空间内求解。
序列二次规划方法是将一个带有等式和不等式约束（可以是非线性）的非线性优化问题转化为二次规划问题求解，二次规划问题类似公式（1）形式。具体转化过程可以参考：
http://www.caam.rice.edu/~adpadu/talks/sqp1.pdf


问题3
对于medium－scale问题，求解二次规划问题涉及到Hessian矩阵。Hessian矩阵的近似计算是通过拟牛顿法得到的，拟牛顿法提供了两个公式可用于Hessian矩阵（或其逆）的迭代：BFGS公式和DFP公式），而初始的Hessian矩阵是任意给的，如给一个单位阵I。
BFGS公式如下：
H(k+1) = H(k) + <q(k),q(k)>/<q(k),s(k)> - <s(k)H(k), s(k)H(k)>/<s(k), H(k)s(k)>
（3）
总结：
fmincon运行首先检查有无梯度表达提供，如有则选则大规模算法（subspace trust region），由此涉及到Hessian阵的近似计算，由于已提供了梯度的公式，则Hessian阵可以直接通过有限差分计算。但是如果用户直接提供了Hessian计算公式，则直接计算。
如果没有梯度表达式提供，fmincon选则SQP算法，算法中Hessian阵可以通过BFGS迭代，初始Hessian阵任给。注意BFGS公式中q项是需要计算目标函数梯度得到的。所以Hessian矩阵的近似计算是需要用到有限差分法。

bie3left

正在理解中 谢谢

liangc

楼主总结的非常好
BFGS公式中q项是需要计算目标函数梯度得到的，
如果没有梯度表达式提供
则q项是用有限差分法计算得来的。
请问以上理解是否正确？

clarkyeah

呵呵，
不用有限差分法计算Hessian除非用户直接提供Hessian的表达式，即使用户提供了梯度表达式，则程序利用梯度表达式来有限差分计算q项,公式如下：
q(k)  = grad(f(x(k+1))) - grad(f(x(k)))

错漏之处还请指教

laughfive

fminunc里如果用medium scale也是用的这个算法吧？
运行首先检查有无梯度表达提供，如有则选则大规模算法（subspace trust region），
如果没有梯度表达式提供，fmincon选择BFGS或DFP 拟牛顿法配合二次三次混合线性搜索。按你的意思，BFGS或DFP就是用来计算HESSIAN的，而优化方法是拟牛顿法，里面用到的搜索方法是二次和三次混合的线性搜索方法。这样理解对吗？我正在看这部分程序，想看懂fminunc（里面的fminusub）里面的内容。能否具体讲下算法？