scipy integrate.quad возвращает неправильное значение - PullRequest
Новая
       148

scipy integrate.quad возвращает неправильное значение

0 голосов
/ 26 июня 2018

я использую scipy integrate.quad для расчета cdf нормального распределения: 

def nor(delta, mu, x):
    return 1 / (math.sqrt(2 * math.pi) * delta) * np.exp(-np.square(x - mu) / (2 * np.square(delta)))


delta = 0.1
mu = 0
t = np.arange(4.0, 10.0, 1)
nor_int = lambda t: integrate.quad(lambda x: nor(delta, mu, x), -np.inf, t)
nor_int_vec = np.vectorize(nor_int)

s = nor_int_vec(t)
for i in zip(s[0],s[1]): 
    print i

пока печатается следующим образом: 

(1.0000000000000002, 1.2506543424265854e-08)
(1.9563704110140217e-11, 3.5403445591955275e-11)
(1.0000000000001916, 1.2616577562700088e-08)
(1.0842532749783998e-34, 1.9621183122960244e-34)
(4.234531567162006e-09, 7.753407284370446e-09)
(1.0000000000001334, 1.757986959115912e-10)

для некоторого x возвращается значение, близкое к нулю, должно быть возвращено 1 может кто-нибудь сказать мне, что не так?

1 Ответ

0 голосов
/ 27 июня 2018

По той же причине, что и в , почему квадр возвращает оба ноля при интегрировании простого гауссовского pdf с очень малой дисперсией? , но, учитывая, что я не могу пометить его как дубликат, вот так:

Вы интегрируете функцию с точной локализацией (в дельте шкалы) в течение очень большого (фактически бесконечного) интервала. Процедура интегрирования может просто пропустить часть интервала, в которой функция существенно отличается от 0, вместо этого оценивая ее как 0. Некоторое руководство требуется. Параметр points может использоваться для этого (см. Связанный вопрос), но поскольку quad в течение бесконечного интервала не поддерживает его, интервал должен быть разделен вручную, например: 

for t in range(4, 10):
    int1 = integrate.quad(lambda x: nor(delta, mu, x), -np.inf, mu - 10*delta)[0] 
    int2 = integrate.quad(lambda x: nor(delta, mu, x), mu - 10*delta, t)[0] 
    print(int1 + int2)

Это печатает 1 или почти 1 каждый раз. Я выбрал mu-10*delta в качестве точки для разделения, полагая, что большая часть функции находится справа от нее, независимо от того, что такое mu и delta.

Примечания:

  1. Используйте np.sqrt и т. Д .; обычно нет смысла помещать math функции в код NumPy. Версии NumPy доступны и векторизованы.
  2. Применение np.vectorize к quad ничего не делает, кроме того, что делает код длиннее и немного сложнее для чтения. Используйте обычный цикл Python или понимание списка. См. NumPy векторизация с интеграцией
...