1 线性规划(lpSolve, linprog)

lpSolve 包中的lp函数,linprog包中的 solveLP函数

2 无约束或区间约束的优化问题

stat::optim/optimize函数

可以利用stat包中的optim()optimize函数

optim函数的例子

#求解局部最优值
f <- function(x) (x[1]^2 + x[2] - 11)^2 + (x[1] + x[2]^2 - 7)^2
optim(par=初值,fn,method=c(),lower=下界,upper=上界)$par
#这个函数需要提供初值,只能实现最小化,可以处理高维
optimize(f, interval, lower = min(interval), upper = max(interval),maximum = FALSE,tol = 误差精度)
#这个函数需要提供区间

result=optim(c(0.2,0.3),f)

https://blog.csdn.net/shen19920619/article/details/55094994
https://segmentfault.com/a/1190000011761522

nlminb函数

nlminb是R语言中用于非线性优化的函数,主要用于在给定约束条件下最小化目标函数。它可以处理有界约束问题,即变量可以有上下界。以下是nlminb函数的基本用法和参数说明:

基本用法

result <- nlminb(start, objective, gradient = NULL, hessian = NULL,  ..., 
						lower = -Inf, upper = Inf, control = list())  

参数说明:

  • start: 初始参数值的向量。
  • objective: 目标函数,它是需要最小化的函数。该函数应该接受一个参数向量并返回一个标量值。
  • gradient: (可选)目标函数的梯度函数。如果提供,应该返回一个与参数向量长度相同的向量。
  • hessian: (可选)目标函数的Hessian矩阵函数。如果提供,应该返回一个对称矩阵。
  • …: 传递给目标函数的其他参数。
  • lower: 参数的下界,可以是单个值或与参数向量长度相同的向量。
  • upper: 参数的上界,可以是单个值或与参数向量长度相同的向量。
  • control: 控制参数的列表,用于调整算法的行为,如收敛标准、最大迭代次数等。

返回值:
nlminb返回一个列表,包含以下元素:

  • par: 优化后的参数值。
  • objective: 在优化参数下的目标函数值。
  • convergence: 收敛状态,0表示成功。
  • iterations: 迭代次数。
  • evaluations: 目标函数和梯度的评估次数。
  • message: 关于优化过程的消息。

示例:
以下是一个简单的示例,演示如何使用nlminb来最小化一个二次函数:

# 定义目标函数  
objective_function <- function(x) {  
  return((x[1] - 1)^2 + (x[2] - 2)^2)  
}  

# 初始参数  
start_params <- c(0, 0)  

# 调用nlminb进行优化  
result <- nlminb(start = start_params, objective = objective_function)  

# 打印结果  
print(result$par)        # 优化后的参数  
print(result$objective)  # 目标函数值  
在这个示例中,目标函数是一个简单的二次函数,nlminb用于找到使该函数最小化的参数值。

3 线性约束的二次规划问题

二次规划的一般形式

在这里插入图片描述

solve.QP函数

solve.QP()是该函数只解决严格的凸二次规划问题,求极小值。

solve.QP(Dmat, dvec, Amat, bvec, meq=0, factorized=FALSE) 
函数参数: 
- Dmat: 为Hessian矩阵 
- dvec: 为向量,和Dmat相对应 
- Amat: 约束的系数矩阵,默认约束为“≥” 
- bvec: 为向量,和Amat相对应 
- meq: 表示从哪一行开始Amat矩阵中的约束是等式,默认为0 
VALUE: 
- solution: 二次规划问题中向量的取值 
- value: 标量,二次规划目标函数的取值 
- unconstrained.solution: 没有约束条件下向量的取值

例子
在这里插入图片描述
代码如下:

D <- matrix(c(1,0,0,2),nr=2)
d <- c(3,1)
A <- t(matrix(c(1,-1,2,1,2,-1),nr=3))##函数表达式中为A^T,而在程序中需要输入A。
b0 <- c(2,1,-4)
solve.QP(Dmat = D,dvec = d, Amat = A,bvec = b0)

