python - 在计算正常样本落在 python 中某个区间内的概率的对数时防止下溢

Question

给定两个实数 x 和 y，我想在 python 中计算以下函数:

log Pr [ x <= t <= y ],

其中 t 从正态分布中采样。

一个天真的实现是使用 scipy.stats.norm .

np.log(scipy.stats.norm.cdf(y) - scipy.stats.norm.cdf(x))

不幸的是，当 x 和 y 远离 0 时，这会导致下溢。如何防止这样的数字错误？

Answer 1

如果在 logspace 中完成，这个问题会稳定得多。

技巧是对小于零的值使用 scipy.stats.norm.logcdf，对大于零的值使用 scipy.stats.norm.logsf。

这与用于计算 log(exp(y) - exp(x)) 的稳定算法相结合给出了合理的结果

import numpy as np
from scipy.stats import norm

def log_subtract(x, y):
    return x + np.log1p(-np.exp(y-x))

def lnprob(x, y):
    if x < 0:
        return log_subtract(norm.logcdf(y), norm.logcdf(x))
    else:
        return log_subtract(norm.logsf(x), norm.logsf(y))