Truncated/Bounded Multivariate Normal?

Thanks a lot. In principle, this works.

However, I realized some undesired effect the introduction of a potential has on the posterior distribution. Even if I specify a trivial potential

muBounded = pm.Potential("muBounded", tt.switch(1, 0, 0))

the sampled posterior distribution changes. Is this expected, and if so, why?

Here’s the model specification with potential:

with pm.Model() as update_model:
    # non-centered parametrization
    mean_error = pm.Normal('mean_error', 0, 1, shape=dimension)
    mu = pm.Deterministic('mu', mu_marginal_mean + tt.dot(mu_cholesky, mean_error))
    muBounded = pm.Potential("muBounded", tt.switch(1, 0, 0))

    L_packed_error = pm.Normal('L_packed_error', 0, 1, shape=packedMatrixDimension[dimension])
    L_packed = pm.Deterministic('L_packed', L_packed_marginal_mean + tt.dot(L_packed_cholesky, L_packed_error))

    obs = pm.MvNormal('obs', muBounded, chol=pm.expand_packed_triangular(dimension, L_packed), observed=observations)

This is the original model without potential:

with pm.Model() as update_model:
    # non-centered parametrization
    mean_error = pm.Normal('mean_error', 0, 1, shape=dimension)
    mu = pm.Deterministic('mu', mu_marginal_mean + tt.dot(mu_cholesky, mean_error))

    L_packed_error = pm.Normal('L_packed_error', 0, 1, shape=packedMatrixDimension[dimension])
    L_packed = pm.Deterministic('L_packed', L_packed_marginal_mean + tt.dot(L_packed_cholesky, L_packed_error))

    obs = pm.MvNormal('obs', mu, chol=pm.expand_packed_triangular(dimension, L_packed), observed=observations)

As you can see from the posterior distribution shown below (see second plot), the sampled distribution probabilities no longer center around the “true” value, which are represented by black vertical lines. Interestingly, it only affects the variable that is not constrained directly (L_packed, and not mu to which I apply the bound).

Figure_1_bounded

The full code of the MWE, as well as the corresponding plot of the model without potential is available in this thread: Updating multivariate priors

I’d be very curious on any opinion on this!