> One can do better yet by using the logspace_add(logx, logy) function
> available in the C API to R, but there is currently not an R-language
> interface
> to that (or logspace_subtract(logx,logy)). They calculate log(exp(logx) +-
> exp(logy))
> with minimal roundoff error.
>
> You could probably also use the built-in plogis() function, with its log=TRUE
> and
> perhaps lower.tail=FALSE arguments.
Here are functions, "f2" and "f3", which compute log(f()), via the above
algorithms,
along with your original "f" and my first suggestion "f1". f2 and f3 give
almost
identical results over a wide range of inputs (pos. and neg.). f1 agrees with
them
for positive inputs but falls apart for larger negative inputs. Any of f1, f2,
or f3 will
do your optimization, although I don't know why you want to optimize a monotonic
function.
f <-
function (k, T_s = 20) {
2 - 2/(1 + exp(-2*T_s*k))
}
f1 <-
function(k, T_s = 20) {
# log(f(k))
log(2) - 2*T_s*k - log1p(exp(-2*T_s*k))
}
f2 <-
function(k, T_s = 20) {
# log(f(k))
log(2) - logspace_add(2*T_s*k, 0)
}
logspace_add <-
function(logx, logy)
{
# log(exp(logx) + exp(logy))
mn <- pmin(logx, logy)
mx <- pmax(logx, logy)
mx + log1p(exp(mn - mx))
}
f3 <-
function(k, T_s = 20) {
# log(f(k))
log(2) + plogis(k, lower.tail=FALSE, log=TRUE, scale=1/(T_s*2))
}
Bill Dunlap
Spotfire, TIBCO Software
wdunlap tibco.com
> -----Original Message-----
> From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On
> Behalf
> Of William Dunlap
> Sent: Monday, October 28, 2013 3:03 PM
> To: Rolf Turner; Shantanu MULLICK
> Cc: r-help@r-project.org
> Subject: Re: [R] Optimize function in R: unable to find maximum of logistic
> function
>
> Note that f(x) returns exactly zero for all x above 1, hence its estimated
> derivative is 0 everywhere in that region. You are asking optimize to find
> the non-zero part of a function which is 0 in more than 80% of its domain.
>
> You can put a trace on f() to see where optimize() looks:
> > trace(f, quote(cat("k=", deparse(k), "\n")), print=FALSE)
> [1] "f"
> > optimize(f, c(0, 5), tol = 1e-14, maximum= TRUE)
> k= 1.90983005625053
> k= 3.09016994374947
> k= 3.81966011250105
> k= 4.27050983124842
> k= 4.54915028125263
> k= 4.72135954999579
> k= 4.82779073125683
> k= 4.89356881873896
> ... eventually homing in on 5 ...
> It looks in a few places, but f() is 0 in all of them so it figures that
> is the mininum and the maximum value it ever takes. I suppose
> you could ask that optimize look in more places for a place with
> a non-zero derivative but it would be very expensive to look in all
> (c. 2^50) places.
>
> You probably should have it maximize the log of your objective
> function, which gives you more range to play with, and you will
> have to do some numerical analysis to minimize underflow and
> roundoff error. E.g., one could use
> f1 <- function (k) {
> # log(f(k))
> T_s <- 20
> log(2) - 2 * T_s * k - log1p(exp(-2 * T_s * k))
> }
> where log1p(x) calculates log(1+x) more accurately than a computer's
> 'log' and '+' can. (Remember that 1+10^-17 exactly equals 1 here.)
>
> One can do better yet by using the logspace_add(logx, logy) function
> available in the C API to R, but there is currently not an R-language
> interface
> to that (or logspace_subtract(logx,logy)). They calculate log(exp(logx) +-
> exp(logy))
> with minimal roundoff error.
>
> You could probably also use the built-in plogis() function, with its log=TRUE
> and
> perhaps lower.tail=FALSE arguments.
>
> Bill Dunlap
> Spotfire, TIBCO Software
> wdunlap tibco.com
>
>
> > -----Original Message-----
> > From: r-help-boun...@r-project.org [mailto:r-help-boun...@r-project.org] On
> > Behalf
> > Of Rolf Turner
> > Sent: Monday, October 28, 2013 2:01 PM
> > To: Shantanu MULLICK
> > Cc: r-help@r-project.org
> > Subject: Re: [R] Optimize function in R: unable to find maximum of logistic
> > function
> >
> >
> > This could be described as a bug, perhaps. Or it could be described as
> > an indication that numerical optimization is inevitably tricky. Notice that
> > if you narrow down your search interval from [0,5] to [0,0.5] you get the
> > right answer:
> >
> > > optimize(f, c(0, 0.5), maximum= TRUE,tol=1e-10)
> > $maximum
> > [1] 4.192436e-11
> >
> > $objective
> > [1] 1
> >
> > I guess there's a problem with finding a gradient that is effectively
> > (numerically) zero when "k" is equal to 5.
> >
> >
> > cheers,
> >
> > Rolf Turner
> >
> >
> > On 10/29/13 06:00, Shantanu MULLICK wrote:
> > > Hello Everyone,
> > >
> > > I want to perform a 1-D optimization by using the optimize() function. I
> > > want to find the maximum value of a "logistic" function. The optimize()
> > > function gives the wrong result.
> > >
> > > My code:
> > > f= function (k) {
> > > T_s = 20
> > > result = (2- 2/(1+ exp(-2*T_s*k)))
> > > return(result)
> > > }
> > > optimize(f, c(0, 5), tol = 0.0000000000001, maximum= TRUE)
> > >
> > > The maximum value for the function happens at k=0, and the maximum value
> > > is
> > > 1. Yet the optimise function, says that the maximum value happens at k=
> > > 4.9995, and the maximum value is 0.
> > >
> >
> > ______________________________________________
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained, reproducible code.
>
> ______________________________________________
> R-help@r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
