|
wolffd@0
|
1 function obslik = mk_hmm_obs_lik_vec(engine, evidence)
|
|
wolffd@0
|
2
|
|
wolffd@0
|
3 % P(o1,o2| h) = P(o1|h) * P(o2|h) where h = Q1,Q2,...
|
|
wolffd@0
|
4
|
|
wolffd@0
|
5 bnet = bnet_from_engine(engine);
|
|
wolffd@0
|
6 ss = length(bnet.intra);
|
|
wolffd@0
|
7 onodes = bnet.observed;
|
|
wolffd@0
|
8 hnodes = mysetdiff(1:ss, onodes);
|
|
wolffd@0
|
9 ns = bnet.node_sizes(:);
|
|
wolffd@0
|
10 ns(onodes) = 1;
|
|
wolffd@0
|
11
|
|
wolffd@0
|
12 Q = length(engine.startprob);
|
|
wolffd@0
|
13 obslik = ones(Q, 1);
|
|
wolffd@0
|
14
|
|
wolffd@0
|
15 for i=1:length(onodes)
|
|
wolffd@0
|
16 o = onodes(i);
|
|
wolffd@0
|
17 %data = cell2num(evidence(o,1));
|
|
wolffd@0
|
18 data = evidence{o,1};
|
|
wolffd@0
|
19 if myismember(o, bnet.dnodes)
|
|
wolffd@0
|
20 obslik_i = eval_pdf_cond_multinomial(data, engine.obsprob{i}.CPT);
|
|
wolffd@0
|
21 else
|
|
wolffd@0
|
22 if bnet.auto_regressive(o)
|
|
wolffd@0
|
23 error('can''t handle AR nodes')
|
|
wolffd@0
|
24 end
|
|
wolffd@0
|
25 %% calling mk_ghmm_obs_lik, which calls gaussian_prob, is slow, so we inline it
|
|
wolffd@0
|
26 %% and use the pre-computed inverse matrix
|
|
wolffd@0
|
27 %obslik_i = mk_ghmm_obs_lik(data, engine.obsprob{i}.mu, engine.obsprob{i}.Sigma);
|
|
wolffd@0
|
28 x = data(:);
|
|
wolffd@0
|
29 m = engine.obsprob{i}.mu;
|
|
wolffd@0
|
30 Qi = size(m, 2);
|
|
wolffd@0
|
31 obslik_i = size(Qi, 1);
|
|
wolffd@0
|
32 invC = engine.obsprob{i}.inv_Sigma;
|
|
wolffd@0
|
33 denom = engine.obsprob{i}.denom;
|
|
wolffd@0
|
34 for j=1:Qi
|
|
wolffd@0
|
35 numer = exp(-0.5 * (x-m(:,j))' * invC(:,:,j) * (x-m(:,j)));
|
|
wolffd@0
|
36 obslik_i(j) = numer / denom(j);
|
|
wolffd@0
|
37 end
|
|
wolffd@0
|
38 end
|
|
wolffd@0
|
39 % convert P(o|ps) into P(o|h) by multiplying onto a (h,o) potential of all 1s
|
|
wolffd@0
|
40 ps = bnet.parents{o};
|
|
wolffd@0
|
41 dom = [ps o];
|
|
wolffd@0
|
42 obspot_i = dpot(dom, ns(dom), obslik_i);
|
|
wolffd@0
|
43 dom = [hnodes o];
|
|
wolffd@0
|
44 obspot = dpot(dom, ns(dom));
|
|
wolffd@0
|
45 obspot = multiply_by_pot(obspot, obspot_i);
|
|
wolffd@0
|
46 % compute p(oi|h) * p(oj|h)
|
|
wolffd@0
|
47 S = struct(obspot);
|
|
wolffd@0
|
48 obslik = obslik .* S.T(:);
|
|
wolffd@0
|
49 end
|
|
wolffd@0
|
50
|
|
wolffd@0
|
51
|
|
wolffd@0
|
52
|