$solution
[1] 2.0 1.5

$value
[1] -3.25

$unconstrained.solution
[1] 3.0 0.5

$iterations
[1] 2 0

$Lagrangian
[1] 0 1 0

$iact
[1] 2

在约束条件下,x=2.0,y=1.5x=2.0,y=1.5,目标函数最小值为-3.25;在无约束条件下x=3.0,y=0.5x=3.0,y=0.5。

osqp::osqp函数

这个是OXFORD大学随后提供的优化包。
安装:

install.packages("remotes")
remotes::install_github("r-lib/remotes#103")
remotes::install_git("git://github.com/OxfordControl/osqp-r",submodules = TRUE)

例子:

在这里插入图片描述

library(osqp)
library(Matrix)

# Define problem data
P <- Matrix(c(4., 1.,
              1., 2.), 2, 2, sparse = TRUE)
q <- c(1., 1.)
A <- Matrix(c(1., 1., 0.,
              1., 0., 1.), 3, 2, sparse = TRUE)
l <- c(1., 0., 0.)
u <- c(1., 0.7, 0.7)

# Change alpha parameter and setup workspace
settings <- osqpSettings(alpha = 1.0)
model <- osqp(P, q, A, l, u, settings)

# Solve problem
res <- model$Solve()

来源参考:https://osqp.org/docs/examples/setup-and-solve.html

4 无约束非线性优化问题(nloptr)

nloptr是一个在R中用于非线性优化的包。它提供了一种寻找目标函数的最小值或最大值的方法。nloptr包中的函数可以通过设置不同的参数来使用不同的优化算法。

以下是nloptr包中最常用的函数:

  • nloptr:这是nloptr包中最重要的函数,用于执行非线性优化。它接受初始值、目标函数、约束条件等参数,并返回优化结果。

  • nloptrx:这个函数提供了一个更灵活的接口,可以用于更复杂的优化问题,包括非线性约束和不等式约束。

  • auglag:这个函数可以用于求解带有约束的非线性优化问题,使用增广拉格朗日方法进行求解。

  • slsqp:这个函数是一个特定的优化算法,用于求解带有约束的非线性优化问题。

实例


library(nloptr)

# 定义目标函数
objective <- function(x) {
  return(x^2 + 5*sin(x))
}

# 设置初始值
x0 <- 1

# 进行优化
result <- nloptr(x0 = x0, eval_f = objective, opts = list("algorithm" = "NLOPT_LN_COBYLA"))

# Call:
#   
#   nloptr(x0 = x0, eval_f = objective, opts = list(algorithm = "NLOPT_LN_COBYLA"))
# 
# 
# 
# Minimization using NLopt version 2.7.1 
# 
# NLopt solver status: 4 ( NLOPT_XTOL_REACHED: Optimization stopped 
#                          because xtol_rel or xtol_abs (above) was reached. )
# 
# Number of Iterations....: 24 
# Termination conditions:  relative x-tolerance = 1e-04 (DEFAULT) 
# Number of inequality constraints:  0 
# Number of equality constraints:    0 
# Optimal value of objective function:  -3.24639426722972 
# Optimal value of controls: -1.110552


# 定义目标函数
objective <- function(x) {
  return((x[1]-2)^2 + (x[2]-3)^2)
}

# 设置初始值
x0 <- c(0, 0)

# 进行优化
result <- nloptr(x0 = x0, eval_f = objective, opts = list("algorithm" = "NLOPT_LN_COBYLA"))

# Call:
#   
#   nloptr(x0 = x0, eval_f = objective, opts = list(algorithm = "NLOPT_LN_COBYLA"))
# 
# 
# 
# Minimization using NLopt version 2.7.1 
# 
# NLopt solver status: 4 ( NLOPT_XTOL_REACHED: Optimization stopped 
#                          because xtol_rel or xtol_abs (above) was reached. )
# 
# Number of Iterations....: 39 
# Termination conditions:  relative x-tolerance = 1e-04 (DEFAULT) 
# Number of inequality constraints:  0 
# Number of equality constraints:    0 
# Optimal value of objective function:  1.26010111855728e-08 
# Optimal value of controls: 1.999952 3.000101

