The Gibbs Sampler

In the previous blog post the details that are needed for an implementation of a Collapsed CRP Gibbs sampler were explicitly presented.
In this post, an implementation in c++ will be presented.
The post will be relatively short as all the technical details that are needed were explicitly given previously.
The code of the c++ implementation can be found here.  Since I will use the package a lot I will regularly commit to fix some bug or add an extra feature.
The videos displays the behaviour of the algorithm for different values of alpha and a different number of iterations.

More specificaly, the hyper parameter alpha is set to 20 the datapoints are given during initialization through a .txt file in the form of [a,b] making it easy to test different datasets. The data were created from 4 different random processes and the task is to find the mixture model that best fits them. The algorithm is set to 100 iterations and the output shows the visualization of cluster assignments to the data until the convergence to 4 clusters. The colormap is used to represent the group each particle belongs and the black circles represent the mean of every cluster. As can be seen, the algorithm converges to three or four clusters which is a reasonable classification given that the data came from 4 different random processes with two of which giving neighbor datapoints. It must be noted that octave uses playback data to correctly visualize the cluster allocation and the sampling even on 100 iterations takes about 1.5 seconds.

Model implementation technical details

Introduction

In this post I will present the details of the resampling filter that was presented in my previous post. More specifically, I will describe the Gibbs sampling algorithm #4 as it was described in Radford M. Neal's paper. This was the basic algorithm that was used in the previous post to augment the resampling process of the FastSLAM algorithm. In Neal's paper the information is presented in a nice and neat but also very abstract manner, so it's difficult for an uninitiated reader to follow all the details. Therefore, information on this post will also take elements from this paper from Herman Kamper as well as this from Xiaodong Yu. In the rest post I will first describe the algorithm as presented by Neal and then I will talk about it's details more explicitly. Finally I will present how these are applied to the aforementioned filter. Again some background on DPMM and MCMC methods is assumed.

Initial description of the collapsed CRP algorithm.

The description of the algorithm taken as is from Neal's paper is the following:

Finally, in a conjugate context, we can often intergrate analytically over the \(\begin{equation} \phi_{c} \end{equation}\), eliminating them from the algorithm. The state of the Markov chain then consists only of the \(\begin{equation} c_{i} \end{equation}\) which we update by Gibbs sampling using the following conditional probabilitites.

\begin{cases}if\  c=c_{j}\ for\ some\ j\neq i:\ P(c_i=c)= b\frac{n_{-i,c}}{n-1+a}\int_{}{}F(y_i,\phi)dH_{-i,c}(\phi) \\\ \ \ \ \ \ \ \ \ \ P(c_i\neq c_j\ for\ all\ i \neq\ j|c_{-i},y_i)=b\frac{a}{n-1+a}\int_{}{}F(y_i,\phi)dG_0(\phi)\end{cases}

Here, \(\begin{equation} H_{-i,c} \end{equation}\) is the posterior distribution of \(\begin{equation} \phi \end{equation}\) based on the prior \(\begin{equation} G_0 \end{equation}\) and all observations \(\begin{equation} y_j \end{equation}\) for which \(\begin{equation} j \neq i \end{equation}\) and \(\begin{equation} c_j =c \end{equation}\).

The third Gibbs sampling method can be summarized as follows:

Algorithm 3: Collapsed Gibbs sampling

---

Let the state of the Markov chain consist of \(\begin{equation} c_1,...,c_n \end{equation}\). Repeatedly sample as follows:

• For \(\begin{equation} i=1,...n \end{equation}\) Draw a new value from \(\begin{equation} c_i|c_{-i},y_i \end{equation}\) as defined by the previous equation.

---

Yes, that's all, no joking here.

A resampling example using matlab.

Introduction

In the previous post, a brief introduction on the Dirichlet distribution, Dirichlet process as well as finite and infinite mixture models was given. Furthermore, methods that can be used to perform inference on such models were enumerated. Finally the research question was revisited and the problem, as well as the approach towards solving it were given.
Before continuing to ROS and a low-level implementation of the model, I wanted to have a working example on matlab taking advantage of its nice vectorized functions.

Current setup

The current setup on matlab is now the following:

• A matalb-based FastSLAM implementation by Grace Lee that was already used in post #1 can be found here
• A matlab-based Gibbs sampler by François Caron can be found here

The particle positions will be used as input to the Gibbs sampler who will cluster the particles as a mixture of Gaussian distributions. The base distribution that will be used as a prior to the Dirichlet process on the sampler will be a Normal-Inverse-Wishart distribution and was chosen due to its conjugacy to the Gaussian. These mixtures output by the sampler will then be used to train a Hierarchical Dirichlet process with the cluster significance proportional to the particle weights that are classified within the cluster and particle significance proportional to the given particle weight.
An example on the single iteration on the resampling step is then the following:

Model description-Warning long post!

Post #3

Introduction

In this post I will describe the basic structure of the filter that I plan on implementing in the resampling step of the FastSLAM .The post is structured as follows:
First a brief introduction regarding Dirichlet processes as well as the basic finite and infinite mixture models using such processes are given. Then the basic inference methods on finite and infinite mixture models are presented. Finally the relation that connects such models to my research topic is analysed.
Reading material regarding DP's can be found here and here. Wikipedia pages have useful but long articles regarding the Dirichlet Distribution and Dirichlet Process and their respective properties(exchangeability, aggregation etc) .

The Dirichlet process, a fast introduction.

Dirichlet distribution: Let :

Be a random pmf.  In addition, suppose that 

First week

Blog update #2.
So, the first week into the intenrship was, like all first weeks, an introductory one; administrative stuff, getting to know the people and the workplace. Really nice working environment and very friendly people.
Now, as far as the internship, I sticked to the internship planning and this week was mostly focused on research regarding the general theories of my topic.
The research done could be divided into two main areas:

1. ROS and the turtlebot simulator
2. The Dirichlet distribution and its properties

Regarding the first one, I will be given a desktop that has Ubuntu 12.04 as the turtlebot simulator is not yet compatible with the 14.04 Ubuntu due to conflicts on .urdf files. As a result, the machine will be running on ROS hydro as it is the most recent working simulator until a patch on 14.04 and ROS indigo is given.

On a sidenote, I also installed Ubuntyu 12.04  and ROS Hydro on my old laptop but the simulator there is too slow so I will use the company's dekstop for now.

Regarding the ROS simulator, I re-did all the tutorials on the ROS wiki to remember how topics work. Furthermore, I read through the code of the gmapping package, which is the built in tool for 2D-SLAM in ROS in order to have a deeper understanding of how it works and what should be changed to implement my alterations when the time comes.

The basic websites I used are the following:

• http://wiki.ros.org/ 
• http://answers.ros.org/questions/

First Post

This blog serves as a log regarding the progress of my internship at DoBots company in Rotterdam.
The blog is scheduled to be updated weekly or bi-weekly if more progress is made and more data needs to be listed.
The internship will start on 22/09/2014.
Cheers.