在 numpy 中使用 exp 和最小正 float128 避免下溢

Question

我正在尝试计算以下比率： w(i) / (sum(w(j)) 其中w 使用指数递减函数更新，即w(i) = w(i) * exp(-k)，k 是一个正参数。所有的数字都是非负数。 然后将该比率用于公式（乘以常数并添加另一个常数）。不出所料，我很快就遇到了下溢问题。

我想这经常发生，但有人可以给我一些关于如何处理这个问题的参考吗？我没有找到合适的转换，所以我尝试做的一件事是将一些最小正数设置为安全阈值，但我没有设法找到最小正浮点数（我在 numpy.float128 中表示数字）。我如何才能在我的机器上实际获得最小 positive 这样的数字？ 代码如下所示：

w = np.ones(n, dtype='float128')
lt = np.ones(n)
for t in range(T):
    p = (1-k) * w / w.sum() + (k/n)
    # Process a subset of the n elements, call it set I, j is some range()
    for i in I: 
        s = p[list(j[i])].sum()
        lt /= s
        w[s] *= np.exp(-k * lt)

其中 k 是 (0,1) 中的某个常数，n 是数组的长度

Answer 1

当处理指数级的小数字时，通常最好在日志空间中工作。例如，log(w*exp(-k)) = log(w) - k，除非 k 本身呈指数级大或 w 为零，否则不会出现任何上溢/下溢问题。而且，如果w 为零，numpy 将正确返回-inf。然后，在求和时，您会考虑出最大的项：

log_w = np.log(w) - k
max_log_w = np.max(log_w)
# Individual terms in the following may underflow, but then they wouldn't
# contribute to the sum anyways.
log_sum_w = max_log_w + np.log(np.sum(np.exp(log_w - max_log_w)))
log_ratio = log_w - log_sum_w

这可能不是您想要的，因为您可以完全排除 k（假设它是一个常量而不是数组），但它应该可以帮助您。

Scikit-learn 实现了与extmath.logsumexp 类似的东西，但与上面的基本相同。