How to bound a Dirichlet?

Questions
1

If I wanted a Dirichlet prior, but I wanted to censor out any samples with values under say 1e-14 or so (to help with numerical stability, as p(x)>infinity messes up find_MAP), how would I do that? This sort of code doesn’t really seem to work (find_MAP gives values outside of the simplex, or just crashes).

bounder = pm.Bound(pm.Dirichlet,lower=np.ones(shape=(arrsize))*1e-5,upper = 1-np.ones(shape=(arrsize))*1e-15)
prior = bounder('a',a = a_vec,shape=(arrsize),testval = testval_arr)

Or is using logratios or something a better approach here?

Thanks,
Abhishek.

2

Something on top of my head is to add a Potential:

    prior = pm.Dirichlet('a', a=a_vec,shape=(arrsize), testval=testval_arr)
    # ensure a is not too small
    bounding_penalty = pm.Potential('p_min_potential', tt.switch(a < 1e-5, -np.inf, 0))
3

It’s giving logp as -inf though, is that ok?
logp = -inf, ||grad|| = 115.79: 100%|██████████| 4/4 [00:00<00:00, 548.83it/s]

using VI is also just oscillation between -inf and NaN.

4

It means that a goes below 1e-5. You should try sampling, which respect the Potential penalty more fatefully.

5

ah ha it works, Thanks!

Pl. close topic.

1 Like