如何使用 nloptr 在 R 中构造一个具有 n 项优化的 objective 函数？
How to construct an objective function with n terms for optimisation in R using nloptr?
我已经更新了这个问题以更清楚地说明问题。
我有以下形式的数据集：
#sample dataset
dt <- read.table(text = 'Plant Type Ownership Capacity VC TargetPLF
1 T Base Pub 1000 1.50 0.80
2 S Base Pub 500 1.62 0.80
3 Y Base Pub 500 2.43 0.75
4 K Base Pub 500 1.93 0.80
5 J Base Pub 500 3.40 0.80', header =T)
demand <- c(1500,2000,2200,3000,2800)
我的 objective 函数如下所示：
eval_f <- function(x) #objective function
{
return((1.50*1000*x[1] + 1.62*500*x[2] + 2.43*500*x[3] + 1.93*500*x[4] + 3.40*500*x[5]) / (1000*x[1] + 500*x[2] + 500*x[3] + 500*x[4] + 500*x[5]))
}
其他参数：
x0 <- rep(0.45,5) #initialisation values
local_opts <- list( "algorithm" = "NLOPT_LD_MMA", "xtol_rel" = 1.0e-15 ) #options for optimisation
opts = list( "algorithm"= "NLOPT_GN_ISRES","xtol_rel"= 1.0e-15,"maxeval"= 160000,"local_opts" = local_opts,"print_level" = 0 ) #optimisation algorithm
# constraint function
eval_g0_eq <- function(x) {
return(1000*x[1] + 500*x[2] + 500*x[3] + 500*x[4] + 500*x[5] - 1500)
}
lb = c(rep(0.40,5)) #lower bounds of controls
ub = dt$Target #upper bounds of controls
nloptr(x0,eval_f,lb=lb,ub=ub,eval_g_eq = eval_g0_eq,opts=opts) #optimisation
如何进行 5 次优化，每次使用不同的约束函数？
所以在第二次迭代中，约束函数应该从向量 'demand' 中取出第二个元素，变成这样：
# constraint function for second iteration
eval_g0_eq <- function(x) {
return(1000*x[1] + 500*x[2] + 500*x[3] + 500*x[4] + 500*x[5] - 2000)
}
nloptr(x0,eval_f,lb=lb,ub=ub,eval_g_eq = eval_g0_eq,opts=opts) #optimisation
此外，有没有一种方法可以使用 dt 中的值构造 eval_f 和 eval_g0_eq，而不是像现在这样手动完成？
我不知道这个 nlopt 函数或库来自哪里，但我认为您不需要重复创建 objective function 或 constraint function。您可以使用 R 的矢量化计算潜力来执行此操作。
dt <- read.table(text = 'Plant Type Ownership Capacity VC TargetPLF
1 T Base Pub 1000 1.50 0.80
2 S Base Pub 500 1.62 0.80
3 Y Base Pub 500 2.43 0.75
4 K Base Pub 500 1.93 0.80
5 J Base Pub 500 3.40 0.80', header =T)
demand <- c(1500,2000,2200,3000,2800)
library(nloptr)
dt
#> Plant Type Ownership Capacity VC TargetPLF
#> 1 T Base Pub 1000 1.50 0.80
#> 2 S Base Pub 500 1.62 0.80
#> 3 Y Base Pub 500 2.43 0.75
#> 4 K Base Pub 500 1.93 0.80
#> 5 J Base Pub 500 3.40 0.80
下面的代码将执行等于 dt
中的行数的迭代
lapply(seq(nrow(dt)), function(.x) nloptr(x0 = rep(0.45,5),
eval_f = function(x) sum(dt$VC * dt$Capacity * x)/sum(dt$Capacity * x),
lb=c(rep(0.40,5)),
ub = dt$Target,
eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x] ,
opts = list( "algorithm"= "NLOPT_GN_ISRES",
"xtol_rel"= 1.0e-15,
"maxeval"= 160000,
"local_opts" = list("algorithm" = "NLOPT_LD_MMA",
"xtol_rel" = 1.0e-15 ),
"print_level" = 0 )
)
)
它的输出将是
#> [[1]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.9506666666645
#> Current value of controls: 0.7 0.4 0.4 0.4 0.4
#>
#>
#>
#> [[2]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.893
#> Current value of controls: 0.8 0.8 0.4 0.8 0.4
#>
#>
#>
#> [[3]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.92181376479112
#> Current value of controls: 0.7883007 0.7815353 0.583928 0.7845948 0.4043761
#>
#>
#>
#> [[4]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.93672658072513
#> Current value of controls: 0.7870349 0.7534322 0.6660887 0.770566 0.4118776
#>
#>
#>
#> [[5]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.93152444911471
#> Current value of controls: 0.795248 0.7687759 0.6712354 0.7692782 0.4034175
我已经更新了这个问题以更清楚地说明问题。
我有以下形式的数据集：
#sample dataset
dt <- read.table(text = 'Plant Type Ownership Capacity VC TargetPLF
1 T Base Pub 1000 1.50 0.80
2 S Base Pub 500 1.62 0.80
3 Y Base Pub 500 2.43 0.75
4 K Base Pub 500 1.93 0.80
5 J Base Pub 500 3.40 0.80', header =T)
demand <- c(1500,2000,2200,3000,2800)
我的 objective 函数如下所示：
eval_f <- function(x) #objective function
{
return((1.50*1000*x[1] + 1.62*500*x[2] + 2.43*500*x[3] + 1.93*500*x[4] + 3.40*500*x[5]) / (1000*x[1] + 500*x[2] + 500*x[3] + 500*x[4] + 500*x[5]))
}
其他参数：
x0 <- rep(0.45,5) #initialisation values
local_opts <- list( "algorithm" = "NLOPT_LD_MMA", "xtol_rel" = 1.0e-15 ) #options for optimisation
opts = list( "algorithm"= "NLOPT_GN_ISRES","xtol_rel"= 1.0e-15,"maxeval"= 160000,"local_opts" = local_opts,"print_level" = 0 ) #optimisation algorithm
# constraint function
eval_g0_eq <- function(x) {
return(1000*x[1] + 500*x[2] + 500*x[3] + 500*x[4] + 500*x[5] - 1500)
}
lb = c(rep(0.40,5)) #lower bounds of controls
ub = dt$Target #upper bounds of controls
nloptr(x0,eval_f,lb=lb,ub=ub,eval_g_eq = eval_g0_eq,opts=opts) #optimisation
如何进行 5 次优化，每次使用不同的约束函数？
所以在第二次迭代中，约束函数应该从向量 'demand' 中取出第二个元素，变成这样：
# constraint function for second iteration
eval_g0_eq <- function(x) {
return(1000*x[1] + 500*x[2] + 500*x[3] + 500*x[4] + 500*x[5] - 2000)
}
nloptr(x0,eval_f,lb=lb,ub=ub,eval_g_eq = eval_g0_eq,opts=opts) #optimisation
此外，有没有一种方法可以使用 dt 中的值构造 eval_f 和 eval_g0_eq，而不是像现在这样手动完成？
我不知道这个 nlopt 函数或库来自哪里，但我认为您不需要重复创建 objective function 或 constraint function。您可以使用 R 的矢量化计算潜力来执行此操作。
dt <- read.table(text = 'Plant Type Ownership Capacity VC TargetPLF
1 T Base Pub 1000 1.50 0.80
2 S Base Pub 500 1.62 0.80
3 Y Base Pub 500 2.43 0.75
4 K Base Pub 500 1.93 0.80
5 J Base Pub 500 3.40 0.80', header =T)
demand <- c(1500,2000,2200,3000,2800)
library(nloptr)
dt
#> Plant Type Ownership Capacity VC TargetPLF
#> 1 T Base Pub 1000 1.50 0.80
#> 2 S Base Pub 500 1.62 0.80
#> 3 Y Base Pub 500 2.43 0.75
#> 4 K Base Pub 500 1.93 0.80
#> 5 J Base Pub 500 3.40 0.80
下面的代码将执行等于 dt
lapply(seq(nrow(dt)), function(.x) nloptr(x0 = rep(0.45,5),
eval_f = function(x) sum(dt$VC * dt$Capacity * x)/sum(dt$Capacity * x),
lb=c(rep(0.40,5)),
ub = dt$Target,
eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x] ,
opts = list( "algorithm"= "NLOPT_GN_ISRES",
"xtol_rel"= 1.0e-15,
"maxeval"= 160000,
"local_opts" = list("algorithm" = "NLOPT_LD_MMA",
"xtol_rel" = 1.0e-15 ),
"print_level" = 0 )
)
)
它的输出将是
#> [[1]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.9506666666645
#> Current value of controls: 0.7 0.4 0.4 0.4 0.4
#>
#>
#>
#> [[2]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.893
#> Current value of controls: 0.8 0.8 0.4 0.8 0.4
#>
#>
#>
#> [[3]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.92181376479112
#> Current value of controls: 0.7883007 0.7815353 0.583928 0.7845948 0.4043761
#>
#>
#>
#> [[4]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.93672658072513
#> Current value of controls: 0.7870349 0.7534322 0.6660887 0.770566 0.4118776
#>
#>
#>
#> [[5]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.93152444911471
#> Current value of controls: 0.795248 0.7687759 0.6712354 0.7692782 0.4034175