[Rd] Optimization in R

[Andrew Clausen]

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?

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.)

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.

> 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