Ошибка точности в методе `scipy.stats.binom`

Question

Я использую scipy.stats.binom для работы с биномиальным распределением. Учитывая n и p, функция вероятности равна

Сумма от k в диапазоне от 0 до n должна (и действительно дает) 1. Зафиксируя точку x_0, мы можем добавить вероятности в обоих направлениях, и две суммы следует добавить к 1. Однако приведенный ниже код выдает два разных ответа, когда x_0 близко к n.

from scipy.stats import binom

n = 9
p = 0.006985
b = binom(n=n, p=p)
x_0 = 8

# Method 1
cprob = 0
for k in range(x_0, n+1):
    cprob += b.pmf(k)
print('cumulative prob with method 1:', cprob)

# Method 2
cprob = 1
for k in range(0, x_0):
    cprob -= b.pmf(k)
print('cumulative prob with method 2:', cprob)

Я ожидаю, что результаты обоих методов будут согласованы. Для x_0 < 7 согласен, но для x_0 >= 8, как указано выше, я получаю

>> cumulative prob with method 1: 5.0683768775504006e-17
>> cumulative prob with method 2: 1.635963929799698e-16

Ошибка точности в двух методах распространяется через мой код (позже) и дает совершенно разные ответы. Любая помощь приветствуется.

Answer 1

Ошибки округления порядка станка эпсилон ожидаются и являются неизбежными.То, что они размножаются, а потом взрываются, означает, что ваша проблема очень плохо обусловлена.Вам нужно переосмыслить алгоритм или реализацию, в зависимости от того, откуда взялась плохая обусловленность.

В вашем конкретном примере вы можете получить, используя либо np.sum (который пытается быть осторожным с округлением), либо даже math.fsum из стандартной библиотеки.