4 非线性约束的非线性优化问题

Rsolnp2::donlp2函数

安装代码:

install.packages("Rdonlp2", repos="http://R-Forge.R-project.org")

注意在使用donlp2时upper和lower成对出现,缺省一个就会报错。

    donlp2( 
        par, fn,
        par.upper = rep(+Inf, length(par)), par.lower = rep(-Inf, length(par)),   
        A = NULL,
        lin.upper = rep(+Inf, length(par)), lin.lower = rep(-Inf, length(par)), 
        nlin = list(),
        nlin.upper = rep(+Inf, length(nlin)), nlin.lower = rep(-Inf, length(nlin)),
        control = donlp2Control(),
        control.fun = function(lst){return(TRUE)},
        env = .GlobalEnv, 
        name = NULL)

使用可参考:
https://blog.csdn.net/qq_27755195/article/details/51458659

安装Rdonlp2包不成功的解决方法:
使用命令按章Rdonlp2时提示如下信息:

WARNING: Rtools is required to build R packages but is not currently installed. Please download and install the appropriate version of Rtools before proceeding:

https://cran.rstudio.com/bin/windows/Rtools/
Package which is only available in source form, and may need
  compilation of C/C++/Fortran: ‘Rdonlp2’
  These will not be installed

表明缺少Rtools,需要下载安装并配置环境变量,然后重启R,下载地址为:https://cran.r-project.org/bin/windows/Rtools/,
具体操作可参考:https://www.cnblogs.com/liugh/p/9937489.html
安装完Rtools还不能正常安装Rdonlp2的话,就从这里https://r-forge.r-project.org/R/?group_id=156下载后本地安装。

Rsolnp::solnp的用法

功能
solnp函数,用于求解光滑的非线性约束非线性优化问题,返回目标函数的局部极小值。
形式
solnp(pars, fun, eqfun = NULL, eqB = NULL, ineqfun = NULL, ineqLB = NULL,
ineqUB = NULL, LB = NULL, UB = NULL, control = list(), …)
参数

  • pars: 初值。
  • fun: 目标函数,返回值是一个常数。
  • eqfun: (可选) 返回等式约束的值,多个等式约束返回一个vector.
  • eqB: (可选) 等式约束的约束值,和eqfun的返回结果相一致。
  • ineqfun: (可选) 不等式约束,多个返回由表达式结果组成的向量。
  • ineqLB: (可选) 不等式约束的下界约束,维数和ineqfun返回结果一致。
  • ineqUB:(可选) 不等式约束的上界约束,维数和ineqfun返回结果一致。
  • LB: (可选)优化变量的下界向量。
  • UB: (可选)优化变量的上界向量。
  • control: (可选) 优化器的控制参数。

注意,有些时候在R上进行优化也需要制定一个较大的上下界,以避免函数出错。
control参数

  • rho: 对跑出可行域的目标函数进行的惩罚权重(默认 1).This is used as a penalty weighting scaler for infeasibility in the augmented objective function.
    The higher its value the more the weighting to bring the solution into the feasible region (default 1). However, very high values might lead to numerical ill conditioning or significantly slow down convergence.
  • outer.iter: 最大的major iterations (默认 400).
  • inner.iter: 最大的minor iterations (默认 800).
  • delta: Relative step size in forward difference evaluation (default 1.0e-7).
  • tol: 精度(默认 1e-8).
  • trace:目标函数值和优化变量在每个major iteration中输出的次数(默认 1).
    函数返回值
  • pars: 优化变量的结果。
  • convergenc: 是否优化成功的指标,0表示优化成功,1或2表示没有优化成功。
  • values: 优化迭代过程中的目标函数值组成的向量,最后一个值是最终的优化结果。
  • lagrange: Lagrange乘子组成的向量。
  • hessian: 最优解时的hessian矩阵。
  • ineqx0: 取得最优解时的,不等式约束的松弛变量(slack variables)。
  • nfuneval: The number of function evaluations.
  • elapsed: 优化所耗费的时间。

