Как пробовать неоднородные пуассоновские процессы в Python быстрее этого?

Question

Я отбираю процесс Пуассона в миллисекундном масштабе времени, где скорость не фиксирована. Я дискретирую процесс выборки, проверяя каждый интервал дельта размера, существует ли там событие или нет на основе средней скорости в этом интервале. Поскольку я использую Python, он работает немного медленнее, чем я надеялся. Код настоящее время я использую следующий:

import numpy 
def generate_times(rate_function,max_t,delta): 
    times = [] 
    for t in numpy.arange(delta,max_t,delta): 
     avg_rate = (rate_function(t)+rate_function(t+delta))/2.0 
     if numpy.random.binomial(1,1-math.exp(-avg_rate*delta/1000.0))>0: 
      times.extend([t]) 
    return times 

функция скорости может быть произвольным, я не ищу для закрытой форме раствора заданной функцией скорости.

Если вы хотите, чтобы некоторые параметры, чтобы играть с вами, можно попробовать:

max_t = 1000.0 
delta = 0.1 
high_rate = 100.0 
low_rate = 0.0 
phase_length = 25.0 
rate_function = (lambda x: low_rate + (high_rate-low_rate)*0.5*(1+math.sin(2*math.pi*x/phase_length))) 

Answer 1

Вот версия, которая работает около 75x быстрее и реализует ту же функцию:

def generate_times_opt(rate_function,max_t,delta): 
    t = numpy.arange(delta,max_t, delta) 
    avg_rate = (rate_function(t) + rate_function(t + delta))/2.0 
    avg_prob = 1 - numpy.exp(-avg_rate * delta/1000.0) 
    rand_throws = numpy.random.uniform(size=t.shape[0]) 

    return t[avg_prob >= rand_throws].tolist() 

Выход:

11:59:07 [70]: %timeit generate_times(rate_function, max_t, delta) 
10 loops, best of 3: 75 ms per loop 

11:59:23 [71]: %timeit generate_times_opt(rate_function, max_t, delta) 
1000 loops, best of 3: 1.15 ms per loop 

---

Sidenote: Это не лучший способ имитировать процесс Inhomogenous Poisson. От Wikipedia:

Неоднородный пуассоновский процесс с функцией интенсивности Х (г) может быть смоделирован отказом выборки из однородного процесса Пуассона с Й фиксированной скоростью: выбрать достаточно большой А так, что А (Т) = Х р (t) и имитировать пуассоновский процесс с параметром скорости λ. Примите событие от моделирования Пуассона в момент времени t с вероятностью p (t).