Time Varying Overdispersed Poisson Process

You are right, we don’t need that pm.Normal(). The below works for me and is much quicker, taking just a few minutes to sample.

with pm.Model() as m2:
    #a = pm.Gamma('a', 5, 1.5)
    #b = pm.Beta('b', 8, 41)
    
    g, A, s = trend_model(m2, df1['t'], changepoints_prior_scale=0.2, n_changepoints=4)
    x_monthly, beta_monthly = seasonality_model(m2, df1, 'yearly')
    x_daily, beta_daily = seasonality_model(m2, df1, 'weekly')
    x_hourly, beta_hourly = seasonality_model(m2, df1, 'daily')

    y = pm.Deterministic('y', g + det_dot(x_monthly, beta_monthly) + det_dot(x_daily, beta_daily) + det_dot(x_hourly, beta_hourly))
    
    #sigma = pm.HalfCauchy('sigma', 10.)
    #seasonality = pm.Normal('seasonality', mu=y, sigma=sigma, shape=len(df1))
    
    #lam = pm.Gamma('lam', a, b)
    lam = pm.HalfNormal('lam', 20.)
    obs = pm.Poisson('obs', lam + y, observed=df1.y)
    
    trace = pm.sample(1000)
1 Like