Reputation: 712
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
Upvotes: 1
Views: 3319
Reputation: 44330
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.
Upvotes: 5