SciPy binom.logsf - 只是一个便利函数？

Question

我认为

scipy.stats.binom.log*

函数系列能够返回更大范围的值，因为不知何故，它们将所有计算保留在对数空间中。

但是这个例子表明，由于相同输入的精度，

binom.logsf

和

binom.sf

都失败了，所以

log*

函数只是方便函数，即相当于代码

log(binom.sf())

问题
这是预期的行为还是我错过了配置步骤？

for n in range(n_min_max, n_min_max + 20):        
    p_trial_log = binom.logsf(n, N, p_let)
    p_trial = binom.sf(n, N, p_let)
    print(f"args = {n}, {N}, {p_let}")
    print(f"p_trial_log={p_trial_log:.1e}")
    print(f"p_trial=    {p_trial:.4e}")
    print(f"p_trial_exp={exp(p_trial_log):.4e}")
    print("")

输出：

args = 42, 100000, 4e-10
p_trial_log=-5.6e+02
p_trial=    1.2691e-242
p_trial_exp=1.2691e-242

args = 43, 100000, 4e-10
p_trial_log=-5.7e+02
p_trial=    1.1532e-248
p_trial_exp=1.1532e-248

args = 44, 100000, 4e-10
p_trial_log=-5.8e+02
p_trial=    1.0246e-254
p_trial_exp=1.0246e-254

args = 45, 100000, 4e-10
p_trial_log=-6.0e+02
p_trial=    8.9059e-261
p_trial_exp=8.9059e-261

args = 46, 100000, 4e-10
p_trial_log=-6.1e+02
p_trial=    7.5760e-267
p_trial_exp=7.5760e-267

args = 47, 100000, 4e-10
p_trial_log=-6.3e+02
p_trial=    6.3104e-273
p_trial_exp=6.3104e-273

args = 48, 100000, 4e-10
p_trial_log=-6.4e+02
p_trial=    5.1488e-279
p_trial_exp=5.1488e-279

args = 49, 100000, 4e-10
p_trial_log=-6.5e+02
p_trial=    4.1171e-285
p_trial_exp=4.1171e-285

args = 50, 100000, 4e-10
p_trial_log=-6.7e+02
p_trial=    3.2274e-291
p_trial_exp=3.2274e-291

args = 51, 100000, 4e-10
p_trial_log=-6.8e+02
p_trial=    2.4814e-297
p_trial_exp=2.4814e-297

args = 52, 100000, 4e-10
p_trial_log=-7.0e+02
p_trial=    1.8718e-303
p_trial_exp=1.8718e-303

args = 53, 100000, 4e-10
p_trial_log=-7.1e+02
p_trial=    1.3858e-309
p_trial_exp=1.3858e-309

args = 54, 100000, 4e-10
p_trial_log=-7.3e+02
p_trial=    1.0073e-315
p_trial_exp=1.0073e-315

args = 55, 100000, 4e-10
p_trial_log=-7.4e+02
p_trial=    7.2134e-322
p_trial_exp=7.2134e-322

args = 56, 100000, 4e-10
p_trial_log=-inf
p_trial=    0.0000e+00
p_trial_exp=0.0000e+00

args = 57, 100000, 4e-10
p_trial_log=-inf
p_trial=    0.0000e+00
p_trial_exp=0.0000e+00

args = 58, 100000, 4e-10
p_trial_log=-inf
p_trial=    0.0000e+00
p_trial_exp=0.0000e+00

args = 59, 100000, 4e-10
p_trial_log=-inf
p_trial=    0.0000e+00
p_trial_exp=0.0000e+00

args = 60, 100000, 4e-10
p_trial_log=-inf
p_trial=    0.0000e+00
p_trial_exp=0.0000e+00

args = 61, 100000, 4e-10
p_trial_log=-inf
p_trial=    0.0000e+00
p_trial_exp=0.0000e+00

Answer 1

追踪文档和[来源]，对于类

rv_generic

，

_logsf

与

_sf

和

相关

def _logsf(self, x, *args):
    with np.errstate(divide='ignore'):
        return log(self._sf(x, *args))

另请参阅

logsf

的源代码，它有一些额外的输入处理，但最终调用了上面的通用

rv_generic

。

Answer 2

_logsf

在

logsf()

或

rv_discrete

中实现，具体取决于您处理的是离散分布还是连续分布。

rv_continuous

进行一些错误检查，然后调用

logsf()

来实际计算概率。在某些 SciPy 发行版中，未实现

_logsf()

。在这些情况下，超类通过调用特定于分布的

_logsf()

（为所有分布定义）并取该数字的对数来实现它。对于非常小的概率，这会导致浮点下溢。对于其他分布，下溢风险已被识别并通过一些替代计算减轻。例如，正态分布有一个

_sf()

实现，可以避免取生存函数的对数。

这就是它的样子。

因为二项式分布的

类缺少

_logsf() 方法，所以它会调用 _logsf()

 并取该值的对数。

最后，欢迎补丁。如果您知道计算任何分布的记录生存函数的方法，该方法比当前所做的更准确，

您应该使用 SciPy 提出问题。