How to sample efficiently from time series data?

Something along the line of:

s0, e0, i0 = 100., 50., 25.
st0, et0, it0 = [theano.shared(x) for x in [s0, e0, i0]]

C = np.array([3, 5, 8, 13, 21, 26, 10, 3], dtype=np.float64)
D = np.array([1, 2, 3, 7, 9, 11, 5, 1], dtype=np.float64)

def seir_one_step(st0, et0, it0, beta, gamma, delta):
    bt0 = st0 * beta
    ct0 = et0 * gamma
    dt0 = it0 * delta
    
    st1 = st0 - bt0
    et1 = et0 + bt0 - ct0
    it1 = it0 + ct0 - dt0
    return st1, et1, it1

with pm.Model() as model:
    beta = pm.Beta('beta', 2, 10)
    gamma = pm.Beta('rho', 2, 10)
    delta = pm.Beta('gamma', 2, 10)

    (st, et, it), updates = theano.scan(
        fn=seir_one_step,
        outputs_info=[st0, et0, it0],
        non_sequences=[beta, gamma, delta],
        n_steps=len(C))

    ct = pm.Binomial('c_t', et, gamma, observed=C)
    dt = pm.Binomial('d_t', it, delta, observed=D)

    trace = pm.sample()

with model:
    bt = pm.Binomial('b_t', st, beta, shape=len(C))
    ppc_trace = pm.sample_posterior_predictive(trace, var_names=['b_t'])

If you are interested in SEIR like model, Prediction of Danish Covid19 cases is a great example.

2 Likes