jak przyspieszyć model Markko PyMC?

Question

Czy istnieje sposób na przyspieszenie tego prostego modelu PyMC? W przypadku 20-40 punktów danych dopasowanie zajmuje około 5-11 sekund.

import pymc 
import time 
import numpy as np 
from collections import OrderedDict 

# prior probability of rain 
p_rain = 0.5 
variables = OrderedDict() 
# rain observations 
data = [True, True, True, True, True, 
     False, False, False, False, False]*4 
num_steps = len(data) 
p_rain_given_rain = 0.9 
p_rain_given_norain = 0.2 
p_umbrella_given_rain = 0.8 
p_umbrella_given_norain = 0.3 
for n in range(num_steps): 
    if n == 0: 
     # Rain node at time t = 0 
     rain = pymc.Bernoulli("rain_%d" %(n), p_rain) 
    else: 
     rain_trans = \ 
      pymc.Lambda("rain_trans", 
         lambda prev_rain=variables["rain_%d" %(n-1)]: \ 
         prev_rain*p_rain_given_rain + (1-prev_rain)*p_rain_given_norain) 
     rain = pymc.Bernoulli("rain_%d" %(n), p=rain_trans) 
    umbrella_obs = \ 
     pymc.Lambda("umbrella_obs", 
        lambda rain=rain: \ 
        rain*p_umbrella_given_rain + (1-rain)*p_umbrella_given_norain) 
    umbrella = pymc.Bernoulli("umbrella_%d" %(n), p=umbrella_obs, 
           observed=True, 
           value=data[n]) 
    variables["rain_%d" %(n)] = rain 
    variables["umbrella_%d" %(n)] = umbrella 

print "running on %d points" %(len(data)) 
all_vars = variables.values() 
t_start = time.time() 
model = pymc.Model(all_vars) 
m = pymc.MCMC(model) 
m.sample(iter=2000) 
t_end = time.time() 
print "\n%.2f secs to run" %(t_end - t_start) 

z punktów tylko 40 danych, potrzeba 11 sekund, aby uruchomić:

running on 40 points 
[-----------------100%-----------------] 2000 of 2000 complete in 11.5 sec 
11.54 secs to run 

(z 80 punktów w ciągu 20 sekund). To jest zabawny przykład. Wyrażenia wewnątrz Lambda(), które określają przejścia, są w praktyce bardziej złożone. Ta podstawowa struktura kodu jest elastyczna (podczas gdy kodowanie modelu z macierzami przejścia jest mniej elastyczne). Czy istnieje sposób na utrzymanie podobnej struktury kodu, ale lepszą wydajność? W razie potrzeby chętnie przełączymy się na PyMC3. Dzięki.

Answer 1

Łańcuch Markov Monte Carlo jest znanym sekwencyjnym problemem.

Oznacza to, że czas jej działania jest proporcjonalny do liczby kroków i czasu wykonywania funkcji fitness.

Istnieje kilka sztuczek można zrobić, jednakże:

• Zastosowanie pypy (wymaga przerobienia, pymc nie są obsługiwane)
• Zastosowanie Gibbsa próbkowanie poprawić następnym krokiem
• wykorzystać punkty wielokrotny start (w równoległy)
• stosowanie wielu oddziałów (równolegle)
• użytku heurystyczne do zatrzymania łańcucha wcześniej
• Zastosowanie przybliżenie f lub punkty, które są zbliżone do już obliczone

Harder podejścia:

• Zastosowanie Numba (kompiluje funkcji zdolności do kodu maszynowego)
• przepisać swoją funkcję biznesowe w C (lub podobny)
• użytku natywny kod MCMC (inny niż Python, wymaga powyższego)

Wreszcie istnieje wiele badań:

http://www.mas.ncl.ac.uk/~ndjw1/docs/pbc.pdf

https://sites.google.com/site/parallelmcmc/

http://pyinsci.blogspot.com/2010/12/efficcient-mcmc-in-python.html (pypy)