From ripley at stats.ox.ac.uk  Wed May 17 18:46:45 2006
From: ripley at stats.ox.ac.uk (ripley at stats.ox.ac.uk)
Date: Wed, 17 May 2006 18:46:45 +0200 (CEST)
Subject: [Rd] optim "CG" bug w/patch proposal (PR#8786)
Message-ID: <20060517164645.59D782477B@slim.kubism.ku.dk>

On Wed, 17 May 2006, maechler at stat.math.ethz.ch wrote:

>
>>>>>> "Duncan" == Duncan Murdoch <murdoch at stats.uwo.ca>
>>>>>>     on Tue, 16 May 2006 08:34:06 -0400 writes:
>
>    Duncan> On 5/16/2006 4:56 AM, westfeld at inf.tu-dresden.de
>    Duncan> wrote:
>    >> Probably I included too much at once in my bug report. I
>    >> can live with an unfulfilled wishlist and thank you for
>    >> thinking about it. The "badly-behaved" function is just
>    >> an example to demonstrate the bug I reported. I think it
>    >> is a bug if optim returns (without any warning) an
>    >> unmatching pair of par and value: f(par) != value. And it
>    >> is easily fixed.
>
>    >>  Andreas
>
>    Duncan> I agree with you that on return f(par) should be
>    Duncan> value.  I agree with Brian that changes to the
>    Duncan> underlying strategy need much more thought.
>
> I agree (to both).
> However, isn't Andreas' patch just fixing the problem
> and not changing the underlying strategy at all?
> [No, I did not study the code in very much detail ...]

The (minor) issue is that x is updated but not f(x).  I think the intended 
stategy was to update neither, so Andreas' patch was a change of stategy. 
In particular, a question is if this should be marked as a convergence 
failure.  But people really need to read the reference before commenting,
and I at least need to find the time to do so in more detail.

> Martin Maechler
>
>    >> Prof Brian Ripley wrote:
>    >>
>    >>> [Sorry for the belated reply: this came in just as I was leaving for a
>    >>> trip.]
>    >>>
>    >>> I've checked the original source, and the C code in optim does
>    >>> accurately reflect the published algorithm.
>    >>>
>    >>> Since your example is a discontinuous function, I don't see why you
>    >>> expect CG to work on it.  John Nash reports on his extensive
>    >>> experience that method 3 is the worst, and I don't think we should let
>    >>> a single 2D example of a badly-behaved function override that.
>    >>>
>    >>> Note that no other optim method copes with the discontiuity here: had
>    >>> your reported that it would have been clear that the problem was with
>    >>> the example.
>    >>>
>    >>> On Fri, 21 Apr 2006, westfeld at inf.tu-dresden.de wrote:
>    >>>
>    >>>> Dear R team,
>    >>>>
>    >>>> when using optim with method "CG" I got the wrong $value for the
>    >>>> reported $par.
>    >>>>
>    >>>> Example:
>    >>>> f<-function(p) {
>    >>>> if (!all(p>-.7)) return(2)
>    >>>> if (!all(p<.7)) return(2)
>    >>>> sin((p[1])^2)*sin(p[2])
>    >>>> }
>    >>>> optim(c(0.1,-0.1),f,method="CG",control=list(trace=0,type=1))
>    >>>> $par 19280.68 -10622.32
>    >>>> $value -0.2346207 # should be 2!
>    >>>>
>    >>>> optim(c(0.1,-0.1),f,method="CG",control=list(trace=0,type=2))
>    >>>> $par 3834.021 -2718.958
>    >>>> $value -0.0009983175 # should be 2!
>    >>>>
>    >>>> Fix:
>    >>>> --- optim.c     (Revision 37878)
>    >>>> +++ optim.c     (Arbeitskopie)
>    >>>> @@ -970,7 +970,8 @@
>    >>>> if (!accpoint) {
>    >>>> steplength *= stepredn;
>    >>>> if (trace) Rprintf("*");
>    >>>> -                           }
>    >>>> +                           } else
>    >>>> +                               *Fmin = f;
>    >>>> }
>    >>>> } while (!(count == n || accpoint));
>    >>>> if (count < n) {
>    >>>>
>    >>>> After fix:
>    >>>> optim(c(0.1,-0.1),f,method="CG",control=list(trace=0,type=1))
>    >>>> $par 0.6993467 -0.4900145
>    >>>> $value -0.2211150
>    >>>> optim(c(0.1,-0.1),f,method="CG",control=list(trace=0,type=2))
>    >>>> $par 3834.021 -2718.958
>    >>>> $value 2
>    >>>>
>    >>>> Wishlist:
>    >>>
>    >> [wishlist deleted]
>    >>
>    >>
>
>    Duncan> ______________________________________________
>    Duncan> R-devel at r-project.org mailing list
>    Duncan> https://stat.ethz.ch/mailman/listinfo/r-devel
>
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
>
>

-- 
Brian D. Ripley,                  ripley at stats.ox.ac.uk
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford,             Tel:  +44 1865 272861 (self)
1 South Parks Road,                     +44 1865 272866 (PA)
Oxford OX1 3TG, UK                Fax:  +44 1865 272595



[Rd] optim "CG" bug w/patch proposal (PR#8786)
ripley at stats.ox.ac.uk