turion / dunai-bayes

BSD 3-Clause "New" or "Revised" License
3 stars 0 forks source link

The implementation of bayesfilter #3

Open reubenharry opened 1 year ago

reubenharry commented 1 year ago

Yes, I think I need a bit of guidance here what typical nomenclature and conventions are. It's certain that we need some kind of importance function (is that the right name?), because exact equality immediately looses way to much probability mass. To make this configurable, I defined the SoftEq class, so one can configure the importance function via newtypes, while having hopefully sensible defaults for Double and so on. But I'm fuzzy on where it belongs. I used it here in the bayes filter, but that's a bit weird. I would have expected to use it in the particle filter directly (somewhere close to runPopulation), but couldn't find the right place.

So this turns out to be a common problem in probabilistic programming languages related to Borel's paradox (which is the issue of conditioning on measure-0 sets). If I write:

example = do
  x <- prior
  o <- normal x 1
  condition (o == 1)
  return x

there's no hope that any inference algorithm will get the true posterior (without doing a source code transformation) because the probability of drawing a given number from a normal (or any continuous distribution) is 0. Instead we could write

example = do
  x <- prior
  factor (normalPdf x 1 o)
  return x

This works fine, and should be (modulo typos/thinkos) denotationally equivalent.

My concern about the SoftEq solution is that it distorts the probabilistic program that specifies the posterior (if indeed bayesFilter is used to define the posterior). As you say above, it should not really be in the probabilistic program, but rather in the inference somehow, because otherwise the program ceases to be a correct mathematical description of the posterior (unless we're lucky).

My solution in general has been to avoid doing something like SoftEq, and instead specify the model with a data structure like:

data Bayesian m z o = Bayesian
  { prior :: m z,
    generative :: z -> m o
    likelihood :: z -> o -> Log Double 
  }

where generative and likelihood correspond to the same conditional distribution, but in two forms. The second form is the one used in writing posterior:

posterior :: (MonadInfer m, Foldable f, Functor f) => Bayesian m z o -> f o -> m z
posterior Bayesian {..} os = do
  z <- prior
  factor $ product $ fmap (likelihood z) os
  return z

which is the non-reactive analog of bayesFilter.

On the other hand, it's sometimes not possible to give an explicit density for the generative model, so we could offer something like bayesFilter as long as it's clear it's not giving the true posterior.

turion commented 1 year ago

@reubenharry I had an idea how to tackle this problem, extending your Bayesian data structure:

https://github.com/turion/monad-bayes/blob/0d0422a2badc58c9315ed8127c74c1d74c96f41e/src/Control/Monad/Bayes/WithPDF.hs

reubenharry commented 1 year ago

Just at a surface level read, I'm wondering whether the ideas in https://github.com/psg-mit/probzelus-haskell/tree/master/haskell/src are relevant. They are a Haskell port of ideas from the Gen probabilistic programming language, which also passes along pdfs in order to deal with the Borel conditioning issue, I think in a similar way to what you are suggesting. See in particular Distributions.hs, Metaprob.hs and MetaGen.hs.

turion commented 1 year ago

Good catch!

https://github.com/psg-mit/probzelus-haskell/blob/a4b66631451b6156938a9c5420cfff2999ecbbc6/haskell/src/Distributions.hs#L20

This is pretty similar to my idea, but a bit purer even. The rest I must admit I don't fully understand, probably there are some good ideas in there that I can't appreciate yet.