例子
在这里插入图片描述
代码:

library(pacman)
p_load(Rsolnp)

fn1=function(x)
{
  z=x[1]^2*sin(x[2])+x[2]^2*cos(x[1])
  return(z)
}
eqn1=function(x){
  return(x[1]*x[2])
}
ineqn1=function(x){
  z1=x[1]+x[2]
  z2=3*x[1]-x[2]
  z3=sin(x[1])*cos(x[2])
  return(c(z1,z2,z3))
}
x0 = c(10,10)
powell=solnp(x0, fun = fn1, eqfun = eqn1, eqB = c(2),
             ineqfun = ineqn1,ineqLB = c(2,1,-Inf),
             ineqUB = c(Inf,3,3),LB = c(-100,-100), UB = c(100,100))

结果

> powell
$pars
[1] 1.403075 1.425440

$convergence
[1] 0

$values
 [1] -138.309264  -40.714718   -5.028974    2.438729    2.315228    2.295481
 [7]    2.295413    2.285524    2.287054    2.287053    2.287053

$lagrange
              [,1]
[1,]  5.424216e-01
[2,] -9.412557e-13
[3,] -7.857526e-13
[4,]  3.065161e-12

$hessian
             [,1]         [,2]        [,3]        [,4]        [,5]
[1,]  0.988384217  0.003371369 -0.02097895 -0.51309600 -0.52773026
[2,]  0.003371369  0.156623138 -0.20311846 -0.05872235  0.35866137
[3,] -0.020978952 -0.203118459  0.90354377 -0.99230133 -0.87784508
[4,] -0.513096002 -0.058722354 -0.99230133  5.80401809 -0.06335659
[5,] -0.527730257  0.358661374 -0.87784508 -0.06335659  1.77786334

$ineqx0
[1] 2.8285155 2.7837853 0.1428122

$nfuneval
[1] 236

$outer.iter
[1] 10

$elapsed
Time difference of 0.04886913 secs

$vscale
[1] 2.28705348 0.00000001 1.00000000 1.00000000 1.00000000 1.00000000
[7] 1.00000000

Rsolnp::gosolnp函数

功能
当目标函数不光滑或者有很多极小值点时,很难判断全局最优解,并且结果往往依赖于初值的选取。该函数可以随机地设置许多初值来运行求解器,以解决前述的问题。

函数形式

gosolnp(pars = NULL,   
        fun,   
        eqfun = NULL,   
        eqB = NULL,   
        ineqfun = NULL,   
        ineqLB = NULL,   
        ineqUB = NULL,   
        LB = NULL,   
        UB = NULL,   
        control = list(),   
        ...)  

参数

  • pars: 初值。
  • fixed: 需要保持固定而不能随机生成的pars索引。
  • fun: 目标函数,返回值是一个常数。
  • eqfun: (可选) 返回等式约束的值,多个等式约束返回一个vector.
  • eqB: (可选) 等式约束的约束值,和eqfun的返回结果相一致。
  • ineqfun: (可选) 不等式约束,多个返回由表达式结果组成的向量。
  • ineqLB: (可选) 不等式约束的下界约束,维数和ineqfun返回结果一致。
  • ineqUB:(可选) 不等式约束的上界约束,维数和ineqfun返回结果一致。
  • LB: 优化变量的下界向量,该函数是不可选的!
  • UB: 优化变量的上界向量,该函数是不可选的!
  • control: (可选) 优化器的控制参数和solnp函数一致。
  • distr: 一个和优化参数相同维数的向量,表示每个元素产生时服从的分布。 1:uniform, 2: truncated nromal, 3: normal.
  • distr.opt: 当参数产生不服从uniform时,需要设置分布的参数mean和sd.
    -** n.restarts**: 重启求解器的次数。
  • n.sim: 在每一重启求解器时随机参数的个数。注意,如果不等约束存在时,往往会被拒绝。当不等约束的区间比较宽时,这个参数才会用到。
  • cluster: 如果你想用并行计算的话,这里会启动并行池,但是一定要记得关闭并行池!!
  • rseed: (可选) 用初始化随机数产生的种子,不设置时就用系统时间作为种子。
    函数返回值
  • pars: 优化变量的结果。
  • convergenc: 是否优化成功的指标,0表示优化成功,1或2表示没有优化成功。
  • values: 优化迭代过程中的目标函数值组成的向量,最后一个值是最终的优化结果。
  • lagrange: Lagrange乘子组成的向量。
  • hessian: 最优解时的hessian矩阵。
  • ineqx0: 取得最优解时的,不等式约束的松弛变量(slack variables)。
  • nfuneval: The number of function evaluations.
  • elapsed: 优化所耗费的时间。
  • start.pars: 用来重启求解器的参数。

