Poisson-gamma regression

Hi all,

I am new to Bayesian inference via Edward. I am implementing a Poisson-gamma regression in Edward. I compiled my Edward’s model. However, I have been unsuccessful with running Variational Bayes inference in Edward. Here is a snippet of code with the error that I’m getting:

N = 120
M = 5

x_train = np.reshape(z[:N], (N, 1))
y_train = np.reshape(w[:N], (N, 1))
t_train = np.reshape(range(N), (N, 1))
s_train = np.array([[12.0], [36.0], [60.0], [84.0], [108.0]])

X = tf.placeholder(tf.float32, [N, 1])
t = tf.placeholder(tf.float32, [N, 1])
s = tf.placeholder(tf.float32, [M, 1])

# Define priors
theta = Normal(mu=tf.zeros(1), sigma=tf.ones(1))
theta_star = Normal(mu=tf.zeros(1), sigma=tf.ones(1))
nu0 = Normal(mu=tf.ones(1)*1.5, sigma=tf.ones(1))
nu1 = Normal(mu=tf.ones(1)*2.0, sigma=tf.ones(1))

alpha_ = tf.exp(nu0)
beta_ = tf.exp(nu1)
gamma = Gamma(alpha=tf.ones([M, 1])*alpha_, beta=tf.ones([M, 1])*beta_)


def Kernel(t, s):
    mat = []
    for i in range(N):
        mat += [[]]
        for j in range(M):
            mat[i] += [1.2*tf.exp(-0.25*tf.pow(t[i] - s[j], 2))]
            
        mat[i] = tf.stack(mat[i])
    return tf.reshape(tf.stack(mat), [N, M])

# Poisson intensity
mu = theta*X + theta_star*tf.matmul(Kernel(t, s), gamma)

# Define likelihood
y = Poisson(lam=mu, value=tf.zeros_like(mu))


qtheta = Normal(mu=tf.Variable(tf.random_normal([1])),
            sigma=tf.nn.softplus(tf.Variable(tf.random_normal([1]))))
qtheta_star = Normal(mu=tf.Variable(tf.random_normal([1])),
            sigma=tf.nn.softplus(tf.Variable(tf.random_normal([1]))))

# qnu0 = Normal(mu=tf.Variable(tf.random_normal([1])), sigma=tf.nn.softplus(tf.Variable(tf.random_normal([1]))))
# qnu1 = Normal(mu=tf.Variable(tf.random_normal([1])), sigma=tf.nn.softplus(tf.Variable(tf.random_normal([1]))))

# qgamma = Gamma(alpha=tf.nn.softplus(tf.Variable(tf.random_normal([M, 1]))), beta=tf.nn.softplus(tf.Variable(tf.random_normal([M, 1]))))


inference = ed.KLqp({theta: qtheta, theta_star: qtheta_star}, data={X: x_train, t: t_train, s: s_train, y: y_train})
---------------------------------------------------------------------------
InvalidArgumentError                      Traceback (most recent call last)
<ipython-input-48-5b104cd13dec> in <module>()
----> 1 inference = ed.KLqp({theta: qtheta, theta_star: qtheta_star}, data={X: x_train, t: t_train, s: s_train, y: y_train})

/usr/local/lib/python2.7/site-packages/edward/inferences/klqp.pyc in __init__(self, *args, **kwargs)
     52   """
     53   def __init__(self, *args, **kwargs):
---> 54     super(KLqp, self).__init__(*args, **kwargs)
     55 
     56   def initialize(self, n_samples=1, kl_scaling=None, *args, **kwargs):

/usr/local/lib/python2.7/site-packages/edward/inferences/variational_inference.pyc in __init__(self, *args, **kwargs)
     21   """
     22   def __init__(self, *args, **kwargs):
