1
$\begingroup$

I'm using the affine-invariant sampler from emcee to draw samples from a $p$ dimensional posterior, using $M$ parallel chains ($M>10$).

Since my model is p-dimensional with $p>1$, I'm also using the multiESS method to estimate the effective sample size.

My question is: since I have $M$ chains, how should I apply multiESS? Should I apply it over each chain separately? Should I concatenate all the chains and apply it on the single flat chain? Is one method more appropriate than the other?

One thing to keep in mind (I'm not sure if this will affect the decision) is that in the affine-invariant sampler the chains are not independent.

Improve this question
$\endgroup$
2
  • $\begingroup$ It depends on what your final Monte Carlo estimator is. As in, how are you combining the $M$ non-indepedent samples to obtain estimates of posterior quantities? $\endgroup$ – Greenparker Jul 28 '18 at 11:11
  • $\begingroup$ I generate the single flat chain (ie: I concatenate all the chains) for each parameter, and from that I obtain the marginalized posterior for each one. $\endgroup$ – Gabriel Jul 29 '18 at 15:39
1
$\begingroup$

I am using notation from the Goodman and Weare (2010) paper. This is how they describe their estimation problem in Section 3.

There are $L$ walkers in the ensemble, and each at each time point $t$, what is observed is $X(t) = (X_1(t), \dots, X_L(t) )$, for $t = 1, \dots, T_e$. The final estimator for $E[f(x)]$ is, $$\hat{A}_e = \dfrac{1}{T_e} \sum_{t=1}^{T_e} \left( \dfrac{1}{L} \sum_{k=1}^{L} f(X_k(t)) \right) \,.$$

If the walkers were all independent, then you could calculate multiESS and ess for each walker and then add them all together (as long as there was proper burn-in. Burn-in in this case more important than burn in for a single long run).

However, the walkers are not independent, so they define $$F(X(t)) = \dfrac{1}{L} \sum_{k=1}^{L} f(X_k(t)) \, \quad \text{ and thus } \quad\hat{A}_e = \dfrac{1}{T_e}\sum_{t=1}^{T_e} F(X_{t})\, $$ That is, $F(X(t))$ is the average of $f(x)$ obtained from $L$ walkers at time $t$.

Now here we make the assumption that $F(X(t))$ has a Markovian structure (or atleast a $\phi$-mixing structure). I can't confirm this assumption because I don't quite understand what the sampler is doing in the paper, but they make a similar assumption in Section 3 when they define their $\tau_e$.

Assuming this, $F(X_1), F(X_2), \dots, F(X_{T_e})$ each has variance $\Omega/L$ and exhibit structured correlation so that for the estimator $\hat{A}_e$, so we have the variance $$\text{Var}(\hat{A}_e) = \dfrac{\Sigma}{T_e} = \dfrac{L\Sigma}{LT_e}$$

Compare this with the estimator for if all observations were IID (let $\Lambda = \text{Var}_{\pi} f(x)$), $$\text{Var}(\hat{A}_s) = \dfrac{\Lambda}{LT_e} $$

So the multivariate effective sample size is then, $$mESS = LT_e \left(\dfrac{\det(\Lambda)}{\det(L\Sigma)} \right)^{1/p}\,. $$

Note that to implement this using multiESS in mcmcse, you must calculate multi.mcse on $F(X(t))$ and not on $f(X(t))$, but to calculate $\Lambda$ you should concatenate all the data, and then calculate the sample covariance. To avoid confusions, I would suggest, don't use the function multiESS, and calculate using the determinant function in R. The code should be similar to the following.

library(mcmcse)
L <- 3
Te <- 1e3
walk1 <- matrix(rnorm(Te*4), nrow = 1000, ncol = 4)
walk2 <- matrix(rnorm(Te*4), nrow = 1000, ncol = 4)
walk3 <- matrix(rnorm(Te*4), nrow = 1000, ncol = 4)

Fx <- (walk1 + walk2 + walk3)/L
Sigma <- mcse.multi(Fx, size = "cuberoot")$cov

concat_full <- rbind(walk1, walk2, walk3)
Lambda <- var(concat_full)

mess <- (L*Te) * (det(Lambda)/det(L*Sigma) )^(1/4)
Improve this answer
$\endgroup$
4
  • $\begingroup$ Thank you for the detailed answer Greenparker! A few questions/cmmts if you don't mind, so I'm sure I understood: 1. given that $n=LT_e$, one difference with the mcse.multiESS formula is that I should use $L\Sigma$ instead of $\Sigma$ in the denominator. Is this correct? 2. The second difference is that I should obtain $\Sigma$ (the output of mcse.multi?) using $F(X(t))$, that is: the average of $f(X(t))$) over my $L$ chains? $\endgroup$ – Gabriel Jul 29 '18 at 15:42
  • $\begingroup$ 3. $\Lambda$ is obtained the same way I would for a single chain (ie: the standard way used by mcse.multiESS), and in this case the single chain would be my $L$ concatenated chains, is this correct? $\endgroup$ – Gabriel Jul 29 '18 at 15:43
  • $\begingroup$ 4. If the above is correct, then I must obtain $\Lambda$ from my concatenated chains of length $LT_e$, but $\Sigma$ from my averaged $L$ chains of length $T_e$. $\endgroup$ – Gabriel Jul 29 '18 at 15:50
  • 1
    $\begingroup$ You got most of it right; just that $\Lambda$ cannot be obtained from multiESS, but from using the "var" function on the concatenated chain. I have put some pseudo code for when you have 3 walkers of length 1000, for a 4 dimensional problem (assuming f is the identity function). $\endgroup$ – Greenparker Jul 30 '18 at 16:27

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.