我不小心发现了 R base 中的 as.integer
或 det
函数的奇怪错误。有谁知道这里发生了什么以及如何预防?
我正在计算以下 3×3 矩阵的行列式:
mat <- matrix(c(15, 6, 116, 10, 13, 16, 14, 23, 56), ncol = 3)
看起来像这样:
[,1] [,2] [,3]
[1,] 15 10 14
[2,] 6 13 23
[3,] 116 16 56
有两件事很容易看出:所有条目都是整数,并且六组三条目中的每一组都至少包含在同一行或列中的两个条目至少一个偶数。因此行列式必须是偶数。
通过键入 det(mat)
向 R 询问此行列式的实际值,它返回看起来像偶数的东西:8952
。但是你瞧:在 R 的内心深处,它实际上是一个非整数或奇数,因为在键入 as.integer(det(mat))
时我们得到 8951
.
这是怎么回事? 8951显然是错误的。此外,不太明显的是,值 8952 是正确的,可以用笔和纸看出。
所以我的问题是:
这是怎么回事?
当被要求计算整数矩阵的行列式时,我如何强制 R 给我正确整数值?
最佳答案
根本原因:is.integer
截断而不是舍入和 float 学记录的中间值来解释第二个结果,结合 print 控制台 REPL 的一部分,用于解释
det(mat)
的初始结果:
print( det(mat), digits =16)
[1] 8951.999999999993
理论上的答案很可能是 8952,但 R 不是符号数学引擎。
您可以使用 Rmpfr 包(如@BenBolker 所建议的那样)来提高精度级别:
library(Rmpfr)
mat <- mpfr(mat, 64)
as.integer( det(mat) )
[1] 8952
as.integer
截断而不是舍入。参见 ?as.integer
。 R 可以在不损失精度的情况下处理整数的加法或乘法,但一旦发生除法,就可能会出现浮点错误。 (实际上,问题的出现是因为 det
的默认设置是使用 determinant
和 log=TRUE
,然后对复数结果的模取幂。)从帮助页面的值部分:
Non-integral numeric values are truncated towards zero (i.e.,
as.integer(x)
equalstrunc(x)
there), and imaginary parts of complex numbers are discarded (with a warning).
https://stackoverflow.com/questions/65401271/