---> 23     super(VariationalInference, self).__init__(*args, **kwargs)
     24 
     25   def initialize(self, optimizer=None, var_list=None, use_prettytensor=False,

/usr/local/lib/python2.7/site-packages/edward/inferences/inference.pyc in __init__(self, latent_vars, data)
    121           var = tf.Variable(ph, trainable=False, collections=[])
    122           self.data[key] = var
--> 123           sess.run(var.initializer, {ph: value})
    124         elif isinstance(value, np.number):
    125           if np.issubdtype(value.dtype, np.float):

/usr/local/lib/python2.7/site-packages/tensorflow/python/client/session.pyc in run(self, fetches, feed_dict, options, run_metadata)
    765     try:
    766       result = self._run(None, fetches, feed_dict, options_ptr,
--> 767                          run_metadata_ptr)
    768       if run_metadata:
    769         proto_data = tf_session.TF_GetBuffer(run_metadata_ptr)

/usr/local/lib/python2.7/site-packages/tensorflow/python/client/session.pyc in _run(self, handle, fetches, feed_dict, options, run_metadata)
    963     if final_fetches or final_targets:
    964       results = self._do_run(handle, final_targets, final_fetches,
--> 965                              feed_dict_string, options, run_metadata)
    966     else:
    967       results = []

/usr/local/lib/python2.7/site-packages/tensorflow/python/client/session.pyc in _do_run(self, handle, target_list, fetch_list, feed_dict, options, run_metadata)
   1013     if handle is None:
   1014       return self._do_call(_run_fn, self._session, feed_dict, fetch_list,
-> 1015                            target_list, options, run_metadata)
   1016     else:
   1017       return self._do_call(_prun_fn, self._session, handle, feed_dict,

/usr/local/lib/python2.7/site-packages/tensorflow/python/client/session.pyc in _do_call(self, fn, *args)
   1033         except KeyError:
   1034           pass
-> 1035       raise type(e)(node_def, op, message)
   1036 
   1037   def _extend_graph(self):

InvalidArgumentError: Node 'inference_15065146640/0/Variable_24/read': Unknown input node 'Variable_24'

I swapped your x_train, y_train, t_train lines with the following in order to run the script

x_train = np.zeros((N, 1))
y_train = np.zeros((N, 1))
t_train = np.zeros((N, 1))

The script terminated successfully.

Could you try it again? My guess is that the iPython kernel was interacting strangely with TensorFlow sessions, so you need only restart the kernel.

Hi Dustin,

I restarted the kernel, and I swapped the few lines of code you suggested. The script ran successfully.

1 Like

Hi Dustin,

I realize the problem I am having is with the input tensors shapes and types. x_train, y_train, and t_train are integer matrices. The input data looks like this:

u = np.reshape(z, (120, 1))

u = u.astype(float)

w = np.reshape(w, (120, 1))

z = np.reshape(range(121)[1:], (120, 1))

z =z.astype(float)

u[:5],w[:5], z[:5]

(array([[ 1.],
        [ 2.],
        [ 3.],
        [ 4.],
        [ 5.]]), array([[4],
        [5],
        [9],
        [6],
        [8]]), array([[ 1.],
        [ 2.],
        [ 3.],
        [ 4.],
        [ 5.]]))

Here is the code:

N = 120
M = 5

x_train = u
y_train = w
t_train = z
s_train = np.array([[12.00], [36.00], [60.00], [84.00], [108.00]])

X = tf.placeholder(tf.float32, [N, 1])
t = tf.placeholder(tf.float32, [N, 1])
s = tf.placeholder(tf.float32, [M, 1])

# Define priors
theta = Normal(mu=tf.zeros(1), sigma=tf.ones(1))
theta_star = Normal(mu=tf.zeros(1), sigma=tf.ones(1))
nu0 = Normal(mu=tf.ones(1)*1.5, sigma=tf.ones(1))
nu1 = Normal(mu=tf.ones(1)*2.0, sigma=tf.ones(1))

alpha_ = tf.exp(nu0)
beta_ = tf.exp(nu1)
gamma = Gamma(alpha=tf.ones([M, 1])*alpha_, beta=tf.ones([M, 1])*beta_)


def Kernel(t, s):
    mat = []
    for i in range(N):
        mat += [[]]
        for j in range(M):
            mat[i] += [1.2*tf.exp(-0.25*tf.pow(t[i] - s[j], 2))]
            
        mat[i] = tf.stack(mat[i])
    return tf.reshape(tf.stack(mat), [N, M])

# Poisson intensity
mu = theta*X + theta_star*tf.matmul(Kernel(t, s), gamma)

# Define likelihood
y = Poisson(lam=mu, value=tf.zeros_like(mu))


qtheta = Normal(mu=tf.Variable(tf.random_normal([1])),
            sigma=tf.nn.softplus(tf.Variable(tf.random_normal([1]))))
qtheta_star = Normal(mu=tf.Variable(tf.random_normal([1])),
            sigma=tf.nn.softplus(tf.Variable(tf.random_normal([1]))))

qnu0 = Normal(mu=tf.Variable(tf.random_normal([1])), sigma=tf.nn.softplus(tf.Variable(tf.random_normal([1]))))
qnu1 = Normal(mu=tf.Variable(tf.random_normal([1])), sigma=tf.nn.softplus(tf.Variable(tf.random_normal([1]))))
qgamma = Gamma(alpha=tf.nn.softplus(tf.Variable(tf.random_normal([M, 1]))), beta=tf.nn.softplus(tf.Variable(tf.random_normal([M, 1]))))

inference = ed.KLqp({theta: qtheta, theta_star: qtheta_star, nu0: qnu0, nu1: qnu1, gamma: qgamma}, data={X: x_train, t: t_train, s: s_train, y: y_train})

inference.run(n_samples=5, n_iter=250)

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-34-ec1ea2c18385> in <module>()
----> 1 inference.run(n_samples=5, n_iter=250)

/usr/local/lib/python2.7/site-packages/edward/inferences/inference.pyc in run(self, variables, use_coordinator, *args, **kwargs)
    188       Passed into ``initialize``.
    189     """
--> 190     self.initialize(*args, **kwargs)
    191 
    192     if variables is None:

/usr/local/lib/python2.7/site-packages/edward/inferences/klqp.pyc in initialize(self, n_samples, kl_scaling, *args, **kwargs)
     76     self.n_samples = n_samples
     77     self.kl_scaling = kl_scaling
---> 78     return super(KLqp, self).initialize(*args, **kwargs)
     79 
     80   def build_loss_and_gradients(self, var_list):

/usr/local/lib/python2.7/site-packages/edward/inferences/variational_inference.pyc in initialize(self, optimizer, var_list, use_prettytensor, *args, **kwargs)
     62       var_list = list(var_list)
     63 
---> 64     self.loss, grads_and_vars = self.build_loss_and_gradients(var_list)
     65 
     66     if optimizer is None:

/usr/local/lib/python2.7/site-packages/edward/inferences/klqp.pyc in build_loss_and_gradients(self, var_list)
    124       #    return build_score_entropy_loss_and_gradients(self, var_list)
    125       else:
--> 126         return build_score_loss_and_gradients(self, var_list)
    127 
    128 

/usr/local/lib/python2.7/site-packages/edward/inferences/klqp.pyc in build_score_loss_and_gradients(inference, var_list)
    531     for x in six.iterkeys(inference.data):
    532       if isinstance(x, RandomVariable):
--> 533         x_copy = copy(x, dict_swap, scope=scope)
    534         p_log_prob[s] += tf.reduce_sum(
    535             inference.scale.get(x, 1.0) * x_copy.log_prob(dict_swap[x]))

/usr/local/lib/python2.7/site-packages/edward/util/random_variables.pyc in copy(org_instance, dict_swap, scope, replace_itself, copy_q)
    174     kwargs['name'] = new_name
    175     # Create new random variable with copied arguments.
--> 176     new_rv = rv.__class__(*args, **kwargs)
    177     return new_rv
    178   elif isinstance(org_instance, tf.Tensor):

/usr/local/lib/python2.7/site-packages/edward/models/random_variable.pyc in __init__(self, *args, **kwargs)
     70 
     71     if value is not None:
---> 72       t_value = tf.convert_to_tensor(value, self.dtype)
     73       expected_shape = (self.get_batch_shape().as_list() +
     74                         self.get_event_shape().as_list())

/usr/local/lib/python2.7/site-packages/tensorflow/python/framework/ops.pyc in convert_to_tensor(value, dtype, name, preferred_dtype)
    635       name=name,
    636       preferred_dtype=preferred_dtype,
--> 637       as_ref=False)
    638 
    639 

/usr/local/lib/python2.7/site-packages/tensorflow/python/framework/ops.pyc in internal_convert_to_tensor(value, dtype, name, as_ref, preferred_dtype)
    700 
    701         if ret is None:
--> 702           ret = conversion_func(value, dtype=dtype, name=name, as_ref=as_ref)
    703 
    704         if ret is NotImplemented:

/usr/local/lib/python2.7/site-packages/tensorflow/python/ops/variables.pyc in _TensorConversionFunction(v, dtype, name, as_ref)
    647       raise ValueError(
    648           "Incompatible type conversion requested to type '%s' for variable "
--> 649           "of type '%s'" % (dtype.name, v.dtype.name))
    650     if as_ref:
    651       return v._ref()  # pylint: disable=protected-access

ValueError: Incompatible type conversion requested to type 'float32' for variable of type 'int32_ref'

Good catch. This is because the Poisson random variable uses a float data type and not int:

>>> from edward.models import Poisson
>>> x = Poisson(lam=5.0, value=0)
>>> print(x)
RandomVariable("Poisson_1/", shape=(), dtype=float32)

This means any data binded to Poisson random variables must be float and not int. Similarly, you defined your placeholders for t and s to be of type float, although you’re binding it to data of type int.

I added a Github issue to raise a more informative error message if there’s a type mismatch (https://github.com/blei-lab/edward/issues/560).

As an aside, I’m not actually sure why Poisson uses float and not int. It is a result of this [line]
(https://github.com/tensorflow/tensorflow/blob/master/tensorflow/contrib/distributions/python/ops/poisson.py#L86) in TensorFlow; maybe because of unbounded support?

Hi Dustin,

I’m a postdoc in Dr. Barbara Engelhardt’s lab at Princeton. I’m working on a spatial regression for a marked point Poisson process with the Poisson mean that is a sum of fixed effects part and a doubly stochastic process represented by a convolution of a Gaussian kernel and gamma process. The placeholders t and s are inputs to the kernel and gamma is gamma process in the code. I am getting nan’s for the loss when I ran the variational inference. I’m going to try to run a simple Poisson regression to understand how Poisson distribution works in Edward. You have helped me a lot.

Great! Another dimension worth experimenting is to not simplify the model but to simplify the inference. That is, I recommend just trying MAP. Suchi Saria’s lab has produce some excellent healthcare applications using Bayesian nonparametric models and MAP.