search

minaskar / zeus

⚡️ zeus: Lightning Fast MCMC ⚡️
https://zeus-mcmc.readthedocs.io/
GNU General Public License v3.0
227 stars 33 forks source link

"mk" autocorrelation calculation #22

Open grahamgower opened 3 years ago

grahamgower commented 3 years ago

Hi @minaskar,

I looked at the autocorrelation stuff you've just added, and it looks really useful. Thanks!

https://github.com/minaskar/zeus/blob/658ad948b6f2e6edff88af8632b064b651e234fc/zeus/autocorr.py#L54-L65

I had a question/comment about the "mk" calculation. It looks like you're wanting to calculate the autocorrelation over each of the chains, which makes more sense to me than averaging over the chains (dfm) or averaging before calculating autocorrelation (Goodman-Weare). But with the mk method, by simply reshaping the 2d y array into a 1d array, the autocorrelation will be artificially lower than it should be, on account of also look at the correlation across the boundary points between chains. I.e. if the chains each have length n, then the autocorrelation of y2 = y.reshape((-1)) will be artifically lower at each position in y2 that is an integral multiple of n.

Wouldn't it be more appropriate to instead take ffts of each chain separately and average over the chains in frequency domain?

minaskar commented 3 years ago

Hi @grahamgower,

Thanks for your comment.

I think I agree with what you say about the "MK" recipe although I believe the effect is very minor. I haven't seen a case in which it actually affects the outcome. Please let me know if you can think of anything where this could be a problem. A nice property of this estimator, apart from its lower variance (see the figure), is that it can help diagnose violations of ergodicity. E.g. If the walkers are unable to jump from mode to mode in a multimodal case the the MK estimate of the autocorrelation time would be equal to the length of the chain.

I've tried what you suggested and it seems to give similar results to the Goodman-Weare estimator as you can see in the figure bellow: compare

grahamgower commented 3 years ago

Thanks for your insights. I guess all these estimators are asymptotically unbiased when the chain length tends to inf.

I've tried what you suggested and it seems to give similar results to the Goodman-Weare estimator as you can see in the figure bellow:

Ah! This is obvious in hindsight because the fourier transform preserves linearity. Sorry, I should admit that I haven't looked at a fourier transform in quite some years (and never used one outside a classroom).

grahamgower commented 3 years ago

Here's another simple alternative that seems reasonable: https://arxiv.org/abs/2009.01799. They just take the ACvF of mean-centered chains, where the mean is calculated globally over all chains.

minaskar commented 3 years ago

Thanks! That looks interesting, I'll have a look.

I'm not 100% sure that this applies to chains like the ones that zeus generates (from interacting walkers) but I'll have to read the paper first.