This week of work proved to be quite productive. I was consistent with my proposed time line. For this week work I proposed to write base class structure for Hamiltonian Monte Carlo (HMC), implement methods for leapfrog, modified Euler and algorithm for finding a reasonable starting value of epsilon. During the start I wrote the leapfrog and modified Euler as methods of HMC class, but my mentor told me to write a base class and using that base class write leapfrog and modified Euler as different classes.
The earlier structure looked something like this:
class HamiltonianMCda(object):
def __init__(self, discretize_time='leapfrog', *args):
# Some arguments and parameters
# chooses discretization algorithm depending upon
# the string passed to discretize_time
pass
def leapfrog(self):
pass
def modifiedEuler(self):
pass
But why did we settle upon base class implementation? With the earlier structure things were not flexible from user point of view. What if user want to plug-in his own implementations. After the changes I created a base class called DiscretizeTime , class inheriting this could then be passed as an argument for discretization of time. Advantage of having a base class is that it provides a basic structure to the things, and adds extensibility. Now things look something like this:
class DiscretizeTime(object):
def __init__(self, *args):
pass
def discretize_time(self):
# returns the initialized values
pass
class LeapFrog(DiscretizeTime):
def __init__(self, *args):
pass
def _discretize_time(self):
# computes the values and initializes the parameters
pass
class HamiltonianMCda(object):
def __init__(self, discretize_time=LeapFrog, *args):
# discretize_time is a subclass of DiscretizeTime
pass
Now using these base class user can pass his/her own implementations as an argument.
I also wrote other base classes for finding Log of probability Distribution and Gradient of the log. During the time of writing my GSoC proposal me and my couldn’t decide how these methods for finding log and its gradient should be implemented. In this week meeting with the mentor, I proposed that I’ll write method in each model classes we will be implementing, and use that, if user doesn’t provide a class inheriting the base class for gradients and we settled upon it.
Though this week work was great, but still there are things which remain unclear, from theoretical point. How to parameterize a continuous model still remains in doubt. Currently I have assumed that model parameterization will be of a matrix or array type structure. This assumption is good enough for the most of the common models I have came across, but things cannot be stated with certainty that it will generalize to all kind of continuous models. Me and my mentor are looking into things more deeply and hopefully we will find some solution soon. For the next week I’ll try to finish the Hamiltonian Monte Carlo sampler.