my question concerns the logp function in the PYMC package in Python. Ultimately I want to calculate a quantity that goes by many names, namely the Bayes-factor/ evidence/ marginal-likelihood of the model, but perhaps first demystifying logp will get me closer.

Let us assume that I have instantiated the model below

mc = MCMC(myModel)

and have then done my sampling.

What does calling the logp exactly return?


Is it the log-posterior for the last sample? and if so, how can I access the log-posterior value for other samples generated by the MCMC process? (Any hints on how to use this for calculating the marginal likelihood welcome!!!)

Thank you all for your time, N.


Calling myModel.logp returns the log probability of the model at the current parameter setting. So if you ran MAP() it would be the logp at the MAP, if you sampled, it would be the logp at the last sampled parameter value. Thus, these are not the values you're looking for :).

Unfortunately, computing the marginal probability of a model is pretty difficult because you have to integrate out all parameters. Often, this can only be done analytically (which works for certain conjugate models).

There are a few ways to approximate this from the trace of a model using the harmonic mean, but there are some serious issues with this approach. See Radford Neal's blog post for a good description of the issues and the comments section for ideas on how to remedy them.

Personally, I think Bayes Factors for model comparison are more trouble than they are worth. I often use the DIC, despite it's limitations, or posterior predictive checks.

| cite | improve this answer | |
  • $\begingroup$ Thanks twiecki, you confirm my suspicion about logp. I am aware that map estimation is a poor man's Bayes estimation as people say. I will try to overcome the problem in some other way... cheers. $\endgroup$ – ngiann Jan 28 '15 at 14:10

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.