How to model very complicated likelihood automatically？

[NFYG]

Let’s think about a likelihood which is non-trivial

y=(x+\mu+\sigma\epsilon)^2,\ \text{with}\ \epsilon\sim\mathcal{N}(0,1)

I know this can be transformed into ChiSquared dist analytically.

But we still will encounter more complicated and non-linear function (maybe deep nerual network!!) involving many different random variables. The complication make it impossible to design a custom DensityDist.

Is there any way to model a likelihood automatically?

[ricardoV94]

Trivial stuff PyMC CustomDist can do it for you. Borrowing your example:

import pymc as pm

def dist(mu, sigma, size):
  return (mu + sigma * pm.Normal.dist(size=size)) ** 2

mu = 2
sigma = 1.5
x1 = pm.CustomDist.dist(mu, sigma, dist=dist)

value = 0.7
print(pm.logp(x1, value).eval())  # -1.9362354720843369

If your complicated likelihood function is still relatively simple (mostly chained invertible functions), PyMC can figure it for you automatically. Otherwise you’ll have to step in.

But to take a step back, many times you don’t really need it. You can define an arbitrarily complex function for the mean / standard deviation and then assume a normal likelihood (or some other simple form) around this. Whether this is appropriate depends on your problem obviously.

[NFYG]

Finally, I choose to use approximate bayesian computation (ABC), like using kernel density estimate (KDE) to calculate likelihood approximately.