以下代码显然是错误的,问题出在哪里?
i <- 0.1
i <- i + 0.05
i
## [1] 0.15
if(i==0.15) cat("i equals 0.15") else cat("i does not equal 0.15")
## i does not equal 0.15
以下代码显然是错误的,问题出在哪里?
i <- 0.1
i <- i + 0.05
i
## [1] 0.15
if(i==0.15) cat("i equals 0.15") else cat("i does not equal 0.15")
## i does not equal 0.15
由于不是所有的数字都可以在IEEE浮点运算(几乎所有计算机用来表示十进制数并进行数学运算的标准)中精确表示,所以你不会总是得到你期望的结果。这一点尤其真实,因为一些简单的有限小数(如0.1和0.05)在计算机中并不能精确表示,所以对它们进行算术运算的结果可能与“已知”答案的直接表示不完全相同。
这是计算机算术的一个众所周知的限制,并在多个地方进行了讨论:
R FAQ有一个专门讨论这个问题的问题:R FAQ 7.31 Patrick Burns的《R地狱》将第一个“圆圈”专门用于解决这个问题(从第9页开始) David Goldberg的“计算机科学家应该了解的浮点算术知识”,《ACM计算调查》23,1(1991-03),5-48 doi>10.1145/103162.103163(也有修订版本) 《浮点数指南-每个程序员都应该了解的浮点算术知识》 0.30000000000000004.com比较了不同编程语言中的浮点算术 包括几个Stack Overflow问题 为什么浮点数不准确? 为什么十进制数不能在二进制中精确表示? 浮点数运算是否有问题? “浮点数不准确”的规范副本(关于此问题的规范答案的元讨论)在R
中,解决这个问题的标准方法不是使用==
,而是使用all.equal
函数。或者更准确地说,由于all.equal
函数在存在差异时会提供很多细节,所以使用isTRUE(all.equal(...))
。
if(isTRUE(all.equal(i,0.15))) cat("i equals 0.15") else cat("i does not equal 0.15")
i equals 0.15
all.equal
而不是==
的更多示例(最后一个示例旨在显示这将正确显示差异)。0.1+0.05==0.15
#[1] FALSE
isTRUE(all.equal(0.1+0.05, 0.15))
#[1] TRUE
1-0.1-0.1-0.1==0.7
#[1] FALSE
isTRUE(all.equal(1-0.1-0.1-0.1, 0.7))
#[1] TRUE
0.3/0.1 == 3
#[1] FALSE
isTRUE(all.equal(0.3/0.1, 3))
#[1] TRUE
0.1+0.1==0.15
#[1] FALSE
isTRUE(all.equal(0.1+0.1, 0.15))
#[1] FALSE
1.1-0.2
#[1] 0.9
0.9
#[1] 0.9
sprintf("%.54f",1.1-0.2)
#[1] "0.900000000000000133226762955018784850835800170898437500"
sprintf("%.54f",0.9)
#[1] "0.900000000000000022204460492503130808472633361816406250"
sprintf("%a",0.9)
#[1] "0x1.ccccccccccccdp-1"
sprintf("%a",1.1-0.2)
#[1] "0x1.ccccccccccccep-1"
sprintf("%a",1.1-0.2-0.9)
#[1] "0x1p-53"
2^-53
,这很重要,因为这个数字是两个接近1的数值之间最小可表示的差异。 ?.Machine
#....
#double.eps the smallest positive floating-point number x
#such that 1 + x != 1. It equals base^ulp.digits if either
#base is 2 or rounding is 0; otherwise, it is
#(base^ulp.digits) / 2. Normally 2.220446e-16.
#....
.Machine$double.eps
#[1] 2.220446e-16
sprintf("%a",.Machine$double.eps)
#[1] "0x1p-52"
all.equal
。?all.equal
#....
#all.equal(x,y) is a utility to compare R objects x and y testing ‘near equality’.
#....
#all.equal(target, current,
# tolerance = .Machine$double.eps ^ 0.5,
# scale = NULL, check.attributes = TRUE, ...)
#....
a <- c(0.1+0.05, 1-0.1-0.1-0.1, 0.3/0.1, 0.1+0.1)
b <- c(0.15, 0.7, 3, 0.15)
==
不会给出“预期”的结果,而all.equal
不会执行逐元素比较。a==b
#[1] FALSE FALSE FALSE FALSE
all.equal(a,b)
#[1] "Mean relative difference: 0.01234568"
isTRUE(all.equal(a,b))
#[1] FALSE
mapply(function(x, y) {isTRUE(all.equal(x, y))}, a, b)
#[1] TRUE TRUE TRUE FALSE
elementwise.all.equal <- Vectorize(function(x, y) {isTRUE(all.equal(x, y))})
elementwise.all.equal(a, b)
#[1] TRUE TRUE TRUE FALSE
all.equal
,你可以直接复制all.equal.numeric
的相关内部,并使用隐式向量化:tolerance = .Machine$double.eps^0.5
# this is the default tolerance used in all.equal,
# but you can pick a different tolerance to match your needs
abs(a - b) < tolerance
#[1] TRUE TRUE TRUE FALSE
==
更安全,因为它内置了容差。dplyr::near(a, b)
#[1] TRUE TRUE TRUE FALSE
如果应用于浮点数值,标准的R函数%in%
也可能遇到相同的问题。例如:
x = seq(0.85, 0.95, 0.01)
# [1] 0.85 0.86 0.87 0.88 0.89 0.90 0.91 0.92 0.93 0.94 0.95
0.92 %in% x
# [1] FALSE
`%.in%` = function(a, b, eps = sqrt(.Machine$double.eps)) {
any(abs(b-a) <= eps)
}
0.92 %.in% x
# [1] TRUE
dplyr::near
包裹在 any
中也可以用于向量化检查
any(dplyr::near(0.92, x))
# [1] TRUE
加上Brian的评论(以下是原因),你可以通过使用all.equal
来克服这个问题:
# i <- 0.1
# i <- i + 0.05
# i
#if(all.equal(i, .15)) cat("i equals 0.15\n") else cat("i does not equal 0.15\n")
#i equals 0.15
根据Joshua的警告,这里是更新后的代码(感谢Joshua):
i <- 0.1
i <- i + 0.05
i
if(isTRUE(all.equal(i, .15))) { #code was getting sloppy &went to multiple lines
cat("i equals 0.15\n")
} else {
cat("i does not equal 0.15\n")
}
#i equals 0.15
if
语句时,all.equal
函数如果存在差异,不会返回 FALSE
,因此需要用 isTRUE
将其包装起来。 - Joshua Ulrich==
,以使dplanet代码能够正常工作?我也遇到了同样的问题。 - Yacine Hajji这个方法有些hackish,但很快速:
if(round(i, 10)==0.15) cat("i equals 0.15") else cat("i does not equal 0.15")
all.equal(..., tolerance)
参数。all.equal(0.147, 0.15, tolerance=0.05)
为真。 - smcidplyr::near()
是测试两个浮点数向量是否相等的选项。这是来自文档的示例:
sqrt(2) ^ 2 == 2
#> [1] FALSE
library(dplyr)
near(sqrt(2) ^ 2, 2)
#> [1] TRUE
该函数具有内置的容差参数:tol = .Machine$double.eps^0.5
,可以进行调整。默认参数与all.equal()
的默认参数相同。
双精度浮点数中的泛化比较(“<=”,“>=”,“=”):
比较 a <= b:
IsSmallerOrEqual <- function(a,b) {
# Control the existence of "Mean relative difference..." in all.equal;
# if exists, it results in character, not logical:
if ( class(all.equal(a, b)) == "logical" && (a<b | all.equal(a, b))) { return(TRUE)
} else if (a < b) { return(TRUE)
} else { return(FALSE) }
}
IsSmallerOrEqual(abs(-2-(-2.2)), 0.2) # TRUE
IsSmallerOrEqual(abs(-2-(-2.2)), 0.3) # TRUE
IsSmallerOrEqual(abs(-2-(-2.2)), 0.1) # FALSE
IsSmallerOrEqual(3,3); IsSmallerOrEqual(3,4); IsSmallerOrEqual(4,3)
# TRUE; TRUE; FALSE
比较a >= b:
IsBiggerOrEqual <- function(a,b) {
# Control the existence of "Mean relative difference..." in all.equal;
# if exists, it results in character, not logical:
if ( class(all.equal(a, b)) == "logical" && (a>b | all.equal(a, b))) { return(TRUE)
} else if (a > b) { return(TRUE)
} else { return(FALSE) }
}
IsBiggerOrEqual(3,3); IsBiggerOrEqual(4,3); IsBiggerOrEqual(3,4)
# TRUE; TRUE; FALSE
比较a = b:
IsEqual <- function(a,b) {
# Control the existence of "Mean relative difference..." in all.equal;
# if exists, it results in character, not logical:
if ( class(all.equal(a, b)) == "logical" ) { return(TRUE)
} else { return(FALSE) }
}
IsEqual(0.1+0.05,0.15) # TRUE
我曾经遇到过类似的问题,我使用了以下解决方案。
@我找到了一个关于不等间隔切割的解决方案。@我在R中使用了round函数。将选项设置为2位小数并不能解决问题。
options(digits = 2)
cbind(
seq( from = 1, to = 9, by = 1 ),
cut( seq( from = 1, to = 9, by = 1), c( 0, 3, 6, 9 ) ),
seq( from = 0.1, to = 0.9, by = 0.1 ),
cut( seq( from = 0.1, to = 0.9, by = 0.1), c( 0, 0.3, 0.6, 0.9 )),
seq( from = 0.01, to = 0.09, by = 0.01 ),
cut( seq( from = 0.01, to = 0.09, by = 0.01), c( 0, 0.03, 0.06, 0.09 ))
)
[,1] [,2] [,3] [,4] [,5] [,6]
[1,] 1 1 0.1 1 0.01 1
[2,] 2 1 0.2 1 0.02 1
[3,] 3 1 0.3 2 0.03 1
[4,] 4 2 0.4 2 0.04 2
[5,] 5 2 0.5 2 0.05 2
[6,] 6 2 0.6 2 0.06 3
[7,] 7 3 0.7 3 0.07 3
[8,] 8 3 0.8 3 0.08 3
[9,] 9 3 0.9 3 0.09 3
options(digits = 200)
cbind(
seq( from = 1, to = 9, by = 1 ),
cut( round(seq( from = 1, to = 9, by = 1), 2), c( 0, 3, 6, 9 ) ),
seq( from = 0.1, to = 0.9, by = 0.1 ),
cut( round(seq( from = 0.1, to = 0.9, by = 0.1), 2), c( 0, 0.3, 0.6, 0.9 )),
seq( from = 0.01, to = 0.09, by = 0.01 ),
cut( round(seq( from = 0.01, to = 0.09, by = 0.01), 2), c( 0, 0.03, 0.06, 0.09 ))
)
[,1] [,2] [,3] [,4] [,5] [,6]
[1,] 1 1 0.1 1 0.01 1
[2,] 2 1 0.2 1 0.02 1
[3,] 3 1 0.3 1 0.03 1
[4,] 4 2 0.4 2 0.04 2
[5,] 5 2 0.5 2 0.05 2
[6,] 6 2 0.6 2 0.06 2
[7,] 7 3 0.7 3 0.07 3
[8,] 8 3 0.8 3 0.08 3
[9,] 9 3 0.9 3 0.09 3