代码

library(pacman)
p_load(Rsolnp)

fn1=function(x)
{
  z=x[1]^2*sin(x[2])+x[2]^2*cos(x[1])
  return(z)
}
eqn1=function(x){
  z1=x[1]*x[2]
}
ineqn1=function(x){
  z1=x[1]+x[2]
  z2=3*x[1]-x[2]
  z3=sin(x[1])*cos(x[2])
  return(c(z1,z2,z3))
}

res=gosolnp(pars = NULL, fixed = NULL, fun=fn1, eqfun = eqn1, eqB = c(2), 
        ineqfun = ineqn1,ineqLB = c(2,1,-Inf), ineqUB = c(Inf,3,3), LB = c(-100,-100), UB = c(100,100), 
        distr = c(1, 1), distr.opt = list(), n.restarts = 10, n.sim = 20000,cluster = NULL, 
        rseed = NULL,control = list(inner.iter=20))

结果

> res
$pars
[1] 1.403075 1.425440

$convergence
[1] 0

$values
 [1] -6100.911889  -537.602906   158.772954   246.463384   -48.645892    -7.284604     2.431605     2.317129     2.295495
[10]     2.286594     2.286974     2.287053     2.287053     2.287053

$lagrange
              [,1]
[1,]  5.424216e-01
[2,] -7.233024e-13
[3,]  2.581390e-13
[4,]  3.923091e-12

$hessian
            [,1]        [,2]        [,3]      [,4]        [,5]
[1,]  0.91620341  0.18268745 -0.08606492 -1.102502 -0.46314376
[2,]  0.18268745  0.24258774  0.09443442 -1.041828  0.02485526
[3,] -0.08606492  0.09443442  0.87122932 -2.169818 -0.89204316
[4,] -1.10250227 -1.04182752 -2.16981847 12.382395  1.35653707
[5,] -0.46314376  0.02485526 -0.89204316  1.356537  1.84242189

$ineqx0
[1] 2.8285155 2.7837853 0.1428122

$nfuneval
[1] 344

$outer.iter
[1] 13

$elapsed
Time difference of 0.05988693 secs

$vscale
[1] 2.28705348 0.00000001 1.00000000 1.00000000 1.00000000 1.00000000 1.00000000

$start.pars
[1] 28.34000 82.92805

$rseed
[1] 1580362550

alabama::auglag函数

功能
主要用于求解带有约束的优化问题。它实现了增广拉格朗日法(Augmented Lagrangian Method),适用于非线性优化问题,特别是当目标函数和约束条件都是非线性时。

用法

auglag(par, fn, gr, hin, hin.jac, heq, heq.jac, 
control.outer=list(), control.optim = list(), ...)

