Logiciels / Software

Chainxem, Valxem, Treexem

Il s'agit d'un ensemble de routines pour l'identification des modèles de chaîne et d'arbre de Markov cachés. Ces routines nécessitent Matlab et reprennent l'interface du logiciel XEMgaus développé par Christophe Biernacki et dont la refonte a abouti à MIXMOD, dédié à l'identification mélanges gaussiens indépendants. Les procédures ci-dessous étendent les fonctionnalités de XEMGaus, puisque Chaixem est dédié à l'identification des chaînes de Markov cachées, Treexem à celle des arbres binaires de Markov cachés et Valxem à celle des chaînes de Markov cachés à observations manquantes.

Ces logiciels permettent le calcul de probabilités, la restauration des états cachés et l'estimation des paramètres pour les modèles ci-dessus. Il est également possible de simuler des processus. L'estimation peut se faire par l'algorithme EM, CEM, SEM et EM à la Gibbs ou suivant toute combinaison séquentielle de ces algorithmes. L'utilisateur a la possibilité de choisir les lois d'émission parmi les 14 modèles gaussiens, la loi exponentielle ou les lois non paramétriques pour des observations discrètes. L'utilisateur peut également imposer des contraintes sur la matrice de transition (matrice de type diagonale par bandes) ou sur la loi de l'état initial (loi stationnaire).

Les logiciels développés reprennent les mêmes fonctionnalités que XEMGaus pour l'estimation, à savoir :

choix de la valeur initiale du paramètre ou détermination par un ou plusieurs tirages aléatoires ;
choix de la condition d'arrêt : croissance relative de la log-vraisemblance, nombre d'itérations ou arrêt lorsque l'une des deux conditions est satisfaite ;
prise en compte des états connus lorsqu'il y en a ;
estimation des paramètres par la moyenne (éventuellement en supprimant les valeurs des premières itérations) ou par la valeur maximisant la vraisemblance quand des algorithmes stochastiques sont utilisés.

De plus, plusieurs réalisations indépendantes du modèle (éventuellement avec un nombre de données différent pour chacune) peuvent être utilisées pour identifier un même modèle.

Exemple d'utilisation

L'exemple suivant illustre l'utilisation du logiciel Chainxem pour la simulation de deux séquences, de longueur 500 et 600, suivant un modèle de chaîne de Markov cachée stationnaire à deux états cachés, de matrice de transition $P = \left[ \begin{array}{cc} 0,8 & 0,2 \\ 0,4 & 0,6 \end{array}\right],$ de loi stationnaire $\pi = \left[\begin{array}{cc} 0,6667 & 0,3333 \end{array} \right],$ de lois d'émission gaussiennes de moyennes respectives $\mu_2 = 3$ et $\mu_1 = -3$ et de variances respectives $\sigma^2_2 = 0,5$ et $\sigma^2_1 = 1$ .

Les paramètres du modèle sont alors estimés par l'algorithme EM à la Gibbs suivi de l'algorithme EM, en utilisant simultanément les deux séquences simulées.

% simulation d'une chaine de Markov cachee a deux etats

» model.fam = 'normal';

% lois d'emissions gaussiennes

» par.A = [ 0.8000 0.2000; 0.4000 0.6000 ];

% matrice de transition

» par.p = [ 0.6667 0.3333];

% distribution stationnaire de A

» par.mu = [-3; 3];

% les deux parametres d'emission : moyenne des lois gaussiennes

» par.S(:,:,1) = 1;

» par.S(:,:,2) = 0.5;

% les deux parametres d'emission : variance des lois gaussiennes

» [x,z] = chainxrnd([500 600], par,model);

» size(x{1})

ans =

500 1

% x{1} est une sequence simulee de longueur 500

% z{1} est la chaine de Markov correspondante

% estimation des parametres

» model.k=2;

% nombre d'etats caches du modele

» model.model='plkI';

% modele stationnaire, gaussien, avec une variance

% dependante de l'etat cache

