From ggrothendieck at gmail.com  Sat Aug  4 19:05:15 2007
From: ggrothendieck at gmail.com (Gabor Grothendieck)
Date: Sat, 4 Aug 2007 13:05:15 -0400
Subject: [Rd] Optimization in R
In-Reply-To: <46B4A200.10109@stats.uwo.ca>
References: <20070804051231.GB3016@econ.upenn.edu>
	<46B47E53.2040801@stats.uwo.ca>
	<20070804143032.GC16431@econ.upenn.edu> <46B4A200.10109@stats.uwo.ca>
Message-ID: <971536df0708041005h55c2488by76709cf4fb989690@mail.gmail.com>

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
>



[Rd] Optimization in R
Gabor Grothendieck