Pm.sample_posterior_predictive() from an AR1 process

Questions
1

I’m trying to fit an AR1 a time-series as following:

with pm.Model() as ar1_mat1:
    theta = pm.Normal('theta',mu= 0, sigma = 3.0)
    tau = pm.HalfNormal('tau',sigma =1)
    
    likelihood = pm.AR1('likelihood', k=theta, tau_e = tau, observed = y_standard.values)
    
    trace = pm.sample(1000, tune=2000, init="advi+adapt_diag", random_seed=2,return_inferencedata=True)

    poster = pm.sample_posterior_predictive(trace)

Question is, the posterior predictive samples all have the same values for likelihood (the observed data)… Is this expected or is there something i’m missing?
Standardized “y” looks like this:

image
image899×683 37.2 KB

and the likelihood is as follows:

{'likelihood': array([[ 0.85319392, -0.14915978,  0.82250963, ..., -0.17984407,
         -0.8855829 , -0.895811  ],
        [ 0.85319392, -0.14915978,  0.82250963, ..., -0.17984407,
         -0.8855829 , -0.895811  ],
        [ 0.85319392, -0.14915978,  0.82250963, ..., -0.17984407,
         -0.8855829 , -0.895811  ],
        ...,
        [ 0.85319392, -0.14915978,  0.82250963, ..., -0.17984407,
         -0.8855829 , -0.895811  ],
        [ 0.85319392, -0.14915978,  0.82250963, ..., -0.17984407,
         -0.8855829 , -0.895811  ],
        [ 0.85319392, -0.14915978,  0.82250963, ..., -0.17984407,
         -0.8855829 , -0.895811  ]])}
2

Also, is there a cleaner way to make future predictions using AR1/AR in pymc3?

3

Where/how did you get likelihood you showed?

4

It’s the only variable that’s there in poster = pm.sample_posterior_predictive(trace)

5

Sorry, I thought you were doing prior prediction.

How is sampling going and what do your parameter estimates look like? Are you getting divergences? I suspect there may be some sampling issues and that you are either just getting your own observed data back or some other degenerate behavior. But that’s just a guess.

6

There actually are divergences (7-8) per chain. I think I might be instantiating the AR1 wrong. In any case, how might I generate y values from the learned parameters? Is it something like:

with ar1_mat1:
    posterior_samples = pm.sample_posterior_predictive(trace,var_names=['tau','theta','center'])

and then y:

y[t] = y[t-1]* theta + np.random.normal(loc = center, scale = np.sqrt(1/tau))
1 Like
7

If you want to generate posterior predictions “by hand”, then you don’t need sample_posterior_predictive(). You can just unpack parameter values from the posterior returned by pm.sample(). Then plug them back into your model as you describe. If you haven’t already, you should check out the AR1 notebook.