## Markov chain Monte Carlo using Metropolis-Hastings algorithm in Java

We will focus on a case in the Markov chain Monte Carlo (MCMC) domain. MCMC contains algorithms for sampling from a probability distribution. The sampling happens by a Markov chain, i.e. it is a random process that changes from one state to another. In this article, we will use the Metropolis-Hastings algorithm to generate a random walk with a proposal density. We will then either accept or reject the proposed moves according to a certain condition.

We will perform some experiments to study the behavior of the Metropolis-Hastings on a test suite of target distributions. The test suite consists in this case of a Gaussian and a uniform distribution in one dimension. The distributions are characterized by ther mean and variance and the target distribution is the mixture

$p\left ( x \right ) = w_{1}\frac{1}{\sigma _{1}\sqrt{2\Pi }}e^{-\frac{1}{2}\frac{\left (x-\mu _{1}\right )^{2}}{\sigma _{2}^{1}}} + w_{2}\frac{1}{\sigma _{2}\sqrt{2\Pi }}e^{-\frac{1}{2}\frac{\left (x-\mu _{2}\right )^{2}}{\sigma _{2}^{2}}}$

The weights $0 \leqslant w_{1}, w_{2} \leq 1$ and $w_{1} + w_{2} = 1$ measure the relative contribution of both Gaussian distributions. Without loss of generality, we can always assume that the first mean $\mu _{1}= 0$

As for a target distribution, we will us $p(x) \propto 0.3e^{-2x^{2}}+ 0.7e^{-0.2(x-10)^{2}}$

For additional notes, I refer to the paper An Introduction to MCMC for Machine Learning by C. Andrieu, N. de Freitas, A. Doucet and M. Jordan.

We will use symmetric proposal distributions in the Metropolis-Hastings algorithm. In that case the acceptance probability reduces to
$A(\boldsymbol{\mathrm{x}},\mathrm{\mathbf{x^{*}}}) = \mathrm{min} \begin{Bmatrix} 1, & \frac{p(\mathrm{\mathbf{x^{*}}})}{p(\mathrm{\mathbf{x}})} \end{Bmatrix}$

where $\mathbf{x} \enspace$ and $\mathbf{x}^{*}$ are respectively the current state and the proposed new state accepted with probability $A(\mathbf{x},\mathbf{x}^{*})$

