Question about sampling from prior predictive with time series distributions

Questions
1

Hello folks,

I’m attempting to sample from the prior predictive distribution of a time series model. The original model is rather complex but here is a basic example of what can go wrong.

Simple model:

T = 10

with pm.Model() as test_model:
    alpha = pm.Normal('alpha', 0, 1, shape=(4,))
    ar = pm.distributions.timeseries.AR('ar', alpha, 0.1, shape=(T,),
                                       init=pm.HalfNormal.dist(0.1))
    
with test_model:
    prior_pred = pm.sample_prior_predictive()

Error:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-97-0052735bd0c9> in <module>
      7 
      8 with test_model:
----> 9     prior_pred = pm.sample_prior_predictive()

/anaconda3/lib/python3.7/site-packages/pymc3/sampling.py in sample_prior_predictive(samples, model, vars, var_names, random_seed)
   1355     names = get_default_varnames(vars_, include_transformed=False)
   1356     # draw_values fails with auto-transformed variables. transform them later!
-> 1357     values = draw_values([model[name] for name in names], size=samples)
   1358 
   1359     data = {k: v for k, v in zip(names, values)}

/anaconda3/lib/python3.7/site-packages/pymc3/distributions/distribution.py in draw_values(params, point, size)
    422         while to_eval or missing_inputs:
    423             if to_eval == missing_inputs:
--> 424                 raise ValueError('Cannot resolve inputs for {}'.format([str(params[j]) for j in to_eval]))
    425             to_eval = set(missing_inputs)
    426             missing_inputs = set()

ValueError: Cannot resolve inputs for ['ar']

Another example:

T = 10

with pm.Model() as test_model:
    alpha = pm.Normal('alpha', 0, 1, shape=(4,))
    grw = pm.distributions.timeseries.GaussianRandomWalk('grw',
                                                        mu=alpha,
                                                        sd=0.1,
                                                        shape=(T, 4))
with test_model:
    prior_pred = pm.sample_prior_predictive()

The identical error:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-116-01b948d7ee4a> in <module>
      8                                                         shape=(T, 4))
      9 with test_model:
---> 10     prior_pred = pm.sample_prior_predictive()

/anaconda3/lib/python3.7/site-packages/pymc3/sampling.py in sample_prior_predictive(samples, model, vars, var_names, random_seed)
   1355     names = get_default_varnames(vars_, include_transformed=False)
   1356     # draw_values fails with auto-transformed variables. transform them later!
-> 1357     values = draw_values([model[name] for name in names], size=samples)
   1358 
   1359     data = {k: v for k, v in zip(names, values)}

/anaconda3/lib/python3.7/site-packages/pymc3/distributions/distribution.py in draw_values(params, point, size)
    422         while to_eval or missing_inputs:
    423             if to_eval == missing_inputs:
--> 424                 raise ValueError('Cannot resolve inputs for {}'.format([str(params[j]) for j in to_eval]))
    425             to_eval = set(missing_inputs)
    426             missing_inputs = set()

ValueError: Cannot resolve inputs for ['grw']

It is rather likely that the error is with me, but I do not know what I am doing wrong. I have read Error with sample_prior_predictive and know that there’s some dev work going on with draw_values but am unsure why this would break the above model but not a non-time-series model. For example,

T = 10

with pm.Model() as test_model:
    alpha = pm.Normal('alpha', 0, 1, shape=(4,))
    not_ar = pm.Normal('not_ar', alpha, 0.1, shape=(T, 4))
    
with test_model:
    prior_pred = pm.sample_prior_predictive()

works just fine.

Another (non-time series) example that works fine:

T = 10

with pm.Model() as test_model:
    alpha = pm.Normal('alpha', 0, 1, shape=(4,))
    not_ar = pm.MvNormal('not_ar', alpha, cov=np.eye(4), shape=(T, 4))
    
with test_model:
    prior_pred = pm.sample_prior_predictive()

Any advice would be quite welcome. I get around this now by writing my own prior predictive sampler in numpy but would like to avoid doing this if possible.

Thank you all in advance!

-Dave

EDIT: silly me, the time series methods don’t define a .random method. I’ll write up at least some of them and submit a pull request eventually. Thanks for letting me rubber-duck!

1 Like