[Rd] Optimization in R

[Duncan Murdoch]

On 04/08/2007 1:05 PM, Gabor Grothendieck wrote:
> 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.

Why?  This would make sense if optim() had a lot of code common to all 
methods, but it doesn't.

Duncan Murdoch

> 
> 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
>>
> 
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel