Mercurial > hg > camir-aes2014

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

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
first hg version after svn
author wolffd
date Tue, 10 Feb 2015 15:05:51 +0000
parents
children
rev   line source
wolffd@0 1 % Like a durational HMM, except we use soft evidence on the observed nodes.
wolffd@0 2 % Should give the same results as HSMM/test_mgram2.
wolffd@0 3
wolffd@0 4 past = 1;
wolffd@0 5 % If past=1, P(Yt|Qt=j,Dt=d) = P(y_{t-d+1:t}|j)
wolffd@0 6 % If past=0, P(Yt|Qt=j,Dt=d) = P(y_{t:t+d-1}|j) - future evidence
wolffd@0 7
wolffd@0 8 words = {'the', 't', 'h', 'e'};
wolffd@0 9 data = 'the';
wolffd@0 10 nwords = length(words);
wolffd@0 11 word_len = zeros(1, nwords);
wolffd@0 12 word_prob = normalise(ones(1,nwords));
wolffd@0 13 word_logprob = log(word_prob);
wolffd@0 14 for wi=1:nwords
wolffd@0 15 word_len(wi)=length(words{wi});
wolffd@0 16 end
wolffd@0 17 D = max(word_len);
wolffd@0 18
wolffd@0 19
wolffd@0 20 alphasize = 26*2;
wolffd@0 21 data = letter2num(data);
wolffd@0 22 T = length(data);
wolffd@0 23
wolffd@0 24 % node numbers
wolffd@0 25 W = 1; % top level state = word id
wolffd@0 26 L = 2; % bottom level state = letter position within word
wolffd@0 27 F = 3;
wolffd@0 28 O = 4;
wolffd@0 29
wolffd@0 30 ss = 4;
wolffd@0 31 intra = zeros(ss,ss);
wolffd@0 32 intra(W,[F L O])=1;
wolffd@0 33 intra(L,[O F])=1;
wolffd@0 34
wolffd@0 35 inter = zeros(ss,ss);
wolffd@0 36 inter(W,W)=1;
wolffd@0 37 inter(L,L)=1;
wolffd@0 38 inter(F,[W L O])=1;
wolffd@0 39
wolffd@0 40 % node sizes
wolffd@0 41 ns = zeros(1,ss);
wolffd@0 42 ns(W) = nwords;
wolffd@0 43 ns(L) = D;
wolffd@0 44 ns(F) = 2;
wolffd@0 45 ns(O) = alphasize;
wolffd@0 46 ns2 = [ns ns];
wolffd@0 47
wolffd@0 48 % Make the DBN
wolffd@0 49 bnet = mk_dbn(intra, inter, ns, 'observed', O);
wolffd@0 50 eclass = bnet.equiv_class;
wolffd@0 51
wolffd@0 52 % uniform start distrib over words, uniform trans mat
wolffd@0 53 Wstart = normalise(ones(1,nwords));
wolffd@0 54 Wtrans = mk_stochastic(ones(nwords,nwords));
wolffd@0 55 %Wtrans = ones(nwords,nwords);
wolffd@0 56
wolffd@0 57 % always start in state d = length(word) for each bottom level HMM
wolffd@0 58 Lstart = zeros(nwords, D);
wolffd@0 59 for i=1:nwords
wolffd@0 60 l = length(words{i});
wolffd@0 61 Lstart(i,l)=1;
wolffd@0 62 end
wolffd@0 63
wolffd@0 64 % make downcounters
wolffd@0 65 RLtrans = mk_rightleft_transmat(D, 0); % 0 self loop prob
wolffd@0 66 Ltrans = repmat(RLtrans, [1 1 nwords]);
wolffd@0 67
wolffd@0 68 % Finish when downcoutner = 1
wolffd@0 69 Fprob = zeros(nwords, D, 2);
wolffd@0 70 Fprob(:,1,2)=1;
wolffd@0 71 Fprob(:,2:end,1)=1;
wolffd@0 72
wolffd@0 73
wolffd@0 74 % Define CPDs for slice 1
wolffd@0 75 bnet.CPD{eclass(W,1)} = tabular_CPD(bnet, W, 'CPT', Wstart);
wolffd@0 76 bnet.CPD{eclass(L,1)} = tabular_CPD(bnet, L, 'CPT', Lstart);
wolffd@0 77 bnet.CPD{eclass(F,1)} = tabular_CPD(bnet, F, 'CPT', Fprob);
wolffd@0 78
wolffd@0 79
wolffd@0 80 % Define CPDs for slice 2
wolffd@0 81 bnet.CPD{eclass(W,2)} = hhmmQ_CPD(bnet, W+ss, 'Fbelow', F, 'startprob', Wstart, 'transprob', Wtrans);
wolffd@0 82 bnet.CPD{eclass(L,2)} = hhmmQ_CPD(bnet, L+ss, 'Fself', F, 'Qps', W+ss, 'startprob', Lstart, 'transprob', Ltrans);
wolffd@0 83
wolffd@0 84
wolffd@0 85 if 0
wolffd@0 86 % To test it is generating correctly, we create an artificial
wolffd@0 87 % observation process that capitalizes at the start of a new segment
wolffd@0 88 % Oprob(Ft-1,Qt,Dt,Yt)
wolffd@0 89 Oprob = zeros(2,nwords,D,alphasize);
wolffd@0 90 Oprob(1,1,3,letter2num('t'),1)=1;
wolffd@0 91 Oprob(1,1,2,letter2num('h'),1)=1;
wolffd@0 92 Oprob(1,1,1,letter2num('e'),1)=1;
wolffd@0 93 Oprob(2,1,3,letter2num('T'),1)=1;
wolffd@0 94 Oprob(2,1,2,letter2num('H'),1)=1;
wolffd@0 95 Oprob(2,1,1,letter2num('E'),1)=1;
wolffd@0 96 Oprob(1,2,1,letter2num('a'),1)=1;
wolffd@0 97 Oprob(2,2,1,letter2num('A'),1)=1;
wolffd@0 98 Oprob(1,3,1,letter2num('b'),1)=1;
wolffd@0 99 Oprob(2,3,1,letter2num('B'),1)=1;
wolffd@0 100 Oprob(1,4,1,letter2num('c'),1)=1;
wolffd@0 101 Oprob(2,4,1,letter2num('C'),1)=1;
wolffd@0 102
wolffd@0 103 % Oprob1(Qt,Dt,Yt)
wolffd@0 104 Oprob1 = zeros(nwords,D,alphasize);
wolffd@0 105 Oprob1(1,3,letter2num('t'),1)=1;
wolffd@0 106 Oprob1(1,2,letter2num('h'),1)=1;
wolffd@0 107 Oprob1(1,1,letter2num('e'),1)=1;
wolffd@0 108 Oprob1(2,1,letter2num('a'),1)=1;
wolffd@0 109 Oprob1(3,1,letter2num('b'),1)=1;
wolffd@0 110 Oprob1(4,1,letter2num('c'),1)=1;
wolffd@0 111
wolffd@0 112 bnet.CPD{eclass(O,2)} = tabular_CPD(bnet, O+ss, 'CPT', Oprob);
wolffd@0 113 bnet.CPD{eclass(O,1)} = tabular_CPD(bnet, O, 'CPT', Oprob1);
wolffd@0 114
wolffd@0 115 evidence = cell(ss,T);
wolffd@0 116 %evidence{W,1}=1;
wolffd@0 117 sample = cell2num(sample_dbn(bnet, 'length', T, 'evidence', evidence));
wolffd@0 118 str = num2letter(sample(4,:))
wolffd@0 119 end
wolffd@0 120
wolffd@0 121
wolffd@0 122 if 1
wolffd@0 123
wolffd@0 124 [log_obslik, obslik, match] = mk_mgram_obslik(lower(data), words, word_len, word_prob);
wolffd@0 125 % obslik(j,t,d)
wolffd@0 126 softCPDpot = cell(ss,T);
wolffd@0 127 ens = ns;
wolffd@0 128 ens(O)=1;
wolffd@0 129 ens2 = [ens ens];
wolffd@0 130 for t=2:T
wolffd@0 131 dom = [F W+ss L+ss O+ss];
wolffd@0 132 % tab(Ft-1, Q2, Dt)
wolffd@0 133 tab = ones(2, nwords, D);
wolffd@0 134 if past
wolffd@0 135 tab(1,:,:)=1; % if haven't finished previous word, likelihood is 1
wolffd@0 136 %tab(2,:,:) = squeeze(obslik(:,t,:)); % otherwise likelihood of this segment
wolffd@0 137 for d=1:min(t,D)
wolffd@0 138 tab(2,:,d) = squeeze(obslik(:,t,d));
wolffd@0 139 end
wolffd@0 140 else
wolffd@0 141 for d=1:max(1,min(D,T+1-t))
wolffd@0 142 tab(2,:,d) = squeeze(obslik(:,t+d-1,d));
wolffd@0 143 end
wolffd@0 144 end
wolffd@0 145 softCPDpot{O,t} = dpot(dom, ens2(dom), tab);
wolffd@0 146 end
wolffd@0 147 t = 1;
wolffd@0 148 dom = [W L O];
wolffd@0 149 % tab(Q2, Dt)
wolffd@0 150 tab = ones(nwords, D);
wolffd@0 151 if past
wolffd@0 152 %tab = squeeze(obslik(:,t,:));
wolffd@0 153 tab(:,1) = squeeze(obslik(:,t,1));
wolffd@0 154 else
wolffd@0 155 for d=1:min(D,T-t)
wolffd@0 156 tab(:,d) = squeeze(obslik(:,t+d-1,d));
wolffd@0 157 end
wolffd@0 158 end
wolffd@0 159 softCPDpot{O,t} = dpot(dom, ens(dom), tab);
wolffd@0 160
wolffd@0 161
wolffd@0 162 %bnet.observed = [];
wolffd@0 163 % uniformative observations
wolffd@0 164 %bnet.CPD{eclass(O,2)} = tabular_CPD(bnet, O+ss, 'CPT', mk_stochastic(ones(2,nwords,D,alphasize)));
wolffd@0 165 %bnet.CPD{eclass(O,1)} = tabular_CPD(bnet, O, 'CPT', mk_stochastic(ones(nwords,D,alphasize)));
wolffd@0 166
wolffd@0 167 engine = jtree_dbn_inf_engine(bnet);
wolffd@0 168 evidence = cell(ss,T);
wolffd@0 169 % we add dummy data to O to force its effective size to be 1.
wolffd@0 170 % The actual values have already been incorporated into softCPDpot
wolffd@0 171 evidence(O,:) = num2cell(ones(1,T));
wolffd@0 172 [engine, ll_dbn] = enter_evidence(engine, evidence, 'softCPDpot', softCPDpot);
wolffd@0 173
wolffd@0 174
wolffd@0 175 %evidence(F,:) = num2cell(2*ones(1,T));
wolffd@0 176 %[engine, ll_dbn] = enter_evidence(engine, evidence);
wolffd@0 177
wolffd@0 178
wolffd@0 179 gamma = zeros(nwords, T);
wolffd@0 180 for t=1:T
wolffd@0 181 m = marginal_nodes(engine, [W F], t);
wolffd@0 182 gamma(:,t) = m.T(:,2);
wolffd@0 183 end
wolffd@0 184
wolffd@0 185 gamma
wolffd@0 186
wolffd@0 187 xidbn = zeros(nwords, nwords);
wolffd@0 188 for t=1:T-1
wolffd@0 189 m = marginal_nodes(engine, [W F W+ss], t);
wolffd@0 190 xidbn = xidbn + squeeze(m.T(:,2,:));
wolffd@0 191 end
wolffd@0 192
wolffd@0 193 % thee
wolffd@0 194 % xidbn(1,4) = 0.9412 the->e
wolffd@0 195 % (2,3)=0.0588 t->h
wolffd@0 196 % (3,4)=0.0588 h-e
wolffd@0 197 % (4,4)=0.0588 e-e
wolffd@0 198
wolffd@0 199
wolffd@0 200 end