Using NUTS sampler for likelihood with matrix exponential

Question

I am trying to infer the generator of a continuous markov process observed at discrete intervals. If the generator of the markov process is $T$, then the stochastic matrix for the discrete time intervals is given by $ P = \exp(T \Delta t)$. To implement this using pymc, I wrote the custom distribution class

import pymc3
from pymc3.distributions import Discrete
from pymc3.distributions.dist_math import bound



class ContinuousMarkovChain(Discrete):

    def __init__(self, t10=None, t01=None, dt=None, *args, **kwargs):
        super(ContinuousMarkovChain, self).__init__(*args, **kwargs)
#         self.p = p
#         self.q = q
        self.p = tt.slicetype

        self.gt0 = (t01 >0) & (t10> 0) 

        T = tt.stacklists([[-t01, t01], [t10,-t10]])
        self.p = ts.expm(T*dt)


    def logp(self, x):

        return bound(tt.log(self.p[x[:-1],x[1:]]).sum(), self.gt0)

I can use find_MAP and the Slice sampler with this class, but it fails with NUTS. The error message is:

AttributeError: 'ExpmGrad' object has no attribute 'grad'

I thought that NUTS only needed information about the gradient, so why is it trying to take the Hessian of expm?

You mentioning "discrete" anything makes me think something will break :) What is `expm`? Is it differentiable? — chris, Jan 17 '17 at 01:32
`expm` is the matrix exponential function, and its derivative is implemented in theano. The distribution is discrete, but it is continuous wrt with respect to its parameters. — nbren12, Jan 17 '17 at 01:41

score 0 · Answer 1 · answered May 08 '17 at 15:50

I thought Pymc3 needs Hessian in parameter space to set the step-size and directionality for the parameters when using NUTS algorithm. Maybe you can define the grad of ExpmGrad yourself. A relative discussing is here https://github.com/pymc-devs/pymc3/issues/1226

Using NUTS sampler for likelihood with matrix exponential

1 Answers1