Computational Intelligence, SS08
2 VO 442.070 + 1 RU 708.070 last updated:
Course Notes (Skriptum)
Online Tutorials
Practical Course Slides
Animated Algorithms
Interactive Tests
Key Definitions
Literature and Links

Homework 8: Hidden Markov Models

[Points: 12.5; Issued: 2008/05/30; Deadline: 2008/06/13; Tutor: Susanne Rexeis; Infohour: 2008/06/06, 15:30-16:30, HS i11; Einsichtnahme: 2008/06/27, 15:30-16:30, HS i11; Download: pdf; ps.gz]

Go through the tutorial ``Hidden Markov Models'' and download the accompanying MATLAB programs and data. The tasks issued in this homework are a subset of the tasks in the tutorial.

You can use a printout of this homework designation as a basis for your elaboration, and fill in your results in the tables. Add sheets of paper with descriptions of your observations and printouts of MATLAB figures. Add to the end of the report the MATLAB script and your Viterbi decoding function vit_ghmm you programmed for producing all the results.

Additionally, you should send the matlab code of the viterbi function for generating the results to the tutor ( with the subject EW. Please provide your Name and Matrikelnummer in the scripts and in the body of the email.

  1. Load the Hidden Markov Models (HMMs) $ \ensuremath\boldsymbol{\Theta}_1,\ldots,\ensuremath\boldsymbol{\Theta}_4$, and make a sketch of each of the models with states and transition probabilities. The parameters of the Markov models and of the Gaussian emission pdfs are stored in the file data.mat. Each HMM $ \ensuremath\boldsymbol{\Theta}_i$ is stored as an object, e.g., hmm1, with fields hmm1.trans, hmm1.pi, hmm1.means, and hmm1.vars. The trans field contains the transition matrix $ \mathbf{A}$, and the pi field the prior probability vector $ \boldsymbol{\pi}$. The means field contains a matrix composed of mean vectors of the Gaussian emission pdfs, where each column of the matrix corresponds to one state of the HMM (to access the mean vector $ \ensuremath\boldsymbol{\mu}$ of, e.g., the second state from hmm1 use: hmm1.means(:,2). The vars field contains a 3 dimensional array of covariance matrices, where the third dimension corresponds to the state (e.g., to access the covariance matrix $ \ensuremath\boldsymbol{\Sigma}$ of state 1 use hmm1.vars(:,:,1)).
  2. Generate samples from the HMMs hmm1, hmm2, hmm3, and hmm4 and plot them with plotseq and plotseq2. In the resulting plots, the obtained sequences are represented by a yellow line where each point is overlaid with a colored dot. The different colors of the dots indicate the state from which a particular sample has been drawn.

    » % Example: generate a sequence of length T (e.g. 80) from HMM1

    » [X,ST] = sample_ghmm(hmm1,T)

    » plotseq(X,ST) % View of both dimensions as separate sequences

    » plotseq2(X,ST,hmm) % 2D view with location of Gaussian states

    Draw several sequences for each HMM and compare. Compare the MATLAB figures with your sketch of the models and add (some of) them to your homework elaboration. Moreover, explain: What is the effect of the different transition matrices of the HMMs on the sequences obtained? Hence, what is the role of the transition probabilities in the HMM?

  3. Pattern recognition with HMM's: Classify the sequences $ X_1, X_2, X_3,X_4$, given in the file Xdata.mat, in a maximum likelihood sense with respect to the four Markov models $ \ensuremath\boldsymbol{\Theta}_1, \ensuremath\boldsymbol{\Theta}_2, \ensuremath\boldsymbol{\Theta}_3$, and $ \ensuremath\boldsymbol{\Theta}_4$. Use the function loglik_ghmm to compute the log-likelihood $ \log p(X_i\vert\ensuremath\boldsymbol{\Theta}_j)$ of a sequence $ X_i$ with respect to a HMM $ \ensuremath\boldsymbol{\Theta}_j$. Store the results in a matrix (they will be used in the next section).

    » load Xdata

    » % Example:

    » logLike(1,1) = loglik_ghmm(X1,hmm1)

    » logLike(1,2) = loglik_ghmm(X1,hmm2)


    » logLike(i,j) = loglik_ghmm(Xi,hmmj)


    Filling the logLike matrix can be done automatically with the help of loops:

    >> for i=1:4,
         for j=1:4,
          stri = num2str(i);
          strj = num2str(j);
          eval([ 'logLike(' , stri , ',' , strj , ')=...
                loglik_ghmm(X' , stri , ',hmm' , strj , ');' ]);
    You can find the maximum of each row in the matrix with the MATLAB function max:
    >> for i=1:4;
        [v,index] = max(logLike(i,:));
        disp(['X',num2str(i),' -> HMM',num2str(index)]);
    Sequence $ \log p(X_i\vert\ensuremath\boldsymbol{\Theta}_1)$ $ \log p(X_i\vert\ensuremath\boldsymbol{\Theta}_2)$ $ \log p(X_i\vert\ensuremath\boldsymbol{\Theta}_3)$ $ \log p(X_i\vert\ensuremath\boldsymbol{\Theta}_4)$ Most likely model
    $ X_1$          
    $ X_2$          
    $ X_3$          
    $ X_4$          
  4. Viterbi decoder: Write a MATLAB function [loglik,path] = vit_ghmm(data,HMM) to implement the Viterbi decoding algorithm to find the most likely state sequence $ Q$ for a given observation sequence $ X_i$ for HMMs with Gaussian emission probabilities. Use the function mk_ghmm_obs_lik to calculate the observation probabilities (Gaussian) for each state and time step. Please perform the calculations in the log domain, where the multiplications of the probabilities for the parameters $ \delta$ and $ \psi$ become additions. So you can prevent numerical issues for long sequences.
    [loglik,path] = vit_ghmm(data,HMM)
    % Compute the path and the log likelihood of a given model and
    % observation sequence
    % INPUT
    %      data ... containing a sequence of observation
    %         size(data)=[Number of features of each obs., length of sequence]
    %      HMM is an object containing
    %       HMM.trans = state transition probability matrix;
    %       HMM.pi = prior probability vector;
    %       HMM.means = mean vectors of Gaussian emission pdf for each state;
    %       HMM.vars = covariance matrices of Gaussian em. pdf for each state;
    % OUTPUT
    %       loglik ... log likelihood of the most likely path for the data
    %       path ... most likely path
  5. Use your function vit_ghmm to compute the most likely paths for the sequences $ X_1,\ldots X_4$ with respect to each model $ \ensuremath\boldsymbol{\Theta}_1,\ldots,\ensuremath\boldsymbol{\Theta}_4$. Also compute the log-likelihoods $ \log p^*(X_i\vert\ensuremath\boldsymbol{\Theta}_j)$ along the most likely paths found by the Viterbi decoder. Note down your results below. Compare with the log-likelihoods $ \log p(X_i\vert\ensuremath\boldsymbol{\Theta}_j)$ obtained in the previous section with the function loglik_ghmm(...):

    » diffL = logLike-logLikeViterbi

    Log-likelihoods along the best path:

    Sequence $ \log p^*(X_i\vert\ensuremath\boldsymbol{\Theta}_1)$ $ \log p^*(X_i\vert\ensuremath\boldsymbol{\Theta}_2)$ $ \log p^*(X_i\vert\ensuremath\boldsymbol{\Theta}_3)$ $ \log p^*(X_i\vert\ensuremath\boldsymbol{\Theta}_4)$ Most likely model
    $ X_1$          
    $ X_2$          
    $ X_3$          
    $ X_4$          

    Difference between log-likelihoods of a sequence given the model and log-likelihoods along the best path found by the Viterbi algorithm, $ \log p(X\vert\ensuremath\boldsymbol{\Theta}_i) - \log p^*(X\vert\ensuremath\boldsymbol{\Theta}_i)$:

    Sequence HMM1 HMM2 HMM3 HMM4
    $ X_1$        
    $ X_2$        
    $ X_3$        
    $ X_4$        

    Is the likelihood along the best path a good approximation of the real likelihood of a sequence given a model ?