kthohr / mcmc

A C++ library of Markov Chain Monte Carlo (MCMC) methods
https://mcmclib.readthedocs.io/en/latest/
Apache License 2.0
166 stars 45 forks source link

about the log density kernel #4

Closed sanjeevmk closed 2 years ago

sanjeevmk commented 5 years ago

Hi, Thanks for the library, its incredibly useful!

I've some doubts about what to pass as my energy functional. I have an energy function with respect to some P, say; lets call it f(P) . For my application f(P) = 0 should correspond to high probability samples, while higher energies should correspond to low probability ones. So I use e^( -1 f(P) ) as my density function. Since f(p) can lie in [0,infinity), its not technically a density function. To use your HMC API, my log posterior kernel returns the log of this density, i.e, -1 f(P) , and the corresponding gradient with respect to P. The HMC run is initialized with some some initial values for P, and the function that computes the log density, receives proposals for P.

I sample a 100k samples, which is a lot for my case; I divide the samples into 100 bins of equally spaced energy values, and count the number of samples in each bin. Then i plot these counts as a histogram against the energy values. X-axis is the energy values, and Y-axis is the count in each bin of energy values. This histogram should represent the "learned" target distribution. Since my kernel is e^(-1 f(p) ), I should get a monotonically decreasing curve as the energy values increase. That is, the histogram should look like the curve of e^(-x). However, I get a bell-curve like a gaussian, with a sharp peak in the middle. But* , the sampled values are as I expect them to be.

On the other hand, if I return negative log likelihood from my kernel and adjust the gradients accordingly, and repeat the above process, the histogram I plot looks like the curve of e^( - x) , just as expected. However, in this case the sampled samples don't look right.

Could you please clarify what's happening here, and what should the function return for HMC?

Thanks!