[Rd] trace in uniroot() ?

[William Dunlap]

To record the value of the function as well as the arguments, you can use
the following

instrumentObjectiveFunction <- function(FUN) {
    newFUN <- local({
        INFO <- list()
        function(...) {
            value <- FUN(...)
            INFO[[length(INFO)+1]] <<- list(args=list(...), value=value)
            value
        }
    })
    newFUN
}

E.g.,
> untrace(ff)
> ff0 <- uniroot(instrumentedFF <- instrumentObjectiveFunction(ff), c(0,
10))
> str(environment(instrumentedFF)$INFO)
List of 13
 $ :List of 2
  ..$ args :List of 1
  .. ..$ : num 0
  ..$ value: num -1
 $ :List of 2
  ..$ args :List of 1
  .. ..$ : num 10
  ..$ value: num 146
 $ :List of 2
  ..$ args :List of 1
  .. ..$ : num 0.0678
  ..$ value: num -0.965
 $ :List of 2
  ..$ args :List of 1
  .. ..$ : num 5.03
  ..$ value: num 10.4
 $ :List of 2
  ..$ args :List of 1
  .. ..$ : num 0.49
  ..$ value: num -0.722
...


Bill Dunlap
TIBCO Software
wdunlap tibco.com

On Mon, Aug 13, 2018 at 3:44 PM, J C Nash <profjcnash using gmail.com> wrote:

> Despite my years with R, I didn't know about trace(). Thanks.
>
> However, my decades in the minimization and root finding game make me like
> having
> a trace that gives some info on the operation, the argument and the
> current function value.
> I've usually found glitches are a result of things like >= rather than >
> in tests etc., and
> knowing what was done is the quickest way to get there.
>
> This is, of course, the numerical software developer view. I know "users"
> (a far too vague
> term) don't like such output. I've sometimes been tempted with my svd or
> optimization codes to
> have a return message in bold-caps "YOUR ANSWER IS WRONG AND THERE'S A
> LAWYER WAITING TO
> MAKE YOU PAY", but I usually just satisfy myself with "Not at a
> minimum/root".
>
> Best, JN
>
> On 2018-08-13 06:00 PM, William Dunlap wrote:
> > I tend to avoid the the trace/verbose arguments for the various root
> finders and optimizers and instead use the trace
> > function or otherwise modify the function handed to the operator.  You
> can print or plot the arguments or save them.  E.g.,
> >
> >> trace(ff, print=FALSE, quote(cat("x=", deparse(x), "\n", sep="")))
> > [1] "ff"
> >> ff0 <- uniroot(ff, c(0, 10))
> > x=0
> > x=10
> > x=0.0678365490630423
> > x=5.03391827453152
> > x=0.490045026724842
> > x=2.76198165062818
> > x=1.09760394309444
> > x=1.92979279686131
> > x=1.34802524899502
> > x=1.38677998493585
> > x=1.3862897003949
> > x=1.38635073555115
> > x=1.3862897003949
> >
> > or
> >
> >> X <- numeric()
> >> trace(ff, print=FALSE, quote(X[[length(X)+1]] <<- x))
> > [1] "ff"
> >> ff0 <- uniroot(ff, c(0, 10))
> >> X
> >  [1]  0.00000000 10.00000000  0.06783655
> >  [4]  5.03391827  0.49004503  2.76198165
> >  [7]  1.09760394  1.92979280  1.34802525
> > [10]  1.38677998  1.38628970  1.38635074
> > [13]  1.38628970
> >
> > This will not tell you why the objective function is being called (e.g.
> in a line search
> > or in derivative estimation), but some plotting or other postprocessing
> can ususally figure that out.
> >
> >
> > Bill Dunlap
> > TIBCO Software
> > wdunlap tibco.com <http://tibco.com>
> >
> > On Mon, Jul 30, 2018 at 11:35 AM, J C Nash <profjcnash using gmail.com
> <mailto:profjcnash using gmail.com>> wrote:
> >
> >     In looking at rootfinding for the histoRicalg project (see
> gitlab.com/nashjc/histoRicalg
> >     <http://gitlab.com/nashjc/histoRicalg>),
> >     I thought I would check how uniroot() solves some problems. The
> following short example
> >
> >     ff <- function(x){ exp(0.5*x) - 2 }
> >     ff(2)
> >     ff(1)
> >     uniroot(ff, 0, 10)
> >     uniroot(ff, c(0, 10), trace=1)
> >     uniroot(ff, c(0, 10), trace=TRUE)
> >
> >
> >     shows that the trace parameter, as described in the Rd file, does
> not seem to
> >     be functional except in limited situations (and it suggests an
> >     integer, then uses a logical for the example, e.g.,
> >      ## numerically, f(-|M|) becomes zero :
> >          u3 <- uniroot(exp, c(0,2), extendInt="yes", trace=TRUE)
> >     )
> >
> >     When extendInt is set, then there is some information output, but
> trace alone
> >     produces nothing.
> >
> >     I looked at the source code -- it is in R-3.5.1/src/library/stats/R/nlm.R
> and
> >     calls zeroin2 code from R-3.5.1/src/library/stats/src/optimize.c as
> far as I
> >     can determing. My code inspection suggests trace does not show the
> iterations
> >     of the rootfinding, and only has effect when the search interval is
> allowed
> >     to be extended. It does not appear that there is any mechanism to ask
> >     the zeroin2 C code to display intermediate work.
> >
> >     This isn't desperately important for me as I wrote an R version of
> the code in
> >     package rootoned on R-forge (which Martin Maechler adapted as
> unirootR.R in
> >     Rmpfr so multi-precision roots can be found). My zeroin.R has
> 'trace' to get
> >     the pattern of different steps. In fact it is a bit excessive. Note
> >     unirootR.R uses 'verbose' rather than 'trace'. However, it would be
> nice to be
> >     able to see what is going on with uniroot() to verify equivalent
> operation at
> >     the same precision level. It is very easy for codes to be very
> slightly
> >     different and give quite widely different output.
> >
> >     Indeed, even without the trace, we see (zeroin from rootoned here)
> >
> >     > zeroin(ff, c(0, 10), trace=FALSE)
> >     $root
> >     [1] 1.386294
> >
> >     $froot
> >     [1] -5.658169e-10
> >
> >     $rtol
> >     [1] 7.450581e-09
> >
> >     $maxit
> >     [1] 9
> >
> >     > uniroot(ff, c(0, 10), trace=FALSE)
> >     $root
> >     [1] 1.38629
> >
> >     $f.root
> >     [1] -4.66072e-06
> >
> >     $iter
> >     [1] 10
> >
> >     $init.it <http://init.it>
> >     [1] NA
> >
> >     $estim.prec
> >     [1] 6.103516e-05
> >
> >     >
> >
> >     Is the lack of trace a bug, or at least an oversight? Being able to
> follow iterations is a
> >     classic approach to checking that computations are proceeding as
> they should.
> >
> >     Best, JN
> >
> >     ______________________________________________
> >     R-devel using r-project.org <mailto:R-devel using r-project.org> mailing list
> >     https://stat.ethz.ch/mailman/listinfo/r-devel <
> https://stat.ethz.ch/mailman/listinfo/r-devel>
> >
> >
>

	[[alternative HTML version deleted]]