locked
Discreet time HMM, with Gaussian Mixture observations RRS feed

  • Question

  • Hi,

    I would like to implement a HMM for a discreet time, gaussian mixture observation model.

    Will support for HMM's (baum welch , viterbi etc.) be included in future releases of infer.net?

    Thank you.

    Thursday, August 9, 2012 1:59 PM

All replies

  • You can already implement such an HMM in Infer.NET by unrolling in time, as explained in a previous forum post. The specific algorithms that you mentioned (Baum-Welch and Viterbi) are unlikely to be included in future releases since you can already use Variational Message Passing for HMMs.
    Tuesday, August 14, 2012 12:45 PM
    Owner
  • Even a simple example that is not a HMM would be very much appreciated, maybe that would already help me to get the basics:

    Lets assume I have one discrete variable with 3 states (A,B,C)

    The state is influenced by the variables state in the last timestep, there is a certain probability for a statechange to either B or C, if the node is currently in state A.


    répertoire de l'article

    Tuesday, August 14, 2012 1:33 PM
  • This question is repeated from an earlier thread, where it was answered.
    Tuesday, August 14, 2012 3:27 PM
    Owner
  • using System;
    using System.Collections.Generic;
    using System.Linq;
    using System.Text;
    using MicrosoftResearch.Infer.Models;
    using MicrosoftResearch.Infer.Maths;
    using MicrosoftResearch.Infer;
    using MicrosoftResearch.Infer.Distributions;
    
    namespace HMM
    {
        class markovmodel
        {
            public Variable<int> NumberOfExamples;
            public Range NoE;
    
            public static int NumofNodes;
            public int NumofChains;
    
            //variables
            public Variable<int> NumofStates;
            public Variable<int> NumofEmission;
            public VariableArray<Vector> Transpro;
            public VariableArray<Vector> Emispro;
    
            public Variable<Vector> no0pro;
    
            public VariableArray<int>[] nodes;
            public VariableArray<int>[] observednodes;
    
            //priors
            public Variable<Dirichlet> no0prior;
    
            public VariableArray<Dirichlet> Transprior; // Observed
            public VariableArray<Dirichlet> Emissionprior; // Observed
    
    
            //posteriors
            public Dirichlet n0posterior;
            public Dirichlet[] TPosterior;
            public Dirichlet[] EPosterior;
    
            //engine
           // public InferenceEngine engine = new InferenceEngine(new VariationalMessagePassing());
            
            public InferenceEngine engine = new InferenceEngine();
            /// <summary>
            /// Creates a HMM model
            /// </summary>
            /// <param name="nodeNum">Number of nodes</param>
            /// <param name="chainsNum">Number of chains</param>
            /// /// <param name="statesN">number of states</param>
            /// /// <param name="EmissionN"> number of appearacnce</param>
            public void CreatModel(int nodeNum, int chainsNum, int statesN, int EmissionN)
            {
    
                NumofNodes = nodeNum;
                NumofChains = chainsNum;
                NumofStates = Variable.New<int>().Named("numberstates");
                NumofEmission = Variable.New<int>().Named("emsissions");
                
                NumberOfExamples = Variable.New<int>().Named("NofE");
                NoE = new Range(NumberOfExamples).Named("N");
                
                
    
    
                NumofStates.ObservedValue = statesN;
                NumofEmission.ObservedValue = EmissionN;
                Range N = new Range(NumofStates);
                Range M = new Range(NumofEmission);
    
                // transition probabilities is one
                Transprior = Variable.Array<Dirichlet>(N).Named("transprior");
                Emissionprior = Variable.Array<Dirichlet>(N).Named("emissionprior");
    
    
                Transpro = Variable.Array<Vector>(N).Named("Trans");
                Transpro[N] = Variable<Vector>.Random<Dirichlet>(Transprior[N]);
                Transpro.SetValueRange(N);
    
    
    
                Emispro = Variable.Array<Vector>(N).Named("Emission");
                Emispro[N] = Variable<Vector>.Random<Dirichlet>(Emissionprior[N]);
                Emispro.SetValueRange(M);
    
    
                nodes = new VariableArray<int>[nodeNum];
                observednodes = new VariableArray<int>[nodeNum];
    
    
                no0prior = Variable.New<Dirichlet>().Named("no0prior");
                no0pro = Variable<Vector>.Random(no0prior).Named("no0pro");
                no0pro.SetValueRange(N);
    
                // prior is a uniform dis
                //defination
                VariableArray<int> node;
                VariableArray<int> emssion_node;
    
                for (int k = 0; k < nodeNum; k++)
                {
                    node = Variable.Array<int>(NoE).Named("node_" + k);
                    emssion_node = Variable.Array<int>(NoE).Named("emission_" + k);
                    node.SetValueRange(N);
                    emssion_node.SetValueRange(M);
                    if (k == 0)
                    {
                        node[NoE] = Variable.Discrete(no0pro).ForEach(NoE);//each observations
                    }
                    else
                    {
                        var v=nodes[k - 1].Range;
                        using(Variable.ForEach(v))
                        using (Variable.Switch(nodes[k - 1][NoE]))
                              node[v].SetTo(Variable.Discrete(Transpro[nodes[k - 1][v]]));
                    }
                    nodes[k] = node;
                    var v2 = nodes[k].Range;
                    using (Variable.ForEach(v2))
                    using (Variable.Switch(nodes[k][NoE]))
                        emssion_node[v2].SetTo(Variable.Discrete(Emispro[nodes[k][v2]]));
    
                    observednodes[k] = emssion_node;
                }
                //Ways to reduce memory consumption
                engine.Compiler.ReturnCopies = false;
                engine.Compiler.OptimiseInferenceCode = true;
                engine.BrowserMode = BrowserMode.Never;
                engine.NumberOfIterations = 50;
                //show factor graph
                //engine.ShowFactorGraph = true;
    
                // specify variables to be infered. avoiding the overhead of computing or caching marginals for any other variables
               // engine.OptimiseForVariables= new VariableArray<int>[] { nodes[2], observednodes[2] };
                
                //more settings
                engine.ShowTimings = true;
                engine.ShowWarnings = true;
            }
    
            /// <summary>
            /// Sample data from a given setting HMM
            /// </summary>
            /// <param name="num">num of samples</param>
            /// <param name="tran">transpo</param>
            /// <param name="n0pro">initial pro</param>
            /// <param name="empro">emssion pro</param>
            /// <returns></returns>
            public int[][] sample(int num, Vector[] transpo, Vector n0pro, Vector[] empro)
            {
                int v = NumofNodes;//numbers of nodes
    
                int[][] sample = new int[num][];
    
                for (int i = 0; i < num; i++)
                    sample[i] = new int[v];
                for (int j = 0; j < num; j++)
                {
                    int n0 = Discrete.Sample(n0pro);
                    int[] temp_observation = new int[v];
                    int current_node = n0;// the value for the hidden node
                    for (int p = 0; p < v; p++)
                    {
                        int apperance = Discrete.Sample(empro[current_node]);
                        int next = Discrete.Sample(transpo[current_node]);
                        current_node = next;
                        temp_observation[p] = apperance;
                    }
                    sample[j] = temp_observation;
                }
                return sample;
    
            }
            public void TrainModel(int[][] sample)
            {
                 NumberOfExamples.ObservedValue = sample.Length;
    
                 Dirichlet no0 = Dirichlet.Uniform(NumofStates.ObservedValue);
                 Dirichlet[] tranP = Enumerable.Repeat(Dirichlet.Uniform(NumofStates.ObservedValue), NumofStates.ObservedValue).ToArray();
    
                 Dirichlet[] emsP = Enumerable.Repeat(Dirichlet.Uniform(NumofEmission.ObservedValue), NumofStates.ObservedValue).ToArray();
    
                 no0prior.ObservedValue = no0;
                 Transprior.ObservedValue = tranP;
                 Emissionprior.ObservedValue = emsP;
    
                for (int i = 0; i < NumofNodes-1; i++)
                {
                    int[] observation = new int[sample.Length];
                    for (int j = 0; j < sample.Length; j++)
                    {
                        observation[j] = sample[j][i];
                    }
                    observednodes[i+1].ObservedValue = observation;//observing data for each node
                }
                
                for (int i = NumofNodes - 1; i >= 0; i--)
                {
                    var t = engine.Infer<Discrete[]>(nodes[i]);
                    Console.WriteLine("infering one by one...." + i + " " + t[0].GetProbs());
                }
    
                n0posterior = engine.Infer<Dirichlet>(no0pro);
                TPosterior = engine.Infer<Dirichlet[]>(Transpro);
                EPosterior = engine.Infer<Dirichlet[]>(Emispro);
    
            }
    
            public void DisplayPosterior()
            {
                Console.WriteLine("node 0 prior " + n0posterior);
                Console.WriteLine("transition matrix");
                for (int i = 0; i < NumofStates.ObservedValue; i++)
                    Console.WriteLine(TPosterior[i].GetMean());
    
                Console.WriteLine("emission matrix");
                for (int i = 0; i < NumofStates.ObservedValue; i++)
                    Console.WriteLine(EPosterior[i].GetMean());
    
            }
    
        }
    
        class Program3
        {
            static void Main(string[] args)
            {
    
                Vector no0pro = Vector.FromArray(1.0, 0, 0.0, 0);
                Vector[] transpro = new Vector[] { Vector.FromArray(0.2, 0.2,0.2,0.4 ) /* */, Vector.FromArray(0.1,0.2,0.3,0.4 ) /* */, 
                        Vector.FromArray(0.1, 0.1,0.2,0.6 ),Vector.FromArray(0.2, 0.2,0.1,0.5 )};
    
                Vector[] emssionpro = new Vector[] { Vector.FromArray(0.3, 0.3,0.4 ) /* */, Vector.FromArray(0.3,0.4,0.3 ) /* */, 
                        Vector.FromArray(0.1, 0.4,0.5 ),Vector.FromArray(0.5,0.4,0.1 )};
                markovmodel model = new markovmodel();
    
                int nodes = 11;
                int chain = 1;
                int interstates = 4;
                int outputs = 3;
                model.CreatModel(nodes, chain, interstates, outputs);
                // 3 nodes, 1 chain, 4 internal states, 3 outputs
                int[] a = { 0 };
                int[] b = { 1 };
                int[] c = { 2 };
                int[] d = { 3 };
                //data sequence
    
                var no0 = Dirichlet.PointMass(no0pro);
                var tranP = transpro.Select(v => Dirichlet.PointMass(v)).ToArray();
                var emsP = emssionpro.Select(v => Dirichlet.PointMass(v)).ToArray();
                model.no0prior.ObservedValue = no0;
                model.Transprior.ObservedValue = tranP;
                model.Emissionprior.ObservedValue = emsP;
                Console.WriteLine(model.Transprior.ObservedValue);
    
                //Dirichlet no0 = Dirichlet.Uniform(4);
                /*
                var no0 = Dirichlet.PointMass(no0pro);
                Dirichlet[] tranP = Enumerable.Repeat(Dirichlet.Uniform(4), 4).ToArray();
                Dirichlet[] emsP = Enumerable.Repeat(Dirichlet.Uniform(3), 4).ToArray();
                */
                model.no0prior.ObservedValue = no0;
                model.Transprior.ObservedValue = tranP;
                model.Emissionprior.ObservedValue = emsP;
                //int[][] ob={a,b};
                
                
                for (int i = 1; i < markovmodel.NumofNodes; i++)
                {
                     if (i == 5)
                    {
                        model.observednodes[i].ObservedValue = c;
                    }
                    else if (i == 8)
                    {
                        model.observednodes[i].ObservedValue = c;
                    }
                    else if (i == 9)
                    {
                        model.observednodes[i].ObservedValue = c;
                    }
                    else if (i % 2 == 0)
                    {
                        model.observednodes[i].ObservedValue = b;
                    }
                    else
                        //model.observednodes[i].ClearObservedValue();
                       //model.nodes[i].ClearObservedValue();
                        model.observednodes[i].ObservedValue = a;
                }
                //model.n0posterior = model.engine.Infer<Dirichlet>(model.no0pro);
                //Console.Write(model.n0posterior.GetMean());
                model.NumberOfExamples.ObservedValue = 1;
    
                System.IO.StreamWriter file = new System.IO.StreamWriter("test.txt");
                for (int i = nodes - 1; i >= 0; i--)
                {
                    var t = model.engine.Infer<Discrete[]>(model.nodes[i]);
                    Console.WriteLine(t[0].GetProbs());
                    file.WriteLine(t[0].GetProbs());
                    
                }
                file.Close();
    
                model.TPosterior = model.engine.Infer<Dirichlet[]>(model.Transpro);
                model.EPosterior = model.engine.Infer<Dirichlet[]>(model.Emispro);
    
                Console.WriteLine(model.TPosterior[0].GetMean());
                Console.WriteLine(model.TPosterior[1].GetMean());
                Console.WriteLine(model.TPosterior[2].GetMean());
                Console.WriteLine(model.TPosterior[3].GetMean());
    
                Console.WriteLine(model.EPosterior[0].GetMean());
                Console.WriteLine(model.EPosterior[1].GetMean());
                Console.WriteLine(model.EPosterior[2].GetMean());
                Console.WriteLine(model.EPosterior[3].GetMean());
    
                Console.WriteLine("Press to continue");
                Console.Read();
    
    
                Console.WriteLine("**************************Train with samples drawn from known HMM");
    
                Vector no0pro1 = Vector.FromArray(0.0, 0, 0.0, 1.0);
                Vector[] transpro1 = new Vector[] { Vector.FromArray(0, 0,0,1 ) /* */, Vector.FromArray(0.1,0.2,0.3,0.4 ) /* */, 
                        Vector.FromArray(0.1, 0.1,0.2,0.6 ),Vector.FromArray(1, 0,0,0 )};
    
                Vector[] emssionpro1 = new Vector[] { Vector.FromArray(1, 0,0 ) /* */, Vector.FromArray(0.3,0.3,0.4 ) /* */, 
                        Vector.FromArray(0.1, 0.4,0.5 ),Vector.FromArray(1,0,0 )};
    
                markovmodel model1 = new markovmodel();
                model1.CreatModel(nodes, chain, interstates, outputs);
                
                //int[][] sample = model1.sample(100, transpro, no0pro, emssionpro);
                // train with uniform prior about n0prior, trans and emission
                
                //read from the same dataset
                int[][] sample = new int[100][];
    
                for (int i = 0; i < 100; i++)
                    sample[i] = new int[10];
    
                System.IO.StreamReader tr = new System.IO.StreamReader("sample.txt");
    
                for (int i = 0; i < 100;i++ ){
                    String s=tr.ReadLine();
                    string[] arrString = s.Split(new char[] { ' ' });
                    for (int j = 0; j < 10; j++)
                        sample[i][j] = int.Parse(arrString[j + 1]);
                }
                tr.Close();
    
                /* write data to  a file
                System.IO.StreamWriter sample_tex = new System.IO.StreamWriter("sample.txt");
                for (int i = 0; i < sample.Length;i++ ){
                    sample_tex.Write("sample"+i+' ');
                    for (int j = 0; j < 10; j++)
                    {
                        sample_tex.Write(sample[i][j]);
                        sample_tex.Write(" ");
                    }
                    sample_tex.Write(" ");
                    sample_tex.WriteLine();
                }
                sample_tex.Close();*/
    
                
    
                model1.TrainModel(sample);
                model1.DisplayPosterior();
    
                Console.Read();
    
            }
        }
    }
    

    This is my implementation. hope it helps.

    just remind that the chain length could not exceed around 60 timesteps in my laptop.

    cheers.

    Tuesday, August 14, 2012 5:37 PM