Multivariate LogNormal distribution

CustomDist should be fine.

If you are not in a hurry this will actually be trivial after Infer logp for elemwise transformations of multivariate variables by ricardoV94 · Pull Request #6607 · pymc-devs/pymc · GitHub get’s merged, probably in the next release.

In that case all you will need to do is something like

def dist(mu, cov, size):
  return pm.math.exp(pm.MvNornal.dist(mu, cov, size=size))

with pm.Model():
  mu = ...
  cov = ...
  x = pm.CustomDist("x", mu, cov, dist=dist)

PyMC will be able to figure out the logp automatically from the dist function. For now you would still implement a logp function.