13
$\begingroup$

Mathematically, it is evident that var = sd^2 (or var^0.5 = sd), but R seems to answer TRUE for the latter and FALSE for the former. To my knowledge, both var() and sd() function use n-1 for degree of freedom. What would be the problem?

data <- c(20, 20, 30, 20, 16, 13, 20)

var(data)^0.5 == sd(data)
var(data) == sd(data)^2
Improve this question
$\endgroup$
8
  • 17
    $\begingroup$ There is no problem. This is about numerical precision. Use the all.equal function instead to acknowledge the fact that your machine represents real numbers as floating-point numbers with a specific numerical precision. The documentation of all.equal also explains .Machine$double.eps. $\endgroup$
    – dipetkov
    Sep 26, 2022 at 18:10
  • $\begingroup$ Thanks! Will find the document you mentioned. $\endgroup$
    – HYL
    Sep 26, 2022 at 18:31
  • 3
    $\begingroup$ Try this: print(x <- 1 + 1e-15, 17); sqrt(x)^2 == x $\endgroup$
    – whuber ♦
    Sep 26, 2022 at 19:39
  • 6
    $\begingroup$ This is R FAQ 7.31. $\endgroup$
    – Stephan Kolassa
    Sep 27, 2022 at 5:43
  • 2
    $\begingroup$ A good way to visualize the problem: calculate sqrt(13), write the result on a piece of paper, now square your intermediary result and see if you get 13 again. This is what R (has to) do - it's a calculator, not a symbolic equation solver, it doesn't retain knowledge that a specific intermediary result is supposed to be sqrt(13). $\endgroup$
    – xLeitix
    Sep 27, 2022 at 13:07

4 Answers 4

Reset to default
29
$\begingroup$

Please, check

0.4^2 == 0.16
0.4 == 0.16^0.5

It's not a statistical issue. It's a rounding issue in floating point arithmetics in the language. The same applies to Python

0.4**2 == 0.16
0.4 == 0.16**0.5

None of these floating point numbers is stored in a computer's memory as two- or three-decimal digits, but rather in a binary format. For instance, 0.16 (8-byte) is stored as 11111111000100011110101110000101000111101011100001010001111011, while 0.4**2 gives us 11111111000100011110101110000101000111101011100001010001111100. You can play with this decimal to binary conversion using the following Python code:

import struct
bin(struct.unpack('!Q',struct.pack('!d',0.16))[0])
bin(struct.unpack('!Q',struct.pack('!d',0.4**2))[0])
Improve this answer
$\endgroup$
4
  • 1
    $\begingroup$ Hey Alex - I understand that the above example from mkt is a issue of floating point, but still cannot quite grasp why the two-digit decimal point example of yours(0.4**2) is not identical with 0.16. Can you please further elaborate? $\endgroup$
    – HYL
    Sep 26, 2022 at 20:06
  • 2
    $\begingroup$ Hey. I added elaboration to my answer above $\endgroup$
    – Alex
    Sep 26, 2022 at 20:34
  • $\begingroup$ Thanks so much! $\endgroup$
    – HYL
    Sep 26, 2022 at 20:38
  • $\begingroup$ @Hwi-youngLee If this answers your question, consider accepting it (or whichever other answer you see as best) by clicking on the green tick mark to the left of the answer. $\endgroup$
    – mkt
    Sep 29, 2022 at 7:58
20
$\begingroup$

Recommended read for you is

Goldberg, D. (1991). What every computer scientist should know about floating-point arithmetic. ACM computing surveys (CSUR), 23(1), 5-48.

TL;DR you should never use == (or its equivalent in any other programming language) to compare floating point numbers because they represent numbers only approximately. As in your example, two different computations can lead to floating point representations that are different though close enough given the guarantees given by such representation. What you should do instead is to compare them with some precision, e.g. $|x - y| < \epsilon$, but this is also not a perfect solution and there are better ones. Happily, in most cases, they would be already implemented for you in the programming language you are using (like the mentioned all.equal in R).

Improve this answer
$\endgroup$
2
  • 2
    $\begingroup$ Dunno why this was downvoted. Just printing sd, one can see that it uses sqrt(var(x)) to perform the calculation. The devil is thus hidden in floating point arithmetic. And one needs to be aware of the problems hidden there. $\endgroup$
    – Colombo
    Sep 27, 2022 at 2:25
  • 2
    $\begingroup$ +1 That paper is a classic. $\endgroup$
    – Galen
    Sep 27, 2022 at 4:41
10
$\begingroup$

Correct to 14 decimal points, it's a floating point number issue:

options(digits = 20)
> var(data)
[1] 27.476190476190474499
sd(data)^2
[1] 27.476190476190470946
Improve this answer
$\endgroup$
5
$\begingroup$

Check the differences between those quantities. They are likely to be small.

Computers are funky about how they store long decimals.

Improve this answer
$\endgroup$
1
  • 1
    $\begingroup$ Thanks, it seems that the difference is indeed very small, but still feeling kind of awkward that it produced FALSE... $\endgroup$
    – HYL
    Sep 26, 2022 at 18:11

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.