`model.initial_point()` returns inf if using `pm.distributions.transforms.Ordered()`

import numpy as np
import pymc as pm

N = 2
y = np.array([0.0, 1.0])
with pm.Model() as pymc_model:
    # Priors
    mu = pm.Normal(
        "mu",
        mu=0,
        sigma=2,
        shape=2,
        transform=pm.distributions.transforms.Ordered(),
    )
    sigma = pm.HalfNormal("sigma", sigma=2, shape=2)
    theta = pm.Beta("theta", alpha=5, beta=5)

    # Mixture likelihood
    mix = pm.NormalMixture(
        "mix", w=[theta, 1 - theta], mu=mu, sigma=sigma, observed=y
    )

seed = 42
np.random.seed(seed)
pymc_model.initial_point(random_seed=seed)

It returns:

{'mu_ordered__': array([  0., -inf]),
 'sigma_log__': array([0.69314718, 0.69314718]),
 'theta_logodds__': array(0.)}

PyMC version: 5.12.0

The default initial point for normal (i guess [0, 0]) is not possible under the ordered transformation. You need to pass an initval that respect the ordering constraint. For example:

    mu = pm.Normal(
        "mu",
        mu=0,
        sigma=2,
        shape=2,
        transform=pm.distributions.transforms.Ordered(),
        initval=[-1, 1]
    )
1 Like

Thanks! pm.sample(1000, model=pymc_model, initvals = {"mu": np.array([0.0, 1.0])}) also works when sampling, but I think specifying initval in the way you suggested is better!

1 Like