R Optimization with Equality and Inequality Constraints

R optimization with equality and inequality constraints

On this occasion optim will not work obviously because you have equality constraints. constrOptim will not work either for the same reason (I tried converting the equality to two inequalities i.e. greater and less than 15 but this didn't work with constrOptim).

However, there is a package dedicated to this kind of problem and that is Rsolnp.

You use it the following way:

#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)

Output:

> 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

So the resulting optimal values are:

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

which means that the first parameter is 15 and the rest zero and zero. This is indeed the global minimum in your function since the x2 is adding to the function and 5 * x1 has a much greater (negative) influence than x3 on the outcome. The choice of 15, 0, 0 is the solution and the global minimum to the function according to the constraints.

The function worked great!

Multi-objective optimization with both equality and inequality constraints in R

Use normalization so you can drop the equality constraint:

Objectives:

 w <- w / sum(w)
 f1<-(w1)^2+4*(w2)^4+3*(w3)^2+5*w4+w5
 f2<-((w1)-(w2)^2+3(w3)^2-4*(w4)^3-w5)^2-8

Constraints:

 0 <= w1,w2,w3,w4,w5 <= 1

quadratic optimization in R with both equality and inequality constraints

Here's a working example:

library('quadprog')

# min  
#     -8 x1 -16 x2 + x1^2 + 4 x2^2 
#
# s.t.
#
# x1 + 2 x2 == 12 # equalities 
# x1 +   x2 <= 10 # inequalities (N.B. you need to turn it into "greater-equal" form )
# 1 <= x1 <= 3 # bounds 
# 1 <= x2 <= 6 # bounds 

H <- rbind(c(2, 0),
           c(0, 8))

f <- c(8,16)

# equalities
A.eq <- rbind(c(1,2))
b.eq <- c(12)

# inequalities
A.ge <- rbind(c(-1,-1))
b.ge <- c(-10)

# lower-bounds 
A.lbs <- rbind(c( 1, 0),
               c( 0, 1))
b.lbs <- c(1, 1)

# upper-bounds on variables
A.ubs <- rbind(c(-1, 0),
               c( 0,-1))
b.ubs <- c(-3, -6)

# solve
sol <- solve.QP(Dmat = H,
                dvec = f,
                Amat = t(rbind(A.eq, A.ge, A.lbs, A.ubs)),
                bvec = c(b.eq, b.ge, b.lbs, b.ubs),
                meq = 1) # this argument says the first "meq" rows of Amat are equalities

sol
> sol
$solution
[1] 3.0 4.5

$value
[1] -6

$unconstrained.solution
[1] 4 2

$iterations
[1] 3 0

$Lagrangian
[1] 10  0  0  0 12  0

$iact
[1] 1 5

optim with inequality constraint

You are trying to solve an equation. Whether or not the constraint is met, can only be decided ex post.
You can use uniroot as follows

f <- function(x,sigma=1,Q=1) {x/sigma^2 - log(x/sigma^2) - 1 - 1/Q}
uniroot(f,c(1,5))

giving

$root
[1] 3.146198

$f.root
[1] 3.552369e-06

$iter
[1] 5

$estim.prec
[1] 6.103516e-05

constrained optimization R: another example

constrOptim() uses linear inequality constraints and defines the feasible region by ui %*% param - ci >= 0. If the constraint is 3 * d1 + d2 <= 280, ui is c(-3, -1) and ci is -280.

constrOptim(); inequality constraint is: 3 * d1 + d2 <= 280

Fd <- function(betas) {
    b1 = betas[1]
    b2 = betas[2]
   (224 * b1 + 84 * b2 + b1 * b2 - 2 * b1^2 - b2^2)
}

theta = c(59.999,100)    # because of needing " ui %*% inital_par - ci > 0 "
ui = c(-3, -1)
ci = -280                # those ui & ci mean " -3*par[1] + -1*par[2] + 280 >= 0 "

constrOptim(theta, Fd, NULL, ui = ui, ci = ci, control=list(fnscale=-1))
  # $par
  # [1] 69.00002 72.99993

[Edited]

If you want not inequality but equality constraints, it would be better to use Rsolnp or alabama package. They can use inequality and/or equality constraints (see Constrained Optimization library for equality and inequality constraints).

solnp(); auglag(); equality constraint is: 3 * d1 + d2 = 280

library(Rsolnp); library(alabama); 

Fd2 <- function(betas) {     #  -1 * Fd
   b1 = betas[1]
   b2 = betas[2]
   -1 * (224 * b1 + 84 * b2 + b1 * b2 - 2 * b1^2 - b2^2)
}

eqFd <- function(betas) {  # the equality constraint
    b1 = betas[1]
    b2 = betas[2]
    (3 * b1 + b2 -280)
}

solnp(pars = c(60, 100), fun = Fd2, eqfun = eqFd, eqB = 0)
auglag(par = c(60, 100), fn = Fd2, heq = eqFd)

---

Related Topics