» algo = {'gibbs','em'};

% algorithmes utilises : EM a la Gibbs puis EM

» cvg = 'xmlORmaxit';

% critere d'arret : stabilisation de la log vraisemblance ou

% depassement du nombre maximal d'iterations

» maxit = 300;

% nombre maximal d'iterations pour l'algorithme EM a la Gibbs

% puis pour EM

» parinit.A = 'random';

» parinit.p = 'random';

» parinit.mu = 'random';

» parinit.S = 'rand. var.';

% initialisation aleatoire de l'algorithme

» nbxem = 3;

% nombre de valeurs initiales pour l'algorithme

» bestpar = chainxem('x',x,'model',model,'algo',algo,'cvg',cvg,...

'maxit',maxit,'nbxem',nbxem,'parinit',parinit);

|------------------------------|
| 2 state(s) - model plkI |
|------------------------------|

| xem 1: [gibbs.9][em.3]

| xem 2: [gibbs.100.147][em.3]

| xem 3: [gibbs.5][em.3]

% estimation des parametres par l'algorithme EM a la Gibbs puis EM

% le nombre d'iterations effectuees est indiquee pour chaque execution

% de EM avec un parametre initial different

» bestpar.A

ans =

0.6173 0.3827

0.1939 0.8061

% estimateur de la matrice de transition

» bestpar.p

ans =

0.3363 0.6637

% estimateur de la matrice de la loi stationnaire

» bestpar.mu

ans =

2.9874

-3.0489

% estimateur des moyennes des parametres d'emission

» bestpar.S

ans(:,:,1) =

0.4586

ans(:,:,2) =

0.9289

% estimateur des variances des parametres d'emission

% modele a lois d'emission "non parametriques" (autrement dit multinomiales) % lignes : etats, colonnes : valeurs

par.O = [0.1 0.1 0.4 0.4; 0.3 0.4 0.2 0.1]; model.fam = 'np'; [x,z] = chainxrnd([500 600], par,model);

Instructions de téléchargement et utilisation

Les sources sont disponibles comme archive Linux compressée. L'aide en ligne est obtenue par les commandes

help pour l'aide hypertexte si les fichiers Contents.m sont présent dans le chemin Matlab;
help chainxem pour l'identification des chaînes de Markov cachées;
help valxem pour l'identification des chaînes de Markov cachées à observations manquantes;
help treexem pour l'identification des arbres de Markov cachées binaires.

Retour en début de page.

Chainxem, Valxem, Treexem

This software consists of a set of routines for the identification of hidden Markov chains and hidden Markov trees. These routines require Matlab and their interface is similar to that of the XEMgaus software, developed by Christophe Biernacki. An improved version of XEMgaus, called MIXMOD, is dedicated to the identification of independent mixture models. The procedures below extend the features of XEMGaus since Chaixem addresses the identification of hidden Markov chains, Treexem that of binary hidden Markov trees and Valxem that of hidden Markov chains with missing observations.

The computation of probabilities (typically the likelihood), the hidden state restoration and the parameter estimation can be performed for the models above. These models can also be simulated. The estimation is achieved by the EM algorithm or any combination of its variants: CEM (Classification EM), SEM (Stochastic EM) and EM à la Gibbs. The following observation distributions can be chosen: the 14 Gaussian models, the exponential distribution or nonparametric (i.e. multinomial) distributions for observed processes with finite values. The user can also define some constraints on the transition probability matrix (band-diagonal matrix) or on the initial state distribution (stationary distribution).

These routines have the same features than XEMGaus concerning the estimation, namely:

choice of the initial value of the parameter or random determination (possibility of using several initial random values);
choice of the stopping criterion: relative increase of the likelihood, number of iterations or combination of both criterion;
known states, if any, are taken into account in the estimation;
concerning stochastic versions of EM, parameter estimation can be based on the mean (possibly with deletion of the values obtained through the first iterations) or on the likelihood maximization.

Moreover, several independent processes (possibly with different sizes) can be used for the identification of a single model.

