# Implementing EM clustering for a mixture markov model

I have a mixture Markov model (containing K clusters, or components) that I am trying to train, e.g perform clustering over a set of varying length sequences. Each component of the model is a first order Markov chain. The model is defined by the following:

where

is the model describing behaviour of sequences in the k-th cluster.

Each kth model component(cluster) has initial probability distribution vector (theta_KI) and a transition matrix (theta_tK) - see equation 2. Each sequence seq_i in the set of sequences has a different length L (so the L in the 2nd equation would be different for each sequence in the set).

There is n possible visible states in the model and m sequences in each cluster, and the total amount of sequences is M.

In this case, how would I define e- and m-steps for my EM algorithm? My understanding is that I have to maximize log likelihood function of the model, which would be a sum of likelihoods across all vectors v_i ( synonym with sequence seq_i)- something like:

However, I get stuck not knowing if this is the right definition for what needs to be maximized and which arguments need to be found for that. I have looked for literature but it is very sparse, almost nonexistent.

Another question: is there actually an R library that has functions that can be extended to include such a model for EM algorithm?

• Please define the notation. – AdamO Aug 2 '13 at 15:36

As you say, the EM is an algorithm aiming at maximizing the log-likelihood of your model. The fact that the sequences have different lengths is not a problem. To see why let us first assume that you have trained your model. What does this mean?. That you have estimates for the transition probabilities, $\hat{p}(v_{i}|v_{i-1};\theta_{K})$ and for the priors $\hat{p}(c_{K};\theta_{K})$. How do you perform classification in such a model? You take your sequence, for example $v_{2}, v_{7}, v_{12}$ and calculate,

$$\mathbf{argmax}_{K} \hat{p}(c_{K};\theta_{K})p(v_{2};\theta_{K}) \hat{p}(v_{7}|v_{2};\theta_{K}) \hat{p}(v_{12}|v_{7};\theta_{K})$$

Hence, whether you have longer or shorter sequences does not change anything substantially. You need the same amount of parameters for each cluster. Which, without further simplification, would be the number of free parameters for the transition matrix, $N*(N-1)$, plus 1 for the prior.

The EM algorithm would be as follows: start with some values for the transition matrix (you could try generating a random matrix such that corresponds to a transition matrix) and for the $\theta_{i}'s$.

1. Expectation: estimate the priors. You may do this by calculating the maximum likelihood estimator, i.e. the number of sequences belonging to each cluster divided by the total number of sequences. A sequence belongs to cluster $K$ if $K = \mathbf{argmax}_{j}p(v|c_{j}, \theta)$.
2. Maximization: based on the resulting assignment of sequences to clusters, estimate the transition probabilities. Here, some smoothing technique, like Laplace smoothing might be necessary in order to ensure numerical stability.

You go on like this until convergence, that is, until the value of the likelihood does not further improve.

Regretfully I am not familiar with R.

• Thank you @juampa for a thorough description. Just to clarify, in step 1 - by priors do you mean the probability that a sequence belongs to cluster c_k? Secondly, how would I calculate argmax in each step? Is by argmax meant the parameters of the given cluster? ( transition probabilities. is there a method to estimate this? – zima Aug 7 '13 at 15:17
• 1. the probability that the sequence is generated by component $k$ of the mixture, 2. for each cluster you have a set of sequences, from which you estimate the transition probabilities. Notice that for each component, you get different estimattes. When presented a sequence, you calculate the probability of that sequence for each of the clusters, i.e. the probability that the sequence has been generated by each of the components (as in the expression above), and pick up the index corresponding to the component yielding highest probability. – jpmuc Aug 7 '13 at 16:10
• thank you again for the answer. Perhaps you could suggest a tutorial on Expectation Maximization with regards to Mixed Markov Models with first order Markov Chains? Also, the article that I am using as a guidance, link here, seems to initialize the model including the assumption that there are k sequences in each cluster. I am not completely understanding how is such a model initialized. Is the first step parameters estimated from the complete dataset as if there is 1 cluster? How would you initialize such a model? – zima Aug 9 '13 at 10:52