Different result when using sqrt(x) and x^0.5 in R

Question

I have an issue where there is a difference in the results when I use sqrt(x) rather than x^0.5.

The calculations are being carried out on floating point numbers such as:

0.002296438

Trouble is this truncated version as displayed in Rstudio does not replicate the problem. However the non-truncated version does (any idea how I can get the non-truncated version to display so I can show a working example)?

The errors are indeed small of the order of e^-18 which are not so worrying in themselves. However over even moderately large data sets (10,000 date points) these errors compound to give errors in the variance estimate at the 4th decimal place which is more concerning!

I realize that R is only accurate to 16 decimal places, see answer from nullglob below but these errors seems to be systematic? Every time you run sqrt(x) and x^0.5 they both produce the same answer each time. However these answers are still different from each other.

Formatting Decimal places in R

Is one version considered to be more accurate than the other?

Baz

OK here is the example

[1] 0.002296437934635199226707
> test4=sqrt(0.002296437934635199226707)
> test5=0.002296437934635199226707^0.5
> test6=test5-test4
test6
[1] 6.938894e-18

Re: printing more decimal places, try `print(number, digits=22)` (maximum `digits` is 22). For a permanent change, set `options(digits=22)`. — mathematical.coffee, Jul 24 '15 at 14:04
`?sqrt` says that `sqrt(z) == z^0.5`, so it would be interesting to see an example where this is not true. — mathematical.coffee, Jul 24 '15 at 14:05
another way to test: `x <- 0.04792116; identical(sqrt(x),x^0.5)` — Ben Bolker, Jul 24 '15 at 14:06
Can you show a concrete example where they are not equal? I can't find one. — joran, Jul 24 '15 at 14:16
Hi there pls see the example added to the bottom of the post. — Bazman, Jul 24 '15 at 14:43
Most likely R is using a different C call for each case [as shown here](http://stackoverflow.com/questions/29098168/why-does-0-5-appear-to-be-more-efficient-than-sqrt). None is more accurate than the other (floating point and accuracy don't go well together). To avoid unwanted effects you should round all results to a meaningful level of precision — Panagiotis Kanavos, Jul 24 '15 at 15:00
Please show the evidence that the errors compound to the 4th decimal place of variance. — cxrodgers, Jul 24 '15 at 16:23

score 1 · Accepted Answer · answered Jul 24 '15 at 16:13

Confirmed (more or less) the suggestion of different C calls with the following C code, which calls the sqrt() and pow() functions from the system math library.

#include <math.h>
#include <stdio.h>

int main(int argc, char **argv) {
    double x = 0.002296437934635199226707;
    long double y = 0.002296437934635199226707;

    printf("%1.22g\n",sqrt(x)-pow(x,0.5));
    printf("%1.22Lg\n",sqrtl(y)-powl(y,0.5));

    return(0);
}

and prints

-6.938893903907228377648e-18
-3.388131789017201356273e-21

on my system (32-bit Ubuntu 12.04), i.e. the results for long double are slightly closer than the results for double. The commenters above are correct, though, that if you're worrying about this level of precision you probably have larger problems; is the rest of your algorithm stable to this level of precision?

In order to dig deeper (which probably isn't worth it except for intellectual curiosity) you'd have to find out more about the implementation of these two functions in the system libraries for your particular OS.

Whereas on my system (64 bit CentOS 6) I can't reproduce this in either C or R. I get a difference of 0. It's probably a combination of OS, compiler and maybe even processor. — Nick Kennedy, Jul 24 '15 at 16:52
Actually I can reproduce the long double result but not the double one. — Nick Kennedy, Jul 24 '15 at 17:41

Different result when using sqrt(x) and x^0.5 in R

1 Answers1