Non-linear optimisation/programming with integer variables in R

Question

I wonder if anyone is able to suggest some packages to solve a non-linear optimisation problem which can provide integer variables for an optimum solution? The problem is to minimise a function with an equality constraint subject to some lower and upper boundary constraints.

I've used the 'nloptr' package in R for a non-linear optimisation problem which worked nicely but would now like to extend the method to have some of the variables as integers. From my use and understanding of nloptr so far, it can only return continuous, and not integer variables for an optimum solution.

I believe this sort of problem needs to be solved using mixed-integer non-linear programming.

One example of the problem in a form for nloptr:

min f(x) (x-y)^2/y + (p-q)^2/q
so that (x-y)^2/y + (p-q)^2/q = 10.2

where 
x and p are positive integers not equal to 0 
and 
y and q may or may not be positive integers not equal to 0

The nloptr code for this in R would look like this

library('nloptr')

x1 <- c(50,25,20,15)

fn <- function(x) {
  (((x[1] - x[2])^2)/x[2]) + (((x[3] - x[4])^2)/x[4])
  }

heq <- function(x) {
  fn(x)-10.2
}

lower_limit <- c(0,0,0,0)
upper_limit <- c(67.314, 78, 76.11, 86)


slsqp(x1, fn, lower = lower_limit, upper = upper_limit,  hin = NULL, heq = heq, control = list(xtol_rel = 1e-8, check_derivatives = FALSE))

This would output the following:

$par
[1] 46.74823 29.72770 18.93794 16.22137

$value
[1] 10.2

$iter
[1] 6

$convergence
[1] 4

$message
[1] "NLOPT_XTOL_REACHED: Optimization stopped because xtol_rel or xtol_abs (above) was reached."

This is the sort of I result I am looking for but as stated above, I need x and p as integers.

I've had a look at https://cran.r-project.org/web/views/Optimization.html which has a really good list of packages for mixed-integer non-linear programming but just wondered if anyone had experience with any of them and what they think might be most appropriate to solve the problem as stated above.

There was a similar question about this posted around 7 years ago on here but it ended up with a link to the cran page so thought it would be worth asking again.

Thanks very much for your input.

Cheers,

Andrew

Answer 1

Here's an example of how it doesn't work with CVXR, without a simpler objective function. I suspect the problem is not convex, even with the constraints, so an alternative option is required.

#base example from https://cvxr.rbind.io/cvxr_examples/cvxr_gentle-intro/
#install.packages("CVXR")
library(CVXR)

#modified for Stackoverflow integer MIQP ####
#Solves, but terms not normalised by y and q respectively

# Variables minimized over
x <- Variable(1, integer=TRUE)
y <- Variable(1)
p <- Variable(1, integer=TRUE)
q <- Variable(1)

# Problem definition (terms not normalised by y and q respectively)
objective <- Minimize((x - y)^2 + (p - q)^2)
constraints <- list(x >= 0, y >= 0, p >= 0, q >= 0, 
                    x <= 67.314, y <= 78, p <= 76.11, q <= 86)
prob2.1 <- Problem(objective, constraints)

# Problem solution
solution2.1 <- solve(prob2.1)
solution2.1$status
solution2.1$value
solution2.1$getValue(x)
solution2.1$getValue(y)
solution2.1$getValue(p)
solution2.1$getValue(q)


#modified for Stackoverflow integer NLP (not integer) ####
#Does not solve, not convex?

# Variables minimized over
x <- Variable(1)
y <- Variable(1)
p <- Variable(1)
q <- Variable(1)

# Problem definition
objective <- Minimize((x - y)^2/y + (p - q)^2/q)
constraints <- list(x >= 0, y >= 0, p >= 0, q >= 0, 
                    x <= 67.314, y <= 78, p <= 76.11, q <= 86)
prob2.1 <- Problem(objective, constraints)

# Problem solution
solution2.1 <- solve(prob2.1, gp = TRUE)
solution2.1 <- solve(prob2.1, gp = FALSE)
# solution2.1$status
# solution2.1$value
# solution2.1$getValue(x)
# solution2.1$getValue(y)
# solution2.1$getValue(p)
# solution2.1$getValue(q)


#modified for Stackoverflow integer MINLP ####
#Does not solve

# Variables minimized over
x <- Variable(1, integer=TRUE)
y <- Variable(1)
p <- Variable(1, integer=TRUE)
q <- Variable(1)

# Problem definition
objective <- Minimize((x - y)^2/y + (p - q)^2/q)
constraints <- list(x >= 0, y >= 0, p >= 0, q >= 0, 
                    x <= 67.314, y <= 78, p <= 76.11, q <= 86)
prob2.1 <- Problem(objective, constraints)

# Problem solution
solution2.1 <- solve(prob2.1, gp = TRUE)
solution2.1 <- solve(prob2.1, gp = FALSE)
# solution2.1$status
# solution2.1$value
# solution2.1$getValue(x)
# solution2.1$getValue(y)
# solution2.1$getValue(p)
# solution2.1$getValue(q)


#modified for Stackoverflow integer NLP (not integer) ####
#objective multiplied by y*q, Does not solve, not convex?

# Variables minimized over
x <- Variable(1)
y <- Variable(1)
p <- Variable(1)
q <- Variable(1)

# Problem definition
objective <- Minimize((x - y)^2*q + (p - q)^2*y)
constraints <- list(x >= 0, y >= 0, p >= 0, q >= 0, 
                    x <= 67.314, y <= 78, p <= 76.11, q <= 86)
prob2.1 <- Problem(objective, constraints)

# Problem solution
solution2.1 <- solve(prob2.1, gp = TRUE)
solution2.1 <- solve(prob2.1, gp = FALSE)
# solution2.1$status
# solution2.1$value
# solution2.1$getValue(x)
# solution2.1$getValue(y)
# solution2.1$getValue(p)
# solution2.1$getValue(q)