Hamiltonian Monte Carlo on TensorFlow and Edward

Overview

On this article, I tried Hamiltonian Monte Carlo algorithm to the simple data by TensorFlow and Edward.
Edward lets us use variational inference, Gibbs sampling and Monte Carlo method. And by relatively small changes, we can switch the methods. So I'll try simple HML model here.
About Hamiltonian Monte Carlo itself, I'll write another article for it.

Data

I'll use the same data as the one shown on the following article.

• Simple Bayesian modeling by Stan

About the way to the data, please check that article.
Roughly, the data has just two columns, speed and dist. I’ll use “speed” as explaining variable and “dist” as explained variable.

Modeling

As I wrote, by Edward, we can switch methods relatively easily. Before, to this data, I made a model with Edward using variational inference. This time, I’ll basically use the same code. Only some points I'll change for switching the method to Hamiltonian Monte Carlo.

• Simple regression by Edward: variational inference

The following code is from the article. I changed some trivial points but it is not important. Anyway, this uses variational inference.

import pandas as pd
import tensorflow as tf
import edward as ed 
from edward.models import Normal
import numpy as np

cars = pd.read_csv('cars.csv')

#the number of data points
N = cars.shape[0]
D = 1

X = tf.placeholder(tf.float32, [N, D])

b = Normal(loc=tf.zeros(D), scale=tf.ones(D))
a = Normal(loc=tf.zeros(1), scale=tf.ones(1))
y = Normal(loc=ed.dot(X, b) + a, scale=tf.ones(N))

q_b = Normal(loc=tf.Variable(tf.zeros(D)),
              scale=tf.nn.softplus(tf.Variable(tf.zeros(D))))
q_a = Normal(loc=tf.Variable(tf.zeros(1)),
              scale=tf.nn.softplus(tf.Variable(tf.zeros(1))))

inference = ed.KLqp({b: q_b, a: q_a}, data={X: np.array(cars['speed']).reshape([N, 1]), y: np.array(cars['dist'])})
inference.run(n_iter=1000)

To switch the code to the one for Hamiltonian Monte Carlo, we need to change two points, meaning import part, inference part.

# import Empirical
from edward.models import Normal, Empirical

# inference part
points = 10000
q_b = Empirical(params=tf.Variable(tf.random_normal([points, D])))
q_a = Empirical(params=tf.Variable(tf.random_normal([points, 1])))

inference = ed.HMC({b: q_b, a: q_a}, data={X: np.array(cars['speed']).reshape([N, 1]), y: np.array(cars['dist'])})
inference.run(step_size=1e-3)

For HMC, Hamiltonian Monte Carlo, we need to import Empirical. On the inference part, we also need to set the Empirical responding to the parameters. The variable, points, means the number of samples.
To check the sampled points, we can plot it.

import matplotlib.pyplot as plt
plt.subplot(2,1,1)
plt.title("b")
plt.hist(q_b.sample(1000).eval())
plt.subplot(2,1,2)
plt.title("a")
plt.hist(q_a.sample(1000).eval())
plt.show()