1.4 优化方法

一、优化问题的数学形式:

min \ f(x_{1},x_{2},..x_{n}) \\ subject \ to \left\{\begin{matrix} g_{i}\leq 0 ,&i=1,2...m \\ h_{j}=0,&j=1,2...l \end{matrix}\right.。根据问题的形式给出不同的分类。

二、分类:(自己的分类,仅供参考)

1.无约束优化:没有不等式约束和等式约束。

无约束优化问题的解法:

(0)撒点法(网格搜索法,随机模拟,智能算法(遗传算法等))。自己选择叫这个名字,就是好记而已,网格相当于均匀撒点,随机模拟是利用随机性,而智能算法在撒完点后,进行特定的选择。

(1)直接法(免导数,直接算f(x)的相关性质):Powell方向加速法,步长加速法(模式搜索法)。

(2)间接法(基于导数信息):

下降法:下降法的主要思想是:利用\mathbf{x_{k+1}=x_{k}+\mu _k\Delta x_k},在给定一个初始值后,通过迭代,寻找最优点\mathbf{x_{opt}}.

在式子中,\mathbf{x_k,\Delta x_k}分别是第k次的迭代点以及下一步要进行搜索的方向,它们都是n维向量,而\mu _k\geqslant 0是一个标量,为第k次迭代的步长。由于最小化算法要求对于每一步迭代,都有:{f(x_{k+1})<f(x_{k}),所有在Xk处展开:\mathbf{f(x_{k+1})\approx f(x_k)+(\triangledown f(x_k))^T\Delta x_k},即(\triangledown f(x_k))^T\Delta x_k}< 0。当下一步的方向是梯度负方向时,成为最速下降法。

进一步的,利用黑塞矩阵,进行二次逼近的相应方法称为牛顿法。算法的准确描述如下:

初始化:选择一个初始点x_{1}\in dom f和允许的精度\varepsilon >0,另k=1;

步骤1:计算目标函数在点x_{k}的梯度\triangledown f(x_k)(以及黑塞矩阵\triangledown^2f(x_k)),并选择下降方向\Delta x_{k}=\left\{\begin{matrix} -\triangledown f(x_{k}),\\-(\triangledown ^2f(x_{k}))^{-1} \triangledownf(x_{k}) \end{matrix}\right.(根据所选择方法不同进行计算)

步骤2:选择步长\mu_{k}> 0.

步骤3:进行更新\mathbf{x_{k+1}=x_{k}+\mu _k\Delta x_k}

步骤4:判断停止准则是否满足:若|f(x_{k+1}-f(x_k))|\leq \varepsilon,则停止迭代,并输出Xk;

若不满足,则另k<-k+1,并返回步骤1,进行下一次迭代,直到停止准则满足为止。

其中,根据步长的选择不同,也有一些变种。(引用自《矩阵分析与应用》第二版4.4梯度分析)

当x是有约束的时候,梯度算法中的更新公式要用投影代替,这时候算法就称为梯度投影法。

详细的知识可以查看《矩阵分析与应用》第二版4.4梯度分析那一章,因为接下来的知识看不懂,就没有细看。

因为最近在看论文,本着不会什么查什么的理念,对于平滑的梯度法就只看了这些。

接着就是,对于不平滑的函数来说,梯度不存在,就用次梯度去代替:

对于平滑函数的梯度有:\mathbf{f(y)\geqslant f(x)+(\triangledown f(x))^T(y-x)},\mathbf{\forall x,y\in dom(f(x))},

当函数不光滑时:when \mathbf{\ x\in R^n},\mathbf{f(y)\geqslant f(x)+g^T(y-x)},\mathbf{\forallyy\in dom(f)},\mathbf{g}就称为函数\mathbf{f}\mathbf{x}处的次梯度,记作:\partial f(x).

需要注意的是:\mathbf{g}是一个向量的集合。

一个例子就是:f(x)=|x|,x\in R.\\\partial |x| = \left\{\begin{matrix} {-1},x<0\\{+1},x>0 \\ [-1,1],x=0\end{matrix}\right.。对于一个不平滑的凸函数,有一条结论是这样的:

若点x^*是凸函数f的一个极小点,当且仅当f在x^*可次微分并且\mathbf{0\in \partial f(x^*)}(一阶条件)。

当上述条件中f为可微分时,一阶条件退化为梯度等于\boldsymbol{0}.

一个比较常用的结论是:

用语言表示一下就是:当我们想求\frac{1}{2}(x-z)^2+\lambda|x|,x\in R^1的极小点时,有闭式解形式:x=sign(z)*max(|z|-\lambda,0)

 

2.有约束优化:在梯度分析中,有约束的优化算法标准思想都是将约束优化问题转化成无约束优化问题。转化的方法主要有三种:①拉格朗日乘子法 ②罚函数法 ③增广拉格朗日乘子(拉格朗日乘子与罚函数相结合的方法)。但是我的分类稍微有所不同:因为这篇博客中线性规划要敲对应的代码,所以分成一类。非线性规划时,再使用上述所提的三种方法,当然这种分类并不正确,因为线性规划也使用后面的三种方法进行求解。

(0)(仅有不等式约束,无等式约束时)撒点法:步骤和1中撒点法一样,只不过要增加约束的条件。

(1) 线性规划类:当目标函数以及两类约束条件都是线性的时候,称为线性规划(linear programme)。

a.线性规划标准形式(linprog):

min\ \mathbf{c^Tx} \ s.t.\left\{\begin{matrix} \mathbf{Ax}\leq \mathbf{b}\\\mathbf{A_{eq}x=b_{eq}} \\ \mathbf{Lb\leq x \leq Ub}\end{matrix}\right.,x是列向量,c是目标函数系数。

matlab调用:[xmin,Fmin]=linprog(\mathbf{x,A,b,Aeq,beq,Lb,Ub});

b.线性规划的一个变种:当我们想优化一个矩阵X时:即min\ \mathbf{tr(C^TX)} \ s.t.\left\{\begin{matrix} \mathbf{AX}\leq \mathbf{B}\\\mathbf{A_{eq}x=B_{eq}} }\end{matrix}\right.。(这里为了方便,省略掉上下界的约束)。因为matlab库中的linprog只能对向量\mathbf{x}进行优化,所以我们将X向量化,即我们对vec(X)进行优化,得到(vec(X))^*后,再转化为X^*。利用\mathbf{tr(C^TX)} =\mathbf{vec(C)^Tvec(X)}这一性质,调用linprog时,系数矩阵要变成vec(C)^T,其余的等式约束和不等式约束也要进行相关的处理。

c.整数规划(intlinprog):标准形式:min\ \mathbf{c^Tx} \ s.t.\left\{\begin{matrix} Ax\leqslant b\\ A_{eq}x=b_{eq} \\ Lb\leqslant x\leqslant Ub \\ x_i\in Z,\in I \end{matrix}\right.,

matlab调用:[xmin,Fmin]=linprog(\mathbf{x,A,I,b,Aeq,beq,Lb,Ub}),其中,\mathbf{I}是一个行向量。

(2)非线性规划:当三者至少有一个不为线性时,称为非线性规划。接下来是相关算法:

①拉格朗日乘子法:

原问题的一般形式:min\ f(\mathbf{x}) \\ s.t.\left\{\begin{matrix} f_i(\mathbf{x})\leqslant 0,i=1,...,m\\ h_{i}(\mathbf{x})=0,i=1,...,p \end{matrix}\right.

记拉格朗日函数:L(\mathbf{x},\lambda_i,\nu _i)=f_0(\mathbf{x})+\sum_{i=1}^{m}\lambda_if_i(\mathbf{x})+\sum_{i=1}^{p}\nu_ih_i(\mathbf{x}),\lambda_i\geq 0

写成向量形式:L(\mathbf{x,\lambda,\nu})=f(\mathbf{x})+\mathbf{\lambda^Tg(x)}+\mathbf{\nu^Th(x)},\lambda^T=[\lambda_1,...\lambda_m]\geq \mathbf{0^T},\nu^T=[\nu_1,...\nu_p].其中,\lambda ,\nu称为对偶变量或是乘子向量。

接着就是定义拉格朗日对偶函数:g(\mathbf{\lambda,\nu})=inf_{x\in D} \ L(\mathbf{x,\lambda,\nu}),其中inf表示下界。

在对偶函数中,先假设lambda与nu是一个参数,对于x在D中,求下界意思就是对x求最小值,因为此时lambda和nu是一个参数,所以对偶函数是\mathbf{g(\lambda,\nu)}是一个关于lambda和nu的函数。

设原问题最小值为:p^*=min \ f_0(x),则对偶函数满足:g(\mathbf{\lambda,\nu})\leq p^*。证明如下:

拉格朗日函数L(\mathbf{x},\lambda_i,\nu _i)=f_0(\mathbf{x})+\sum_{i=1}^{m}\lambda_if_i(\mathbf{x})+\sum_{i=1}^{p}\nu_ih_i(\mathbf{x})\leqslant f_0(\mathbf{x}),因为第二项是负值,而第三项是0。

所以必定有:g(\mathbf{\lambda,\nu})=inf_{x\in D} \ L(\mathbf{x,\lambda,\nu})\leq L(\mathbf{x},\lambda_i,\nu _i)\leqslant p^*

拉格朗日对偶问题为:max \ g(\mathbf{\lambda,\nu})\\ s.t.\mathbf{ \lambda>0},

d^*=sup_{\mathbf{(\lambda\geqslant 0,\nu)}} g(\mathbf{x,\lambda,\nu})=max_{\mathbf{(\lambda\geqslant 0,\nu)}}min_{\mathbf{x}L(x,\lambda,\nu),即d^*是对偶问题的最优解,d^*\leq p^*

无论原问题是何种问题,它的对偶问题都是一个凸问题。

当只满足对偶间隙p^*-d^*\geq 0时,称拉格朗日对偶法为弱对偶的。再进一步,如果对偶间隙为0时,称为强对偶的。

强弱对偶的判断是根据Slater准则。多数问题都是强对偶的。

在强对偶的条件下,如果原问题(无论凹凸)涉及到的函数是可微的,那么有KKT条件:

\left\{\begin{matrix} f_i(x^*)\leq 0,i=1,...,m\\ h_i(x^*)=0,i=1,...,p \\ \lambda_i^* \geq 0,1,...,m \\ \lambda_i^*f_i(x^*)=0,i=1,...,m\\\triangledown f_0(x^*)+\sum_{i=1}^{m} \lambda_i^*\triangledown f_i(x^*)+\sum_{i=1}^{q}\nu_i\triangledown h_i(x^*)=0 \end{matrix}\right.

KKT条件中,第一个原始不等式约束,第二个是原始等式约束,第三个是对偶约束,第四个是互补松弛性,第五个是稳定性。

对于任何的原问题来说,KKT条件是必要的。而当原问题是凸问题时,KKT条件是充要的。

一个例子如下:考虑min f(x) \ s.t. \ Ax=b的优化问题

构造拉格朗日函数:L(x,\nu)=f(x)+\nu^T(Ax-b)。假设最优值点为:(x^*,\nu^*),必有:

(x^*,\nu^*)=arg \ max_{\nu} \ min_x L(x,\nu),且(x^*,\nu^*)=arg \ \min_{x} \max_{\nu} L(x,\nu)

当我们知道x^*,\nu^*其中一个变量时,可以通过带入已知的变量解无约束的拉格朗日函数,从而得到另一变量。

\nu^*已知时,x^* = arg \min_{x}L(x,\nu^*);

假设我们现在已知当前的x_k\nu^*,我们想要求x^*,利用无约束优化中的梯度下降法:求一个x_{k+1}使得L(x,\nu^*)降低,

x_{k+1}=x_{k}+\mu_k(-\triangledown L(x,\nu^*)),在此例子中为:x_{k+1}=x_{k}+\mu_k(-\triangledown L(x,\nu^*))=x_{k}-\mu_k(\triangledown f(x)+A^T\nu^*)

因为v^*是未知的,所以我们以\nu_{k}去替代,所以得到求x_{k+1}的迭代公式:x_{k+1}=x_{k}-\mu_k(\triangledown f(x)+A^T\nu_k).

x^*已知时,\nu^* = arg \max_{\nu}L(x^*,\nu)

同理(注意我们此时求max值,所以直接选择梯度方向即可),我们得到下一次的迭代公式:\nu_{k+1}=v_k+\mu_k\triangledown _{\nu}L(x^*,\nu)

在此例中,为\nu_{k+1}=v_k+\mu_k(Ax_k-b)

综上所述有:x_{k+1}=x_{k}-\mu_k(\triangledown f(x)+A^T\nu_k) \\ \nu_{k+1}=\nu_k+\mu_k(Ax_k-b),其中\mu_k就是下一次迭代的步长。

这个转化后通过迭代求得在无约束下极值点的方法叫做对偶上升法。

但是这个方法有一些缺点(比如说:收敛慢等瑕疵),将算法进行改进,就得到了增广拉格朗日法(见下文)。

接下来的是:拉格朗日乘子法的几何解释:

当只有等式约束时,就退化成了梯度平行的样子:一个很好的几何解释是:梯度共线:

(图片来源:四川大学徐小湛高等数学教学视频)

当有不等式约束时,几何解释较好的是:

如何通俗地讲解对偶问题?尤其是拉格朗日对偶lagrangian duality? - 彭一洋的回答 - 知乎

非线性优化中的 KKT 条件该如何理解? - 王小龙的回答 - 知乎

②)罚函数:

问题形式:min\ f(\mathbf{x}) \\ s.t.\left\{\begin{matrix} f_i(\mathbf{x})\geq 0,i=1,...,m\\ h_{i}(\mathbf{x})=0,i=1,...,p \end{matrix}\right.,注意到标准形式的罚函数中,不等式约束的方向!

构造罚函数:F(x,m_k)=f_0(x)+m_kP(x) \\P(x)=\sum \phi (f_i(x))+\sum \varphi (h_i(x)),mk是惩罚因子,通常是很大的一个数。

构造的罚函数式子中,要满足以下条件:

对于不等式约束的惩罚:y\geq 0,\phi (y)=0;y<0,\phi(y)>0;意思就是:当满足不等式约束时,不予惩罚,所以\phi(y)=0,当不等式约束不满足时,要给与大的惩罚。

对于等式约束的惩罚:y=0,\varphi (y)=0;y\neq 0,\varphi(y)> 0;意思就是:当满足等式约束时,不惩罚,不满足等式时,给予惩罚。

比较常见的选取有:P(x)=(\sum max(0,-f_i(x))^2+\sum [h_i(x)]^2.

罚函数可分为外罚函数和内罚函数。

③增广拉格朗日乘子法(ALM):

a.等式约束:记等式约束为\mathbf{h(x)}=[h_1(x),...h_p(x)]^T.

增广拉格朗日函数:\mathbf{L_{\rho }(x,\nu)=f_0(x)+\nu^Th(x)+\rho\varphi (h(x))},\rho为惩罚因子。

由①中的讨论可以得到:(等式约束下,拉格朗日乘子法的迭代公式)x_{k+1} = arg \min_{x}L_{\rho }(x,\nu_k),\nu_{k+1}=\nu_{k}+\mu_{k}\triangledown _{\nu}L(x_k,\nu_k)

注意到第一个式子(求x_{k+1}的式子)不一定必须用梯度下降法求极小值,用其他方法也可以,所以这里不详细写到下一步。

接下来进行改进:将\lambda_{k+1}=\lambda_{k}+\mu_{k}\triangledown _{\nu}L(x_k,\nu)中的\mu_{k}改成\rho_{k},x_{k}改成x_{k+1}。就得到了求解增广拉格朗日的对偶上升法的更新组成:x_{k+1} = arg \min_{x}L_{\rho }(x,\nu_k)\\\nu_{k+1}=\nu_{k}+\rho _{k}\triangledown _{\nu}L(x_{k+1},\nu_k)

续上文中那个等式约束的例子,使用ALM迭代如下:

L_{\rho }(x,\nu)=f_0(x)+\nu^T(Ax-b)+\frac{\rho}{2}||Ax-b||^2,

对应的更新公式如下:x_{k+1} = arg \min_{x}L_{\rho }(x,\nu_k)\\\nu_{k+1}=\nu_{k}+\rho _{k}(Ax_{k+1}-b).

三、一个常用的技巧:

考虑问题:\min_{x} \ f(x)+g(x),存在f和g分别优化时比较简单,但是组合在一起就比较难直接优化了,比如:\mathbf{||Ax-b||_2^2+||x||_{1}^2}等问题。一个比较可行的解法是:

因为原问题等价于min \ f(x)+g(z)\\s.t. x-z=0,所以记L_{\rho }(x,z,\nu)=f(x)+g(z)+\nu^T(x-z)+\frac{\rho}{2}||x-z||_2^2,根据ALM的对偶上升法

有如下的更新式子:1)\{x_{k+1},z_{k+1}\}=arg\min_{x,z}L_{\rho}(x_{k},z_k,\nu) \\2)\nu_{k+1}=\nu_{k}+\rho(x_{k+1}-z_{k+1})。现在仔细看1)中的内容,问题在于如何更新?

利用坐标轮换法,在第k+1步迭代求\{x_{k+1},z_{k+1}\}时:详细步骤如下:

1)a:\{x_{k+1|t+1}\}=arg \min_{x}(f(x)+\nu^T(x-z_{k+1|t})+\frac{\rho}{2}||x-z_{k+1|t}||_2^2),又因为f(x)+\nu^T(x-z)+\frac{\rho}{2}||x-z||_2^2

=f(x)+\frac{\rho}{2}\{||x-z||_2^2+\frac{2}{\rho}\nu^T(x-z)+\frac{\nu^2}{\rho^2}-\frac{\nu^2}{\rho^2}\}=f(x)+\frac{\rho}{2}\{||x-z||_2^2+\frac{2}{\rho}\nu^T(x-z)+\frac{\nu^2}{\rho^2}\}-\frac{\nu^2}{2\rho}

=f(x)+\frac{\rho}{2}||x-z+\frac{\nu}{\rho}||_2^2-\frac{\nu^2}{2\rho},因为只是求arg\min_{x},所以\{x_{k+1|t+1}\}=arg \min_{x}(f(x)+\frac{\rho}{2}||x-z_{k+1|t}+\frac{\nu}{\rho}||_2^2),同理

1)b:\{z_{k+1|t+1}\}=arg \min_{z}(g(z)+\frac{\rho}{2}||z-x_{k+1|t+1}-\frac{\nu}{\rho}||_2^2).由此算法即可将f(x)与g(z)分离出来,更好算。

在通过内部先迭代算出\{x_{k+1},z_{k+1}\}后,再计算2)中的\nu_{k+1}.继续迭代,直到算法收敛。

 

参考:

1.视频:中国科学技术大学:凌青:最优化理论

2.书籍:矩阵分析与应用 张贤达

3.ppt:电子科技大学:张勇:数学建模系列最优化ppt3:非线性优化

 

 

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值