R оптимизация с ограничениями равенства и неравенства

Question

Я пытаюсь найти локальный минимум функции, а параметры имеют фиксированную сумму. Например,

Рх = 10 - 5x1 + 2x2 - x3

и условия таковы,

х1 + х2 + х3 = 15

(x1, x2, x3)> = 0

Где сумма x1, x2 и x3 имеет известное значение и все они больше нуля. В R, это будет выглядеть примерно так,

Fx = function(x) {10 - (5*x[1] + 2*x[2] + x[3])} 
opt = optim(c(1,1,1), Fx, method = "L-BFGS-B", lower=c(0,0,0), upper=c(15,15,15)) 

Я также пытался использовать неравенство с constrOptim, чтобы заставить сумму быть исправлено. Я все еще думаю, что это может быть правдоподобной работой, но я не смог заставить ее работать. Это упрощенный пример реальной проблемы, но любая помощь будет очень оценена.

Answer 1

В этом случае optim не будет работать, очевидно, потому что у вас есть ограничения на равенство. constrOptim не будет работать по той же причине (я попытался преобразовать равенство в два неравенства, то есть больше и меньше 15, но это не с constrOptim).

Однако есть пакет, посвященный этой проблеме, и это Rsolnp.

Вы можете использовать его следующим образом:

#specify your function 
opt_func <- function(x) { 
    10 - 5*x[1] + 2 * x[2] - x[3] 
} 

#specify the equality function. The number 15 (to which the function is equal) 
#is specified as an additional argument 
equal <- function(x) { 
    x[1] + x[2] + x[3] 
} 

#the optimiser - minimises by default 
solnp(c(5,5,5), #starting values (random - obviously need to be positive and sum to 15) 
     opt_func, #function to optimise 
     eqfun=equal, #equality function 
     eqB=15, #the equality constraint 
     LB=c(0,0,0), #lower bound for parameters i.e. greater than zero 
     UB=c(100,100,100)) #upper bound for parameters (I just chose 100 randomly) 

Выход:

> solnp(c(5,5,5), 
+  opt_func, 
+  eqfun=equal, 
+  eqB=15, 
+  LB=c(0,0,0), 
+  UB=c(100,100,100)) 

Iter: 1 fn: -65.0000  Pars: 14.99999993134 0.00000002235 0.00000004632 
Iter: 2 fn: -65.0000  Pars: 14.999999973563 0.000000005745 0.000000020692 
solnp--> Completed in 2 iterations 
$pars 
[1] 1.500000e+01 5.745236e-09 2.069192e-08 

$convergence 
[1] 0 

$values 
[1] -10 -65 -65 

$lagrange 
    [,1] 
[1,] -5 

$hessian 
      [,1]  [,2]  [,3] 
[1,] 121313076 121313076 121313076 
[2,] 121313076 121313076 121313076 
[3,] 121313076 121313076 121313076 

$ineqx0 
NULL 

$nfuneval 
[1] 126 

$outer.iter 
[1] 2 

$elapsed 
Time difference of 0.1770101 secs 

$vscale 
[1] 6.5e+01 1.0e-08 1.0e+00 1.0e+00 1.0e+00 

Таким образом, полученные оптимальные значения:

$pars 
[1] 1.500000e+01 5.745236e-09 2.069192e-08 

, что означает, что первый параметр 15 и остальные ноль и ноль. Это действительно глобальный минимум в вашей функции, поскольку x2 добавляет функцию, а 5 * x1 имеет гораздо большее (отрицательное) влияние, чем x3 на результат. Выбор 15, 0, 0 - это решение и глобальный минимум функции в соответствии с ограничениями.

Функция работала отлично!

Answer 2

Это на самом деле проблема линейного программирования, поэтому естественным подходом было бы использовать линейный программный решатель, такой как пакет lpSolve. Вы должны обеспечить целевую функцию и матрицу ограничений и решатель сделает все остальное:

library(lpSolve) 
mod <- lp("min", c(-5, 2, -1), matrix(c(1, 1, 1), nrow=1), "=", 15) 

Тогда вы можете получить доступ оптимального решения и объективное значения (добавление постоянного члена 10, которая не предусмотрена в решатель):

mod$solution 
# [1] 15 0 0 
mod$objval + 10 
# [1] -65 

линейное программирование решатель должен быть много быстрее, чем в общей нелинейной оптимизации решателя и не должно быть проблем возвращающегося точное оптимальное решение (вместо соседней точки, которые могут быть подвержены ошибкам округления).