0

I'm translating a package from R to IML, and this will be free online when it will be done :). I gain different results from the decomposition of a big matrix, both results seems the same when you look at them, but, if for example I take first 2 columns of U and do U'*U, my 2x2 matrix will be quite different (U_11 = 1.1e-17 and U_11 =1.4e-17). The difference is very small (3e-18) and this lead me to think that it could be something related to the number of decimals that each software use, SAS IML and R. Anyone knows something more about this topic? How can I test this? Thank you.

stat
  • 699
  • 3
  • 10
  • SAS uses double floating points. Do you have more than 16 (or so) digits of precision here? – Joe Jun 29 '15 at 21:30
  • In R those values are effectively zero. R has 16 decimal digits of precision, but after you start multiplying numbers together there is considerable diminution in precision and the `all.equal` function use this as its tolerance: `.Machine$double.eps ^ 0.5` == 1.490116e-08 – IRTFM Jun 29 '15 at 21:44
  • When you finish the project, please consider posting to the SAS/IML File Exchange: http://blogs.sas.com/content/iml/2014/07/16/sasiml-file-exchange.html – Rick Jun 30 '15 at 10:00

1 Answers1

1

In statistics, we describe very small differences as "statistically insignificant." To a numerical analyst, differences that are less than "machine epsilon" (.Machine$double.eps in R or constant("maceps") in SAS) are almost always "numerically insignificant."

Both SAS and R use double precision computations and are probably calling similar numerical libraries. For a difference that small, I would conjecture that the reason is not algorithmic but is because of different compiler flags and optimization flags that each software uses.

Even within a single product, computing a result in two different orders can result in small differences like this. For example, run the following DATA step:

data _null_;
x = (1 + 1 + 1 + 1 + 1 + 1 + 1) / 7;
y = (1/7 + 1/7 + 1/7 + 1/7 + 1/7 + 1/7 + 1/7);
diff = x - y;
put diff=;
run;

My sugggestion is to ignore "numerically insignificant" results when comparing different software. For more information about floating point computations, see The Floating Point Guide. For the real nitty gritty, see "What every computer scientist should know about floating-point arithmetic"

Rick
  • 1,210
  • 6
  • 11
  • thank you Rick, I accept your answer as correct, I'll surely post my work on IML file exchange – stat Jul 07 '15 at 12:32