How to use custom Hamiltonian?

[rzu512]

Hi.

The parameters of my model have no constraint.
I implemented the following two functions and want to use HMC algorithm of zhusuan.
How to do that? If I can't use the function directly, how to implement custom StochasticTensor?

def cal_log_likelihood(p1, p2):
    """Calculate log-likelihood of the model

    Parameters
    ----------
    p1 : tf.Tensor, shape=(n_chains, n_states)
    p2 : tf.Tensor, shape=(n_chains, n_states, n_states)

    Returns
    -------    
    log_likelihood : tf.Tensor, shape=(n_chains,)
    """
    return log_likelihood

def cal_prior(p1, p2):
    """Calculate logarithm of the prior probability

    Parameters
    ----------
    p1 : tf.Tensor, shape=(n_chains, n_states)
    p2 : tf.Tensor, shape=(n_chains, n_states, n_states)

    Returns
    -------    
    log_prior_probability : tf.Tensor, shape=(n_chains,)
    """
    return log_prior_probability

Thank you.

[csy530216]

Hi, for an introductory example of using HMC, you can refer to the toy example here, where HMC is used to sample a Gaussian-distributed variable: https://github.com/thu-ml/zhusuan/blob/master/examples/toy_examples/gaussian.py. You can also look at other examples listed in our README with HMC or SGMCMC as the inference method, which are more practical and complicated than the toy gaussian example.

For your question, it is best to implement your probabilistic model using BayesianNet in ZhuSuan, then your model will automatically have a log_joint property corresponding to the posterior distribution, then you can pass the MetaBayesianNet instance to HMC.

A more direct approach for your question is to implement a log_joint function by yourself. A minimal description is as follows:

# Suppose the latent variable you want to sample is p.

def log_prior(p):
    # Calculate the log prior density...
    return log_prior

def log_likelihood(p):
    # Calculate the log-likelihood...
    return log_likelihood

def log_joint(observed):
    return log_prior(observed['p']) + log_likelihood(observed['p'])

hmc = zs.HMC(step_size=1e-3)
# A tf.Variable should be created to store the current samples.
p = tf.Variable(...)
sample_op, _ = hmc.sample(log_joint, observed={}, latent={'p': p})

with tf.Session() as sess:
    sess.run(tf.global_variables_initializer())
    sess.run(sample_op)
    # Now the value of the tf.Variable p will be updated; p stores the current samples.

For this kind of usage (passing log_joint callable function instead of a MetaBayesianNet instance to hmc.sample), please refer to some older versions of our examples, such as https://github.com/thu-ml/zhusuan/blob/fe8b7cba1e0ec31b37eeb6652d0e87d28fad18fd/examples/topic_models/lntm_mcem.py.

Feel free to ask if you have more questions :)

[thjashin]

just the log joint probability

[rzu512]

What will happen if the log-joint probability is -inf or nan? Would the HMC automatically reject the proposal in these two cases?

[rzu512]

The HMC works well! Thanks a lot. The mass and step-size adaptations really help to fit my model.

[csy530216]

Very glad to hear that!
Regarding the issue of nan, HMC should deal with such cases correctly (see #78).