Example

The following example illustrates a possible use of the Chainxem software. We show how to simulate two hidden Markov chains of respective lengths 500 and 600 in the case of a two-state stationary model with transition probability matrix $P = \left[ \begin{array}{cc} 0.8 & 0.2 \\ 0.4 & 0.6 \end{array}\right],$ with stationary distribution $\pi = \left[\begin{array}{cc} 0.6667 & 0.3333 \end{array} \right],$ and Gaussian emission distributions with means $\mu_2 = 3$ and $\mu_1 = -3$ and variances $\sigma^2_2 = 0.5$ and $\sigma^2_1 = 1$ .

Then the model parameter are estimated using the EM algorithm à la Gibbs followed by the classical EM algorithm from the two simulated sequences:

% simulation of a two-state hidden Markov chain

» model.fam = 'normal';

% Gaussian emission distributions

» par.A = [ 0.8000 0.2000; 0.4000 0.6000 ];

% transition probability matrix

» par.p = [ 0.6667 0.3333];

% stationary distribution associated with A

» par.mu = [-3; 3];

% parameters of the emission distributions: mean of the Gaussian distributions

» par.S(:,:,1) = 1;

» par.S(:,:,2) = 0.5;

% parameters of the emission distributions: variance of the Gaussian distributions

» [x,z] = chainxrnd([500 600], par,model);

» size(x{1})

ans =

500 1

% x{1} is a simulated sequence with length 500

% z{1} is the associated Markov chain

% parameter estimation

» model.k=2;

% number of hidden states

» model.model='plkI';

% stationary Gaussian model with state-dependent variance

» algo = {'gibbs','em'};

% chosen algorithms: EM a la Gibbs and EM

» cvg = 'xmlORmaxit';

% stopping criterion: stabilization of the log-likelihood or

% maximal number of iterations

» maxit = 300;

% maximal number of iterations concerning the EM algorithm a la Gibbs

% and then the EM algorithm

» parinit.A = 'random';

» parinit.p = 'random';

» parinit.mu = 'random';

» parinit.S = 'rand. var.';

% random determination of the initial parameters

» nbxem = 3;

% number of initial values of the parameter

» bestpar = chainxem('x',x,'model',model,'algo',algo,'cvg',cvg,...

'maxit',maxit,'nbxem',nbxem,'parinit',parinit);

|------------------------------|
| 2 state(s) - model plkI |
|------------------------------|

| xem 1: [gibbs.9][em.3]

| xem 2: [gibbs.100.147][em.3]

| xem 3: [gibbs.5][em.3]

% parameter estimation using EM a la Gibbs followed by EM

% for each run of the algorithm, the effective number of iterations is printed

% (each run corresponds to a different initial parameter)

» bestpar.A

ans =

0.6173 0.3827

0.1939 0.8061

% estimate of the transition probability matrix

» bestpar.p

ans =

0.3363 0.6637

% estimate of the stationary distribution

» bestpar.mu

ans =

2.9874

-3.0489

% estimate of the emission parameters: means of the Gaussian distributions

» bestpar.S

ans(:,:,1) =

0.4586

ans(:,:,2) =

0.9289

% estimate of the emission parameters: variances of the Gaussian distributions

% nonparametric (i.e. multinomial) emission distributions % rows: states, columns: values

par.O = [0.1 0.1 0.4 0.4; 0.3 0.4 0.2 0.1]; model.fam = 'np'; [x,z] = chainxrnd([500 600], par,model);

Instructions for download and use

The sources are available as a zipped Linux archive. The online Matlab help is obtained by the commands

help for hypertext help, if the Contents.m files are contained in the Matlab path variable;
help chainxem for the identification of hidden Markov chains;
help valxem for the identification of hidden Markov chains with missing observations;
help treexem for the identification of binary hidden Markov trees.

Back to the top of the page (english part).

Back to my homepage (in French).

Email :

Laboratoire Jean Kuntzmann