How to deal with optim when some parameters become known?

Question

I have faced a problem with passing arguments to optim. Suppose I want to do box constraint minimization on a multivariate function, for example

fr <- function(x) {   ## Rosenbrock  function
  x1 <- x[1]
  x2 <- x[2]
  x3 <- x[3]
  x4 <- x[4]
  100 * (x2 - x1 * x1)^2 + (1 - x1)^2 +
  100 * (x3 - x2 * x2)^2 + (1 - x2)^2 +
  100 * (x4 - x3 * x3)^2 + (1 - x3)^2
}

As usual optim can be used as following:

optim(par = c(0, 1, 1, 2), fr, method = "L-BFGS-B", lower = c(0, 0, 0, 0), upper = c(3, 3, 3, 3))

Now, suppose this procedure repeated in an algorithm which changes lower and upper (box constraints), followed by par, such that in some iterations one, two or three value of parameters become known, for example x1 = 1. in this case I expect optim to handle this by setting the initial value, lower and upper bounds of x1 to 1:

optim(par = c(1, 1, 1, 2), fr, method = "L-BFGS-B", lower = c(1, 0, 0, 0), upper =    c(1, 3, 3, 3))

But by runnig this line I got an error:

Error in optim(par = c(1, 1, 1, 2), fr, method = "L-BFGS-B", lower = c(1,  : non-finite finite-difference value [1]

Now, the question is how can I deal with this feature of optim without defining many new functions when one or some of the parameters become known?

Thank you in advance

Answer 1

It sounds like optim is not able to handle the upper and lower matching. I suppose you could parameterize your function with the known values and use some simple ifelse statements to check if you should be using the passed value from optim or the known value:

# Slightly redefined function to optimize
fr2 <- function(opt.x, known.x) {
  x <- ifelse(is.na(known.x), opt.x, known.x)
  100 * (x[2] - x[1] * x[1])^2 + (1 - x[1])^2 +
  100 * (x[3] - x[2] * x[2])^2 + (1 - x[2])^2 +
  100 * (x[4] - x[3] * x[3])^2 + (1 - x[3])^2
}

# Optimize, and then replace the appropriate indices of the result with known vals
known.x <- c(NA, 1, NA, 1)
opt.result <- optim(par = c(0, 1, 1, 2), fr2, method = "L-BFGS-B",
                    lower = c(0, 0, 0, 0), upper = c(3, 3, 3, 3), known.x=known.x)
opt.result$par <- ifelse(is.na(known.x), opt.result$par, known.x)
opt.result
# $par
# [1] 0.9999995 1.0000000 0.9999996 1.0000000
# 
# $value
# [1] 1.795791e-10
# 
# $counts
# function gradient 
#       13       13 
# 
# $convergence
# [1] 0
# 
# $message
# [1] "CONVERGENCE: REL_REDUCTION_OF_F <= FACTR*EPSMCH"

This code basically ignores the indices passed from optim if they are already known, and just uses the known values in those cases.