Extreme scaling quite often ruins optimization calculations. If you think
available methods
are capable of doing this, there's a bridge I can sell you in NYC.
I've been trying for some years to develop a good check on scaling so I can
tell users
who provide functions like this to send (lots of) money and I'll give them the
best answer
there is (generally no answer at all). Or, more seriously, to inform them that
they should
not expect results unless they scale. Richard Varga once said some decades ago
that any
problem was trivially solvable in the right scale, and he was mostly right.
Scaling is
important.
To see the range of answers from a number of methods, the script below is
helpful. I had
to remove lbfgsb3c from the mix as it stopped mid-calculation in unrecoverable
way. Note
that I use my development version of optimx, so some methods might not be
included in
CRAN offering. Just remove the methods from the ameth and bmeth lists if
necessary.
Cheers, John Nash
# CErickson221223.R
# optim(c(0,0), function(x) {x[1]*1e-317}, lower=c(-1,-1), upper=c(1,1),
# method='L-BFGS-B')
tfun <- function(x, xpnt=317){
if ((length(x)) != 2) {stop("Must have length 2")}
scl <- 10^(-xpnt)
val <- x[1]*scl # note that x[2] unused. May be an issue!
val
}
gtfun <- function(x, xpnt=317){ # gradient
scl <- 10^(-xpnt)
gg <- c(scl, 0.0)
gg
}
xx <- c(0,0)
lo <- c(-1,-1)
up <- c(1,1)
print(tfun(xx))
library(optimx)
ameth <- c("BFGS", "CG", "Nelder-Mead", "L-BFGS-B", "nlm", "nlminb",
"Rcgmin", "Rtnmin", "Rvmmin", "spg", "ucminf", "newuoa", "bobyqa",
"nmkb", "hjkb", "hjn", "lbfgs", "subplex", "ncg", "nvm", "mla",
"slsqp", "anms")
bmeth <- c("L-BFGS-B", "nlminb", "Rcgmin", "Rtnmin", "nvm",
"bobyqa", "nmkb", "hjkb", "hjn", "ncg", "slsqp")
tstu <- opm(x<-c(0,0), fn=tfun, gr=gtfun, method=ameth, control=list(trace=0))
summary(tstu, order=value)
tstb <- opm(x<-c(0,0), fn=tfun, gr=gtfun, method=bmeth, lower=lo, upper=up,
control=list(trace=0))
summary(tstb, order=value)
On 2022-12-23 13:41, Rui Barradas wrote:
Às 17:30 de 23/12/2022, Collin Erickson escreveu:
Hello,
I've come across what seems to be a bug in optim that has become a nuisance
for me.
To recreate the bug, run:
optim(c(0,0), function(x) {x[1]*1e-317}, lower=c(-1,-1), upper=c(1,1),
method='L-BFGS-B')
The error message says:
Error in optim(c(0, 0), function(x) { :
non-finite value supplied by optim
What makes this particularly treacherous is that this error only occurs for
specific powers. By running the following code you will find that the error
only occurs when the power is between -309 and -320; above and below that
work fine.
p <- 1:1000
giveserror <- rep(NA, length(p))
for (i in seq_along(p)) {
tryout <- try({
optim(c(0,0), function(x) {x[1]*10^-p[i]}, lower=c(-1,-1),
upper=c(1,1), method='L-BFGS-B')
})
giveserror[i] <- inherits(tryout, "try-error")
}
p[giveserror]
Obviously my function is much more complex than this and usually doesn't
fail, but this reprex demonstrates that this is a problem. To avoid the
error I may multiply by a factor or take the log, but it seems like a
legitimate bug that should be fixed.
I tried to look inside of optim to track down the error, but the error lies
within the external C code:
.External2(C_optim, par, fn1, gr1, method, con, lower,
upper)
For reference, I am running R 4.2.2, but was also able to recreate this bug
on another device running R 4.1.2 and another running 4.0.3.
Thanks,
Collin Erickson
[[alternative HTML version deleted]]
______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
Hello,
See if this R-Help thread [1] earlier this year is relevant.
In particular, the post by R Core team member Luke Tierney [2], that answers
much better than what I could.
The very small numbers in your question seem to have hit a limit and this limit
is not R related.
[1] https://stat.ethz.ch/pipermail/r-help/2022-February/473840.html
[2] https://stat.ethz.ch/pipermail/r-help/2022-February/473844.html
Hope this helps,
Rui Barradas
______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
______________________________________________
R-devel@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel
