## Riemann, Langevin & Hamilton [reply]

Posted in R, Statistics, University life with tags , , , , , on September 27, 2010 by xi'an

Here is a (prompt!) reply from Mark Girolami corresponding to the earlier post:

In preparation for the Read Paper session next month at the RSS, our research group at CREST has collectively read the Girolami and Calderhead paper on Riemann manifold Langevin and Hamiltonian Monte Carlo methods and I hope we will again produce a joint arXiv preprint out of our comments. (The above picture is reproduced from Radford Neal’s talk at JSM 1999 in Baltimore, talk that I remember attending!) Although this only represents my preliminary/early impression on the paper, I have trouble with the Physics connection. Because it involves continuous time events that are not transcribed directly into the simulation process.

It might be easier to look at the two attached figures rather than the isocontours in the complete phase space (q, p) of the 1-d Gaussian that is in the pic you include in your post. Both figures relate to sampling from a zero mean bivariate Gaussian with covariance matrix $\Sigma$ — marginal variance 1 and cross-covariance $\rho = 0.6$ (example taken from Neal 2010). In the first figure there are three panes showing the 20 integration steps obtained by standard HMC where the metric tensor for the space of bivariate Gaussians is NOT employed and so the leftmost plot shows the proposal effectively moving from -1.5, -1 to 1.5, 1.5 basically a large traversal over the 2d sampling space. The middle plot shows the corresponding momentum variable steps and the right plot shows the total Hamiltonian value (energy or negative joint likelihood) where the difference at the start and end of the integration is very small indicating a high probability of acceptance. That is all fine – large proposal step accepted with high probability.

Now lets look at the second figure. Lay aside the Physics interpretation and let’s try to adopt a geometric one in this case to see if this helps the understanding. So our task is to simulate from this bivariate Gaussian we know that the metric tensor $\Sigma^{-1}$ defines a flat manifold of bivariate densities. Now if we wanted to move from one point to another we would wish to take the most direct route – the shortest path in the space — the geodesic in other words. The geodesics on a manifold are defined by the second order differential equation in terms of the coordinates (our sample space) and the Christofell symbols of the manifold — so to make a proposed move along the geodesic we need to solve the differential equations describing the geodesic. Now by rewriting the geodesic equation in a different form we end up with Hamiltons equations -— in other words the solution flow of the geodesic equations coincides with the Hamiltonian flows — so solving Hamiltons equations is just another way to solve the geodesic equation. So the variable p just emerges from the rewriting of the geodesic equation and nothing more. In this case as the manifold corresponding to bivariate Gaussians is flat and the metric tensor is the inverse of the covariance matrix — the geodesics will be elliptical paths defined by $\Sigma$ — in other words the isocontours of the bivariate Gaussian. So the first panel in this figure shows 15 integration steps and as can be observed the path follows the elliptical isocontour of the target density — as the p variable is the dual in the geodesic equations then it also maps out an isocontour defined by the metric tensor. Continue reading