[Rd] Optimization in R
Gabor Grothendieck
ggrothendieck at gmail.com
Sat Aug 4 19:05:15 CEST 2007
I think it would be desirable for optim to have a dispatching mechanism
that allows users to add their own optimization techniques to those
provided without having to modify optim and without having to come
up with a new visible function. For example, if we call optim as
optim(...whatever..., method = "X") then that might dispatch
optim.X consistent with S3 except here X is explicitly given rather
than being a class.
On 8/4/07, Duncan Murdoch <murdoch at stats.uwo.ca> wrote:
> On 04/08/2007 10:30 AM, Andrew Clausen wrote:
> > Hi Duncan,
> >
> > On Sat, Aug 04, 2007 at 09:25:39AM -0400, Duncan Murdoch wrote:
> >> This is interesting work; thanks for doing it. Could I make a
> >> suggestion? Why not put together a package containing those test
> >> optimization problems, and offer to include other interesting ones as
> >> they arise?
> >
> > Good idea.
> >
> >> You could also include your wrapper for the gsl function
> >> and your own improvements to optim.
> >
> > I have submitted my gsl multimin() wrapper for inclusion into the Rgsl package,
> > which seems to be the "right" place for it. (Its maintainer is currently
> > enjoying a vacation :)
> >
> > It would be nice if all of these methods could be accessed with the existing
> > optim() interface, so that users could easily try different algorithms.
> > That is, users could specify method="BFGS" or "BFGS-R" or "BFGS-GSL", etc.
> > Is there a convenient mechanism for packages registering new methods?
>
> No there isn't, but I don't really see it as necessary. The R optim()
> function is reasonably short, mainly setting up defaults specific to
> each of the supported optimization methods. The C do_optim() function
> has a bit more code common to the methods, but still only about 50
> lines. You could copy this common part into your own function following
> the same argument and return value conventions and a user would just
> need to change the function name in addition to specifying your method,
> not really that much harder. (You could call your function optim and
> use it as a wrapper for stats::optim if you wanted to avoid even this,
> but I wouldn't recommend it, since then behaviour would depend on the
> search list ordering.)
>
> There's a small risk that the argument list to optim() will change
> incompatibly sometime in the future, but I think that's unlikely.
>
> > One incompatibility with my BFGS implementation is that it returns the
> > *inverse* Hessian, which is a natural by-product of the BFGS algorithm.
> > Indeed, R's existing BFGS implementation also calculates the inverse Hessian,
> > and discards it! (Disclaimer: as far as I know, there are no theorems that say
> > that BFGS's inverse Hessians are any good. In practice, they seem to be.)
>
> If you return more than optim() returns it shouldn't have any serious
> ill effects.
>
> > The inverse Hessian is more useful than the Hessian for statistics because it
> > gives the variance-covariance matrix for maximum likelihood estimators. When
> > the Hessian is close to being singular (aka "computationally singular"),
> > solve() sometimes fails to invert it.
> >
> > I think this means we should change the optim() interface. For example, an
> > extra logical parameter, "inv.hessian" could control whether an inv.hessian
> > matrix is returned.
>
> I'd suggest always returning it if it's useful and the calculation is
> reliable and cheap. Adding "inv.hessian" as a general parameter would
> be troublesome with some of the other methods, where the inverse Hessian
> isn't already calculated, because of the inversion problem you mention
> above.
>
> Duncan Murdoch
>
> >
> >> On your first point: I agree that a prototype implementation in R makes
> >> sense, but I suspect there exist problems where the overhead would not
> >> be negligible (e.g. ones where the user has written the objective
> >> function in C for speed). So I think you should keep in mind the
> >> possibility of moving the core of your improvements to C once you are
> >> happy with them.
> >
> > Fair enough.
> >
> > Cheers,
> > Andrew
> >
> > ______________________________________________
> > R-devel at r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-devel
>
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
>
More information about the R-devel
mailing list