Re: [Rd] NaN and linear algebra

stefano iacus Wed, 23 Mar 2005 09:59:43 -0800


On 23/mar/05, at 09:04, Martin Maechler wrote:

"Bill" == Bill Northcott <[EMAIL PROTECTED]>
    on Wed, 23 Mar 2005 10:19:22 +1100 writes:
Bill> On 23/03/2005, at 12:55 AM, Simon Urbanek wrote:
As I see it, the MacOS X behaviour is not IEEE-754 compliant.
I had a quick look at the IEEE web site and it seems quite clear that NaNs should not cause errors, but propagate through calculations to be tested at some appropriate (not too frequent) point.
This is not quite correct and in fact irrelevant to the problem you describe. NaNs may or may not signal, depending on how they are used. Certain operations on NaN must signal by the IEEE-754 standard. The error you get is not a trap, and it's not a result of a signal check, either. The whole problem is that depending on which algorithm is used, the NaNs will be used different ways and thus may or may not use signaling operations.
Bill> It may not violate the letter of IEEE-754 because matrix calculations Bill> are not covered, but it certainly violates the spirit that arithmetic Bill> should be robust and programs should not halt on these sorts of errors.
I don't consider the `solve' error a bug - in fact I would rather get an error telling me that something is wrong (although I agree that the error is misleading - the error given in Linux is a bit more helpful) than getting a wrong result.
Bill> You may prefer the error, but it is not in the sprit of robust Bill> arithmetic. ie
d<-matrix(NaN,3,3)
f<-solve(d)
Bill> Error in solve.default(d) : Lapack routine dgesv: system is exactly Bill> singular
f
    Bill> Error: Object "f" not found
If I would mark something in your example as a bug that would be
det(m)=0, because it should return NaN (remember, NaN==NaN is FALSE;
furthermore if det was calculated inefficiently using Laplace
expansion, the result would be NaN according to IEEE rules). det=0 is
consistent with the error given, though. Should we check this in R
before calling Lapack - if the vector contains NaNs, det/determinant
should return NaN right away?
Bill> Clearly det(d) returning 0 is wrong. As a result based on a Bill> computation including a NaN, it should return NaN. The spirit of Bill> IEEE-754 is that the programmer should choose the appropriate point at Bill> which to check for NaNs. I would interpret this to mean the R Bill> programmer not the R library developer. Surely that is why R provides Bill> the is.nan function.

d
    Bill> [,1] [,2] [,3]
    Bill> [1,]  NaN  NaN  NaN
    Bill> [2,]  NaN  NaN  NaN
    Bill> [3,]  NaN  NaN  NaN
is.nan(solve(d))
Bill> Error in solve.default(d) : Lapack routine dgesv: system is exactly Bill> singular

Bill> This is against the spirit of IEEE-754 because it halts the program.

is.nan(det(d))
    Bill> [1] FALSE
    Bill> That is plain wrong.
Many functions in R will actually bark at NaN inputs (e.g. qr, eigen,
...) - maybe you're saying that we should check for NaNs in solve
before proceeding and raising an error?
    Bill> However, this problem is in the Apple library not R.
    Bill> Bill Northcott
Indeed!
I pretty much entirely agree with your points, Bill, and would
tend to declare that this Apple library is ``broken''
for building a correctly running R.
Let me ask one question I've been wondering about now for a
while:
  Did you run "make check" after building R,
  and "make check" ran to completion without an error?

Martin, the use vecLib was (since now?) the suggested way to build R on OS X with blas support. This is the way OS X binaries are distributed and of course they pass make check (I've checked again just now and both 2.0.0 and r-devel pass make check) so probably the NaN is a special case with don't test for LAPACK? stefano


If yes (which I doubt quite a bit), there *is* a bug in R's
quality control / quality assurance tools -- and I would want to
add a check for the misbehavior you've mentioned.

Martin Maechler, ETH Zurich

______________________________________________
R-devel@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel


______________________________________________
R-devel@stat.math.ethz.ch mailing list
https://stat.ethz.ch/mailman/listinfo/r-devel

Re: [Rd] NaN and linear algebra

Reply via email to