comparison toolboxes/FullBNT-1.0.7/bnt/examples/dynamic/HHMM/Mgram/mgram1.m @ 0:e9a9cd732c1e tip

first hg version after svn
author wolffd
date Tue, 10 Feb 2015 15:05:51 +0000
parents
children
comparison
equal deleted inserted replaced
-1:000000000000 0:e9a9cd732c1e
1 % a multigram is a degenerate 2HHMM where the bottom level HMMs emit deterministic strings
2 % and the the top level abstract states are independent of each other
3 % cf. HSMM/test_mgram2
4
5 words = {'the', 't', 'h', 'e'};
6 data = 'the';
7 nwords = length(words);
8 word_len = zeros(1, nwords);
9 word_prob = normalise(ones(1,nwords));
10 word_logprob = log(word_prob);
11 for wi=1:nwords
12 word_len(wi)=length(words{wi});
13 end
14 D = max(word_len);
15
16 alphasize = 26;
17 data = letter2num(data);
18 T = length(data);
19
20 % node numbers
21 W = 1; % top level state = word id
22 L = 2; % bottom level state = letter position within word
23 F = 3;
24 O = 4;
25
26 ss = 4;
27 intra = zeros(ss,ss);
28 intra(W,[F L O])=1;
29 intra(L,[O F])=1;
30
31 inter = zeros(ss,ss);
32 inter(W,W)=1;
33 inter(L,L)=1;
34 inter(F,[W L])=1;
35
36 % node sizes
37 ns = zeros(1,ss);
38 ns(W) = nwords;
39 ns(L) = D;
40 ns(F) = 2;
41 ns(O) = alphasize;
42
43
44 % Make the DBN
45 bnet = mk_dbn(intra, inter, ns, 'observed', O);
46 eclass = bnet.equiv_class;
47
48
49
50 % uniform start distrib over words, uniform trans mat
51 Wstart = normalise(ones(1,nwords));
52 Wtrans = mk_stochastic(ones(nwords,nwords));
53
54 % always start in state 1 for each bottom level HMM
55 delta1_start = zeros(1, D);
56 delta1_start(1) = 1;
57 Lstart = repmat(delta1_start, nwords, 1);
58 LRtrans = mk_leftright_transmat(D, 0); % 0 self loop prob
59 Ltrans = repmat(LRtrans, [1 1 nwords]);
60
61 % Finish in the last letter of each word
62 Fprob = zeros(nwords, D, 2);
63 Fprob(:,:,1)=1;
64 for i=1:nwords
65 Fprob(i,length(words{i}),2)=1;
66 Fprob(i,length(words{i}),1)=0;
67 end
68
69 % Each state uniquely emits a letter
70 Oprob = zeros(nwords, D, alphasize);
71 for i=1:nwords
72 for l=1:length(words{i})
73 a = double(words{i}(l))-96;
74 Oprob(i,l,a)=1;
75 end
76 end
77
78
79 % Define CPDs for slice
80 bnet.CPD{eclass(W,1)} = tabular_CPD(bnet, W, 'CPT', Wstart);
81 bnet.CPD{eclass(L,1)} = tabular_CPD(bnet, L, 'CPT', Lstart);
82 bnet.CPD{eclass(F,1)} = tabular_CPD(bnet, F, 'CPT', Fprob);
83 bnet.CPD{eclass(O,1)} = tabular_CPD(bnet, O, 'CPT', Oprob);
84
85 % Define CPDs for slice 2
86 bnet.CPD{eclass(W,2)} = hhmmQ_CPD(bnet, W+ss, 'Fbelow', F, 'startprob', Wstart, 'transprob', Wtrans);
87 bnet.CPD{eclass(L,2)} = hhmmQ_CPD(bnet, L+ss, 'Fself', F, 'Qps', W+ss, 'startprob', Lstart, 'transprob', Ltrans);
88
89 evidence = cell(ss,T);
90 evidence{W,1}=1;
91 sample = cell2num(sample_dbn(bnet, 'length', T, 'evidence', evidence));
92 str = lower(sample(4,:))
93
94 engine = jtree_dbn_inf_engine(bnet);
95 evidence = cell(ss,T);
96 evidence(O,:) = num2cell(data);
97 [engine, ll_dbn] = enter_evidence(engine, evidence);
98
99 gamma = zeros(nwords, T);
100 for t=1:T
101 m = marginal_nodes(engine, [W F], t);
102 gamma(:,t) = m.T(:,2);
103 end
104 gamma
105
106 xidbn = zeros(nwords, nwords);
107 for t=1:T-1
108 m = marginal_nodes(engine, [W F W+ss], t);
109 xidbn = xidbn + squeeze(m.T(:,2,:));
110 end
111
112 % thee
113 % xidbn(1,4) = 0.9412 the->e
114 % (2,3)=0.0588 t->h
115 % (3,4)=0.0588 h-e
116 % (4,4)=0.0588 e-e