[Rd] optim "CG" bug w/patch proposal (PR#8786)

[maechler at stat.math.ethz.ch]

>>>>> "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 ...]

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