[R] Rounding error in seq(...)

[David Airey]

 >> can you explain that a little more detailed?
 >> Perhaps I miss the background knowledge - but it seems just absurd  
to me.
 >>
 >> 0.1+0.1+0.1 is 0.3 - there is no rounding involved, is there?
 >>
 >> why is
 >> x <- 0.1 + 0.1 +0.1
 >> not equal to
 >> y <- 0.3
 >
 > Remember that this is in BINARY arithmetic. It's really not any  
stranger
 > than the fact that 1/3 + 1/3 != 2/3 in finite accuracy decimal  
arithmetic
 > (0.33333 + 0.33333 = 0.66666 != 0.66667).

A nice description of this is in section 2.2.4 of Braun and Murdoch's  
book.

They give a striking example.

 > x <- 1:11
 > mean(x)
[1] 6
 > var(x)
[1] 11
 > sum( (x - mean(x))^2) / 10 # two-pass formula for mean
[1] 11
 > (sum(x^2) - 11 * mean(x)^2) / 10 # one-pass formula for mean
[1] 11
 > A <- 1.e10 # add large value to each x_i
 > x <- 1:11 + A
 > var(x) # R function survives
[1] 11
 > sum( (x - mean(x))^2) / 10 # two-pass formula survives
[1] 11
 > (sum(x^2) - 11 * mean(x)^2) / 10 # one-pass formula suffers  
catastrophic loss of precision!
[1] 0