参数

  • par:参数值的起始向量。可以指定任何初始向量,即使这些向量违反不等式约束。
  • fn:要优化的非线性目标函数。一个标量函数。
  • gr:目标函数 fn 在参数处的梯度。它是一个向量函数,接受一个实数向量作为参数,并返回一个与输入向量长度相同的实数向量。默认值为“NULL”,这意味着梯度将通过数值方法计算。在高维问题中,如果提供确切的梯度,计算速度会显著提高。见 示例。
  • hin: 一个向量函数,指定不等式约束,使得 hin[j] > 0 对于所有 j 成立。
  • hin.jac:hin 的雅可比矩阵。如果未指定,将使用有限差分法计算,但如果指定,将会更快。
  • heq:一个向量函数,指定等式约束,使得 heq[j] = 0 对于所有 j 成立。
  • heq.jac: heq 的雅可比矩阵。如果未指定,将使用有限差分法计算,但如果指定,将会更快。
  • control.outer:一个控制参数列表,用于 constrOptim.nl 中的外部循环。有关更多信息,请参见 详细信息。
  • control.optim:一个控制参数列表,用于内部循环中的无约束优化算法。与 optim 或 nlminb 中使用的参数相同。
  • …:传递给 fn、gr、hin、heq 的附加参数。所有这些参数必须接受任何指定的参数,无论是显式地还是通过 … 参数,但它们不需要全部使用。

输出值

  • par: 优化非线性目标函数的参数,满足约束条件,如果收敛成功。
  • value:在终止时目标函数的值。
  • counts: 一个长度为2的向量,表示目标函数 fn 和梯度 gr 被评估的次数。
  • convergence: 一个整数代码,指示收敛的类型。0 表示成功收敛。正整数代码表示未能收敛。
  • outer.iterations: 外部迭代的次数。
  • lambda: 拉格朗日参数的值。它是一个与不等式和等式的总数相同长度的向量。对于非活动的不等式,它必须为零;对于活动的不等式,它必须为非负;对于等式,它可以是任何符号。
  • sigma: 增强罚参数的值,用于二次项。
  • gradient: 在收敛时增强拉格朗日函数的梯度。应该很小。
  • hessian: 在收敛时增强拉格朗日函数的海森矩阵。对于最小化(最大化),它应该是正定(负定)的。
  • ineq: 在收敛时不等式约束的值。所有值必须为非负。
  • equal: 在收敛时等式约束的值。所有值必须接近于零。
  • kkt1:一个逻辑变量,指示一阶KKT条件是否得到满足。
  • kkt2: 一个逻辑变量,指示二阶KKT条件是否得到满足。
    示例
library(alabama)  

# 目标函数  
objective_function <- function(x) {  
  return((x[1] - 1)^2 + (x[2] - 2)^2)  # 最小化的目标  
}  

# 不等式约束函数  
ineq_constraint <- function(x) {  
  return(c(x[1] + x[2]-3 ))  # g(x)-3 > 0
}  

# 等式约束函数  
eq_constraint <- function(x) {  
  return(c(x[1] - x[2]))  # h(x) = 0  
}  

# 初始值  
start <- c(0, 0)  

# 优化  
result <- auglag(par = start,   
                 fn = objective_function,   
                 hin = ineq_constraint,   
                 heq = eq_constraint
                )
                 
# 查看结果  
print(result)  
# $par
# [1] 1.5 1.5
# 
# $value
# [1] 0.5000001
# 
# $counts
# function gradient 
# 199       54 
# 
# $convergence
# [1] 0
# 
# $message
# NULL
# 
# $outer.iterations
# [1] 11
# 
# $lambda
# [1] 0.0005086384 1.0016110362
# 
# $sigma
# [1] 333333.3
# 
# $gradient
# [1]  0.01924964 -0.01904958
# 
# $ineq
# [1] 3.267964e-08
# 
# $equal
# [1] 6.228175e-08
# 
# $kkt1
# [1] FALSE
# 
# $kkt2
# [1] TRUE
# 
# $hessian
# [,1]      [,2]
# [1,]  499665.5 -167003.2
# [2,] -167003.2  499665.5

其他函数

kernlab包中的ipop函数可以求解二次规划问题,利用内点法,计算复杂度是 O ( n 3 ) O(n^3) O(n3)
osqp包中的OSQP在N特别大的时候效率会更高。

Logo

有“AI”的1024 = 2048,欢迎大家加入2048 AI社区

更多推荐