R中最优化函数optim

郭岚 2019-06-25

最优化函数optim

目标函数:

$$f(x_1,x_2)=(1-x_1)^2+100(x_2-x_1^2)^2$$

该函数全局最小值在($x_1=1,x_2=1$)时取到。

下面这种写法是因为有多个自变量函数,传入一个参数x,每个自变量用向量x的分量来表示,从而定义出目标函数。

obj <- function(x){
  x1 <- x[1]
  x2 <- x[2]
  y <- 100*(x2-x1^2)^2+(1-x1)^2
  return(y)
}

x1梯度:$-400*x_1*(x_2-x_1^2)-2*(1-x_1)$

x2梯度:$200*(x_2-x_1^2)$

梯度

gradient <- function(x){
  x1 <- x[1]
  x2 <- x[2]
 c(x1_gradient = -400*x1*(x2-x1^2)-2*(1-x1),
   x2_gradient =200*(x2-x1^2))
}

tip:对于多元函数,需要用向量的形式来输出各个变量上的梯度。

下面以$x_1=0,x_2=3$作为起始点,其他参数采用默认设置,梯度也不输入。(默认时优化算法为单纯型法)

optim(par = c(0,3),fn=obj)

R中最优化函数optim

  • par元素表示最优解取值
  • value值表示目标函数值
  • counts代表调用目标函数与梯度函数的数目,可认为是迭代数目
  • convergence是收敛的代码

    • 0表示成功完成了优化任务
    • 1表示达到了达到了迭代上限而退出,关于迭代上限会在control参数的maxit元素进行说明
    • 10表示退化(退化表示单纯型无法继续移动)的单纯型
    • 51专指优化算法为L-BFGS-B的时候输出的警告信息。(L-BFGS-B为拟牛顿改进方法)
    • 52专指优化算法为L-BFGS-B的时候输出的错误信息。

optim函数参数详解

  • par:表示各变量的初始值,通过向量传入。对于很多算法,初始值非常重要,因此在具体的优化问题中如果能找到最优解附近的初始值(比如之前的经验)传入优化函数,将能大大的增强优化的效率。如果没有合适的初始值,可以选择任意的向量传入。
  • fn:目标函数
  • gr:梯度函数
  • method:优化方法,通过字符串形式传入,表示优化求解的算法,

    • Nelder-Mead:单纯型法,为optim默认优化算法。

      • 思想:通过单纯型的方式不断替换函数的最差的顶点从而得到最优值。因为没有用的梯度故不是非常有效,但方法十分稳健,效率也不低,因此被作为默认算法。
    • CG:是一种共轭梯度法,这种方法充分利用函数的梯度信息,在每一点都能找到一个最合适的方向(类似的算法是最速下降法,直接按照目标函数值最快的下降方向来搜索)来搜索,该方法的方向不一定是下降的方向,因此能避免陷入局部最优的困境。
    • BFGS:是一种拟牛顿法,也称变尺度法。该算法改进了牛顿法中容易受初值的影响的弱点,但是又不需要在每一步优化的过程中计算精确的hessian矩阵及其逆矩阵,在具备牛顿法搜索快的特性的基础上又能有效的搜索全局最优解,一次使用十分广泛,是optim函数中应用最广的算法。
    • L-BFGS-B:是对BFGS算法的一个优化,能够在优化的同时增加箱型约束条件,一定程度上增强了这些无约束的非线性规划方法的功能。
    • SANN是一种模拟退火的方法,与通常的数学函数的算法不同,该算法是一种概率算法,对于各种复杂的情况,尤其是很多不可微的函数,该算法可以找到最优解,但效率上比不上其他数学算法。
    • Brent算法是一种简单的一维搜索算法,通常被其他函数调用,实际使用中几乎不用。
  • lower:当算法选择为L-BFGS-B时,该函数允许传入简单的箱式约束,也就是说变量大于某个实数且小于某个实数,lower表示下界,upper表示上界,都是通过这个传入的。
  • upper:参见lower
  • control:该参数是一个列表,包含优化中的各种设置,很多其他第三方的优化函数也遵循这样的设置方式。常见的设置元素包括最大迭代次数maxit、绝对收敛容忍度abstol、相对收敛容忍度reltol等。详情可以通过?optim来得到。
  • hessian:表示十分返回hessian矩阵,默认为FALSE,但是hessian矩阵对其他的运算还是非常重要的,比如估计参数的置信区间。

例子

使用CG共轭梯度法在默认的迭代次数下求解:

optim(par = c(0,3),fn=obj)

R中最优化函数optim

发现算法到达迭代次数上限而退出。

我们增加迭代次数:

optim(par = c(0,3),fn=obj,method = "CG",control = list(maxit=500))

R中最优化函数optim

可以发现目标函数的值可以减少,进一步增加迭代的次数可以得到最优解,但是效率上比默认的Nelder-Mead单纯型法差太多。注意:CG已经使用梯度信息,不需要我们自己传入自定义的gradient函数

如果我们使用拟牛顿法BFGS,下面分别为不传入梯度信息与自己传入梯度信息

optim(par = c(0,3),fn=obj,method = "BFGS")

R中最优化函数optim

optim(par = c(0,3),fn=obj,gr=gradient,method = "BFGS")

R中最优化函数optim

R中最优化函数optim