Matlab 日志与 Numpy np.log 不匹配

Question

在将旧的 Matlab 代码重写为 NumPy 时，我注意到对数计算方面的差异。在 NumPy 中，我使用np.log，Matlab 使用log函数。

b = [1 1 2 3 5 1 1];
p = b ./ sum(b);
sprintf('log(%.20f) = %.20f', p(5), log(p(5)))

import numpy as np
b = np.array([1, 1, 2, 3, 5, 1, 1])
p = b.astype('float64') / np.sum(b)
print(f'log({p[4]:.20f}) = {np.log(p[4]):.20f}')

对于配备 M1 芯片的 MacBook Pro 2020，我发现小数点后第 16 位不匹配。

log(0.35714285714285715079) = -1.02961941718115834732  # Matlab
log(0.35714285714285715079) = -1.02961941718115812527  # NumPy

我想得到完全相同的结果。知道如何修改我的 Python 代码吗？

Answer 1

默认情况下，MATLAB 和 numpy 都使用具有 52 位尾数的64 位浮点数。这意味着两个 float64 数字之间的最小相对步长是2**-52 = 2.2e-16。这意味着 16 位之后的任何小数都没有意义。您看到的差异可能是由于实现略有不同造成的。您可以使用以下方法检查这一点

np.nextafter(a, 1)-a

为a = np.log(0.35714285714285715079)你得到的大概就是机器精度2.2e-16的大小。 np.finfo(np.float64).eps

即使您查看输入：您提供的小数位数也超出了完全定义 64 位浮点所需的位数。我们可以将显示的小数位数设置为 100，但它仍然只打印 17 位数字，原因如下：

>>> np.set_printoptions(precision=100)
>>> np.array([0.35714285714285715079])
array([0.35714285714285715])                       

MATLAB 和 numpy 之间的差异甚至可能是由于对总和重新排序而引起的，因为浮点加法不具有关联性。如果您确实依赖于小数点后 16 位，那么您应该使用 64 位浮点数以外的其他值。我建议您熟悉浮点类型的实现方式，因为这在使用科学软件时至关重要。如果您愿意，我建议您查看 numpy 的源代码，了解它是如何实现的，并将其与其他开放库进行比较。