In [1]:
!date
Wed Dec 18 13:54:39 PST 2013
In [5]:
import pymc as pm, theano.tensor as T
In [6]:
x = array([  0,   1,   5,  10,  15,  20,  25,  30,  35,  40,  45,  50,  55,
            60,  65,  70,  75,  80,  85,  90,  95, 100, 105])

y = array([-5.31126213, -6.88284349, -7.28148079, -7.20912457, -6.06006241,
           -5.69987917, -5.72478151, -5.62202549, -5.36570549, -4.96331167,
           -4.50282001, -3.99181652, -3.44459009, -2.88168406, -2.35241652,
           -1.82025242, -1.25903034, -0.66321015,  0.06458783,  0.87678754,
            1.70916784,  2.56996703,  3.38351035])
In [7]:
with pm.Model() as m:
    a = pm.Flat('a')
    b = pm.Flat('b')
    c = pm.Flat('c')
    d = pm.Flat('d')
    e = pm.Flat('e')
    f = pm.Flat('f')
    g = pm.Flat('g')
    h = pm.Flat('h')
    
    t1 = a**((x+b)**c)
    t2 = d * T.exp(-e * T.log(x/f)**2)
    t3 = g*h**x
    
    y_pred = t1 + t2 + t3
    y_obs = pm.Normal('y_obs', mu=y_pred/y, sd=1., observed=ones_like(y))
In [8]:
with m: trace = pm.sample(20000, pm.Metropolis())
 [-----------------100%-----------------] 20000 of 20000 complete in 12.9 sec
In [ ]:
with m: trace = pm.sample(20000, pm.NUTS())
In [ ]: