Why does the error occur in the prior prediction but not in the posterior prediction?

Questions v5
1

I was working on a program to estimate the prior and posterior predictive distributions with the following model, but the error occurred only in the former!

# === Prior ===
model = pm.Model() # define a context
with model:
    omega = pm.Beta("omega", alpha=1, beta=1)
    phi = pm.HalfCauchy("phi", beta=1)

    alpha = pm.Deterministic("alpha", omega * (phi - 2) + 1)
    beta = pm.Deterministic("beta", (1 - omega) * (phi - 2) + 1)

    theta = pm.Beta("theta", alpha=alpha, beta=beta)

    prior_samples = pm.sample_prior_predictive(random_seed=42)

print(prior_samples["prior"]["beta"].values)

# === Posterior ===
# prepare dummy data
np.random.seed(42)
data = np.random.binomial(n=1, p=0.4, size=10)

model = pm.Model() # define a context
with model:
    X = pm.Data("X", data)

    omega = pm.Beta("omega", alpha=1, beta=1)
    phi = pm.HalfCauchy("phi", beta=1)

    alpha = pm.Deterministic("alpha", omega * (phi - 2) + 1)
    beta = pm.Deterministic("beta", (1 - omega) * (phi - 2) + 1)

    theta = pm.Beta("theta", alpha=alpha, beta=beta)

    obs = pm.Binomial("obs", n=1, p=theta, observed=X)
    
# sampling
with model:
    idata = pm.sample(chains=2, random_seed=42)

print(idata["posterior"]["theta"].values)

and got error,

Traceback (most recent call last):
  File "/home/user/bayes_inference/test.py", line 38, in <module>
    prior_samples = pm.sample_prior_predictive(random_seed=42)
  File "/home/user/bayes_inference/venv/lib/python3.10/site-packages/pymc/sampling/forward.py", line 417, in sample_prior_predictive
    values = zip(*(sampler_fn() for i in range(samples)))
  File "/home/user/bayes_inference/venv/lib/python3.10/site-packages/pymc/sampling/forward.py", line 417, in <genexpr>
    values = zip(*(sampler_fn() for i in range(samples)))
  File "/home/user/bayes_inference/venv/lib/python3.10/site-packages/pytensor/compile/function/types.py", line 983, in __call__
    raise_with_op(
  File "/home/user/bayes_inference/venv/lib/python3.10/site-packages/pytensor/link/utils.py", line 523, in raise_with_op
    raise exc_value.with_traceback(exc_trace)
  File "/home/user/bayes_inference/venv/lib/python3.10/site-packages/pytensor/compile/function/types.py", line 970, in __call__
    self.vm()
  File "/home/user/bayes_inference/venv/lib/python3.10/site-packages/pytensor/graph/op.py", line 515, in rval
    r = p(n, [x[0] for x in i], o)
  File "/home/user/bayes_inference/venv/lib/python3.10/site-packages/pytensor/tensor/random/op.py", line 330, in perform
    smpl_val = self.rng_fn(rng, *([*args, size]))
  File "/home/user/bayes_inference/venv/lib/python3.10/site-packages/pytensor/tensor/random/op.py", line 127, in rng_fn
    return getattr(rng, self.name)(*args, **kwargs)
  File "numpy/random/_generator.pyx", line 400, in numpy.random._generator.Generator.beta
  File "_common.pyx", line 612, in numpy.random._common.cont
  File "_common.pyx", line 427, in numpy.random._common.check_constraint
ValueError: a <= 0
Apply node that caused the error: beta_rv{0, (0, 0), floatX, True}(RandomGeneratorSharedVariable(<Generator(PCG64) at 0x7F45331FA340>), [], 11, alpha, beta)
Toposort index: 3
Inputs types: [RandomGeneratorType, TensorType(int64, shape=(0,)), TensorType(int64, shape=()), TensorType(float64, shape=()), TensorType(float64, shape=())]
Inputs shapes: ['No shapes', (0,), (), (), ()]
Inputs strides: ['No strides', (0,), (), (), ()]
Inputs values: [Generator(PCG64) at 0x7F45331FA340, array([], dtype=int64), array(11), array(-0.717), array(0.829)]
Outputs clients: [['output'], ['output']]

HINT: Re-running with most PyTensor optimizations disabled could provide a back-trace showing when this node was created. This can be done by setting the PyTensor flag 'optimizer=fast_compile'. If that does not work, PyTensor optimizations can be disabled with 'optimizer=None'.
HINT: Use the PyTensor flag `exception_verbosity=high` for a debug print-out and storage map footprint of this Apply node.

When I read this, I was convinced that the parameters of the Beta distribution must be non-negative, so the error would occur in the current code where omega and phi are generated from Beta and HalfCauchy distributions, respectively.

But why is it that the posterior predictive distribution does not have a similar error and outputs the end result?

Thanks for answering my question.

version-related information is here

Python==3.10.12
pymc==5.14.0
numpy==1.26.4
2

My guess is that theta is the problem. There’s nothing to prevent the quantity phi - 2 from being negative, since a Half-Cauchy can generate values on the interval [0, 2) just fine.

3

Yeah, I have learned and understood that, but why does the above error occur only in the case of a prior predictive distribution and not in the case of a posterior predictive distribution?

4

The negative values imply a logp of -inf, so the sampler never accepts those draws

1 Like
5

A self contained example of this is:

import pymc as pm

with pm.Model() as m:
  sigma = pm.Normal("sigma", 0.001, 1)  # bad prior for sigma
  y = pm.Normal("y", mu=0, sigma=sigma, observed=[-1, 0, 1])

with m:
  try:
      pm.sample_prior_predictive()
  except Exception as exc:
     print(f"{type(exc).__name__}: {exc}")  # ValueError: scale < 0

  idata = pm.sample()
  pm.sample_posterior_predictive(idata)  # Fine

Posterior predictive works because there are only positive draws of sigma, whereas in the unconstrained prior predictive, negative draws are suggested (roughly half the time)

2 Likes