The proposal distributions are the Gaussian
$q(x^{*},x) = \frac{1}{\sigma \sqrt{2\pi }}e^{-\frac{(x^{*}-x)^{2}}{2\sigma ^{2}}}$
and the uniform distribution over the interval $[x-r,x+r]$
$q(x^{*},x) = \left\{\begin{matrix} \frac{1}{2r} & \mathrm{for }\enspace x^{*}\in [x-r,x+r]\\ 0 & \mathrm{otherwise} \end{matrix}\right.$

Now, we can begin with the experiments. We will collect 5,000 samples after the chain has converged to the invariant distribution, i.e. after the burn-in period. We will then construct a histogram based on the collected samples and compared with the target distribution. For these experiments, we will use Java.

In order to be able to write code for sampling and constructing a histogram, we will use two Java libraries: The Commons Math library and the JFreeChart library.

First, we will create two arrays of size 5000. One array is for the generated numbers according to the Gaussian distribution and the other array will hold the generated numbers according to the uniform distribution.

double[] arrayG;
double[] arrayU;

arrayG = new double[5000];
arrayU = new double[5000];


After creating these arrays for the numbers of both distributions, we will now use the distribution classes in order to be able to do the sampling. In this example, the sigma value for the normal distribution will be 0.5, whereas the mean is 0.2. In the case of the uniform distribution, the sigma value is 0.0 and the mean is 100.

NormalDistribution nor = new NormalDistribution(0.5, 0.2);

UniformRealDistribution uni = new UniformRealDistribution(0.0, 100);


For the sampling, we will use a for-lus for each distribution. We then create an array for each probability.

for(int i = 0; i < 5000; i++)
{

arrayG[i] = nor.sample(); }

for(int i = 0; i < 5000; i++)
{
arrayU[i] = uni.sample(); }

double[] probG;
prob = new double[5000];

double[] probU;
prob2 = new double[5000];


After the arrays for the probabilities are initialized, we will now calculate the probabilities by using the density function. For the accepted values using the Metropolis-Hastings algorithm, we will use two arraylists. This means that the accepted values will be added to these lists according to their distributions.

for(int j = 0; j < 5000; j++)
{
probG[j] = nor.density(arrayG[j]);
probU[j] = uni.density(arrayU[j]); }

List arraylistG = new ArrayList();
List arraylistU = new ArrayList();

for(int i = 0; i < 4999; i++)
{
double check = probG[i+1]/prob[i];
if(check >= 1)
{
}

for(int i = 0; i < 4999 ; i++)
{
double check = probU[i+1]/prob2[i];
if(check >= 1)
{
}


So far, we sampled 5000 times from two distributions and added the accepted values from the Metropolis-Hastings algorithm to an arraylist. The only thing we need to do now is plotting our experiments in order to have a visual understanding of the behaviour of the samples within a certain range. As mentioned earlier, we used the JFreeChart to construct a histogram with our numbers. Since we worked with two distributions, we will therefore create two histograms. The code of these histograms is almost similar.

I wrote also some code to print the result in the console. We will use this code to construct our histograms, i.e. the dataset class and its function addSeries will use the print code to draw the results.

double[] printG;
printG = new double[arraylistG.size()];
for(int i = 0; i < arraylistG.size(); i++)
{
printG[i] = (double)arraylistG.get(i); }

double[] printU;
printU = new double[arraylistU.size()];
for(int i = 0; i < arraylistU.size(); i++)
{
printU[i] = (double)arraylistU.get(i); }


Now, we are close to the end. We will use the code listed above to construct our histograms, both for the Gaussian distribution and the uniform distribution.

int number = 50;
HistogramDataset dataset = new HistogramDataset();
dataset.setType(HistogramType.RELATIVE_FREQUENCY);
String plotTitle = "M-H algorithm for Gaussian distribution";
String xaxis = "number";
String yaxis = "value";
PlotOrientation orientation = PlotOrientation.VERTICAL;
boolean show = false;
boolean toolTips = false;
boolean urls = false;
JFreeChart chart = ChartFactory.createHistogram( plotTitle, xaxis, yaxis,
dataset, orientation, show, toolTips, urls);
int width = 500;
int height = 300;
try {
ChartUtilities.saveChartAsPNG(new File("histogram.PNG"), chart, width, height);
} catch (IOException e) {}

HistogramDataset dataset2 = new HistogramDataset();
dataset2.setType(HistogramType.RELATIVE_FREQUENCY);
String plotTitle2 = "M-H algorithm for uniform distribution";
String xaxis2 = "number";
String yaxis2 = "value";
PlotOrientation orientation2 = PlotOrientation.VERTICAL;
boolean show2 = false;
boolean toolTips2 = false;
boolean urls2 = false;
JFreeChart chart2 = ChartFactory.createHistogram( plotTitle2, xaxis2, yaxis2,
dataset2, orientation2, show2, toolTips2, urls2);
int width2 = 500;
int height2 = 500;
try {
ChartUtilities.saveChartAsPNG(new File("histogram2.PNG"), chart2, width2, height2);
} catch (IOException e) {}


Running the code will result in the creation of two histograms, one for each distribution we used in our experiments.

Histogram of Gaussian distribution sampling

Histogram of uniform distribution sampling

## An introduction to Windows Runtime and Windows Runtime objects

Windows Runtime (WinRT) is the principal provider of system services for applications in Windows 8 that will be designed to incorporate the new “Metro” style. It is Microsoft’s new programming model that supports development in C++, C#, VB.NET as well as JavaScript. WinRT components are designed with an eye to interoperability between multiple languages and APIs, including native (C++), managed (C#) and scripting (JavaScript).

WinRT is not a replacement for any part of Windows we have come to know. Desktop applications written in earlier versions of Windows, will run properly in the desktop mode on Windows 8.

The following picture shows an overview of the traditional desktop applications and the new metro style applications. We also see the APIs that supports the programming languages in the new platform.

Through the above APIs, the developer can communicate with other APIs on the basis of numerous objects. There are more than 800 objects available for approaching APIs. Every object implements several interfaces. These interfaces are located in the registry.

These are the API’s of Windows Runtime:

Every WinRT object supports interfaces, and two of these interfaces are very important: IUknown and IInspectable. The first one is a sort of a parent interface where every other interface inherit from it, directly or indirectly. IIspectable is an interface to inspect an object and to know what the members are. If you want to create an object in WinRT, for example, to communicate with the webcam, you call a certain class. This class looks then in the Activation Store, which is located in the registry, to implement the necessary interfaces. Once the interfaces are implemented, this information is stored in the Windows metadata.

Below there is an example of the FileInformation class which has 3 interfaces. As a developer you will not see these interfaces, but they will be implemented as soon as you use the FileInformation class. What exactly happens here is that the FileInformation class will obtain the interfaces from the Activation Store. The object is then activated from the Registry, and stored in the Windows metadata.

So once we have created an object, we can communicate and interact with it. On the left of the image below, we have our object, and at the right of it there are  our applications written in C + +, C # and JavaScript. Now the question is: when will our applications communicate with the object to use its services?

If it’s a C + + application, the compiler goes at compile time to the metadata and takes all the information it needs out of it, and puts it in the executable at runtime so you no longer need to look at the metadata. So it then can communicate directly with the object. In C # and VB, it is a combination of both at compile time and at runtime exploring the metadata and taking the information necessary to communicate with the object. The same happens with JavaScript, only the metadata is generated dynamically at runtime, which is actually intended for dynamic programming languages such as JavaScript.

This was a short introduction about Microsoft’s new platform Windows Runtime and a description of what happens when you create an object in it’s environment.

## Delegates in C#

A delegate is a type that references a method. Once a delegate is assigned a method, it behaves exactly like that method. The delegate method can be used like any other method, with parameters and a return value.

Any method that matches the delegate’s signature, which consists of the return type and parameters, can be assigned to the delegate. This makes it possible to programmatically change method calls, and also plug new code into existing classes. As long as you know the delegate’s signature, you can assign your own delegated method.

The first step is defining the delegate:

public delegate int Calculate (int value1, int value2);

Second, we create a method which will be assigned to the delegate object:

public int add(int value1, int value2){
return value1 + value2;
}


Now, we create the delegate object and assigning the method to it:

MyClass mc = new MyClass();


Console.WriteLine("Adding two values: " + add(10, 6));