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.
