Dynamic Bayesian Network Inference

You could try something on these lines:

with pm.Model():
    n = len(observed)
    v1 = pm.Normal('v1', 0, 1, shape=(n,))
    v2 = pm.Normal('v2', v1, 1, shape=(n,))
    p = pm.Normal('p', 0, 1)
    v3 = pm.AR('v3', [v2[1:], p*tt.ones(n - 1)], 1, init=v2[0], observed=observed)

You would need to parametrize the model better to avoid divergences, but as a toy example to build on, it could maybe work.