Semi-Separable Hamiltonian Monte Carlo for Inference in Bayesian Hierarchical Models
Sampling from hierarchical Bayesian models is often difficult for MCMC methods, because of the strong correlations between the model parameters and the hyperparameters. Recent Riemannian manifold Hamiltonian Monte Carlo (RMHMC) methods have significant potential advantages in this setting, but are computationally expensive. We introduce a new RMHMC method, which we call semi-separable Hamiltonian Monte Carlo, which uses a specially designed mass matrix that allows the joint Hamiltonian over model parameters and hyperparameters to decompose into two simpler Hamiltonians. This structure is exploited by a new integrator which we call the alternating blockwise leapfrog algorithm. The resulting method can mix faster than simpler Gibbs sampling while being simpler and more efficient than previous instances of RMHMC.
The authors introduce a new Hamiltonian and claim that it is better than existing Hamiltonians for heirarchical models. The technical material is very well-written, and the Hamiltonian they present is (as far as I know) novel. It is also very plausible that the new Hamiltonian leads to substantial improvements. However, I found it to be most confusing (and least convincing) in the few places that it argues the new method is better than the old one.Major Issues:1. The authors say that they are measuring "effective sample size" but don't say how they are computing this; they reference the reader to [6]. I haven't read [6] carefully, but certainly it doesn't define effective sample size (indeed the article doesn't contain the word "effective" and uses "size" only to refer to actual counts). Normally I wouldn't be so picky, but I think this is critical to their claims. In particular, for 2 of the 3 tables, their proposed method is MUCH slower than the competitors as measured by clock time; they rely on normalized effective sample size to argue that their method is better despite this. I actually think that, almost regardless of the normalization method they use for computing ESS, this is going to be a pretty unconvincing argument (the same table could be used to show that almost any adaptive algorithm is better than almost any other, since adaptive algorithms care so much about initial burn-in periods). But it is especially unconvincing if the normalization is not defined.This leads me into a bigger reason that their numerics aren't as convincing as I might hope. It is clear from the last column of table 1 (and figure 1) that some of the competitor algorithms are still quite far from stationarity after the number of steps that they run. If each step had a similar cost, that would be a fair comparison: they would be showing that their algorithm converges more quickly, and that's that. The whole variance thing would be a bit of a red herring in this situation (most of the error for the competitors is coming from bias, not variance), but it would be a useful thing to discuss. However, from their own tables, this doesn't seem to be the case. Their algorithm seems to take much more time per step than the competitor algorithms, and they have chosen a NUMBER of STEPS such that their algorithm has converged but the competitors haven't. This really doesn't seem to be a fair comparison, since effective sample size, bias, etc. can all get better as the algorithm runs.To describe the unfairness a bit more graphically, I'll describe an artifical example that might get something like their table 1. The original algorithm the paper discusses is MCMC driven by kernel $K$ with mixing time $T$. The paper proposes a new algorithm: run $K^{7}$ (that is, take every 7'th point from a Markov chain driven by $K$). We then look at a sample obtained from running both algorithms for $\frac{T}{3}$ steps. Just like Table 3, the new algorithm would take about 7 times as long to run as measured by wall clock. Also just like table 3, we would expect the old algorithm to have enormous bias (it hasn't mixed yet) while the new algorithm would have very tiny bias (it has mixed well). Its clear in my artificial example that the comparison is unreasonable, and I think that for the same reason it is a little unreasonable in this paper.Minor Issues:1. Bottom of pp.2, the authors claim that some dynamics are not `valid' because they aren't reversible. They probably meant to say something about how they don't converge to the target.2. There seems to be a bit of a bait-and-switch at the bottom of page 3. Up until now, the new method has been motivated by pointing out that competitor algorithms scale poorly in high dimensions. When the scaling of this paper's algorithm is revealed, it looks pretty similar for the example under consideration, but this is brushed away as being unimportant. This isn't such a big deal for the applicability of the paper; it just makes the note feel a little dishonest. Alternatively, maybe this is really a big difference, in which case they should explain why.3. I wish section 4.3 were longer! This seems like the most important section, but it wasn't enough for me to figure out what was going on. This is the sort of section that would be *vastly* improved by *any* toy example where this heuristic could be made precise. In particular, it would be nice to see that there is some optimal value of total energy in the Gibbs version of HMC, and that their method of passing energy back and forth between components (rather than keeping it constant) actually does better than that optimal value. Even just stating the algorithm and values would be enough, which could be done in a few lines. If the claim really is obvious, I hope that they have such an example.4. In figure 1, I don't think the author's explain the coloring on pictures 2 and 4.