[R] Interesting quirk with fractions and rounding / using == for floating point
J C Nash
profjcnash at gmail.com
Sun Apr 23 18:12:38 CEST 2017
I've some stuff too, but I need to look it up. A few years ago I built
a small test spreadsheet for Gnumeric when working with Jody Goldberg.
In the early 2000s, Jody contacted R (I think Duncan Murdoch) to ask if
it was OK for Gnumeric to use R's distribution function approximations (Yes).
Someone can correct me if wrong, but apparently a few weeks later Jody
sent some improvements, and I believe these were incorporated. Good
open-source improvements like that deserve to be recorded in a better way
than I have done here. But in any event, I can dig out that material, as
well as some of Kahan's PARANOIA materials.
Suggest taking the discussion off-line from here on unless we come up
with some important bugs or improvements or ...
On 2017-04-23 11:42 AM, Richard M. Heiberger wrote:
> I would be happy to participate in designing the test suite you suggest.
> About a year ago I revised FAQ 7.31, based on my talk at the Aalberg R
> conference. It now points, in addition to the Goldberg paper that has
> been referenced there for a long time, to my appendix on precision.
> Here is the link from the FAQ
> A discussion with many easily followed examples is in Appendix G
> "Computational Precision and Floating Point Arithmetic", pages 753-771
> of _Statistical Analysis and Data Display: An Intermediate Course with
> Examples in R_, Richard M. Heiberger and Burt Holland (Springer 2015,
> second edition). This appendix is a free download from
> The appendix is based on the paper I gave at the Aalberg R
> Conference. It uses the Rmpfr package to illustrate exactly what
> is happening at the level of the machine arithmetic. The
> investigation and discussion of the effect of optimization on
> floating point arithmetic should use the Rmpfr tools.
> On Sun, Apr 23, 2017 at 10:06 AM, J C Nash <profjcnash at gmail.com> wrote:
>> Yes. I should have mentioned "optimizing" compilers, and I can agree with
>> trusting exact equality", though I consider conscious use of equality tests
>> Optimizing compilers have bitten me once or twice. Unfortunately, a lot of
>> floating-point work requires attention to detail. In the situation of
>> for convergence, the alternatives to the test I propose require quite a lot
>> code, with variants for different levels of precision e.g., single, double,
>> There's no universal answer and we do have to look "under the hood". A
>> nasty instance (now fortunately long gone) was the Tektronix 4051 graphics
>> where the comparisons automatically included a "fuzz". There was a FUZZ
>> command to
>> set this. Sometimes the "good old days" weren't! Today's equivalent is when
>> is an upstream change to an "optimization" that changes the manner of
>> as in Peter's examples. If we specify x <- y * (a / b), then we should not
>> x <- (y * a) / b.
>> A slightly related case concerns eigenvectors / singular vectors when there
>> degenerate values (i.e., two or more equal eigenvalues). The vectors are
>> determined only to lie in a (hyper)plane. A large computer contracting firm
>> two weeks of high-priced but non-numerical help trying to find the "error"
>> in either
>> an IBM or Univac program because they gave very different eigenvectors.
>> And in my own field of function minimization / nonlinear least squares, it
>> is quite
>> common to have multiple minima or a plateau.
>> Does anyone know if R has a test suite to check some of these situations? If
>> I'll be happy to participate in generating some. They would not need to be
>> very often, and could be useful as a didactic tool as well as checking for
>> differences in platforms.
>> On 2017-04-23 09:37 AM, peter dalgaard wrote:
>>>> On 23 Apr 2017, at 14:49 , J C Nash <profjcnash at gmail.com> wrote:
>>>> So equality in floating point is not always "wrong", though it should be
>>>> with some attention to what is going on.
>>>> Apologies to those (e.g., Peter D.) who have heard this all before. I
>>>> there are many to whom it is new.
>>> Peter D. still insists on never trusting exact equality, though. There was
>>> at least one case in the R sources where age-old code got itself into a
>>> condition where a residual terme that provably should decrease on every
>>> iteration oscillated between two values of 1-2 ulp in magnitude without ever
>>> reaching 0. The main thing is that you cannot trust optimising compilers
>>> these days. There is, e.g., no guarantee that a compiler will not transform
>>> (x_new + offset) == (x_old + offset)
>>> (x_new + offset) - (x_old + offset) == 0
>>> (x_new - x_old) + (offset - offset) == 0
>>> to.... well, you get the point.
>> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
More information about the R-help