|
wolffd@0
|
1 function CPD = update_ess(CPD, fmarginal, evidence, ns, cnodes, hidden_bitv)
|
|
wolffd@0
|
2 % UPDATE_ESS Update the Expected Sufficient Statistics of a softmax node
|
|
wolffd@0
|
3 % function CPD = update_ess(CPD, fmarginal, evidence, ns, cnodes, hidden_bitv)
|
|
wolffd@0
|
4 %
|
|
wolffd@0
|
5 % fmarginal = overall posterior distribution of self and its parents
|
|
wolffd@0
|
6 % fmarginal(i1,i2...,ik,s)=prob(Pa1=i1,...,Pak=ik, self=s| X)
|
|
wolffd@0
|
7 %
|
|
wolffd@0
|
8 % => 1) prob(self|Pa1,...,Pak)=fmarginal/prob(Pa1,...,Pak) with prob(Pa1,...,Pak)=sum{s,fmarginal}
|
|
wolffd@0
|
9 % [self estimation -> CPD.self_vals]
|
|
wolffd@0
|
10 % 2) prob(Pa1,...,Pak) [WIRLS weights -> CPD.eso_weights]
|
|
wolffd@0
|
11 %
|
|
wolffd@0
|
12 % Hidden_bitv is ignored
|
|
wolffd@0
|
13
|
|
wolffd@0
|
14 % Written by Pierpaolo Brutti
|
|
wolffd@0
|
15
|
|
wolffd@0
|
16 if ~adjustable_CPD(CPD), return; end
|
|
wolffd@0
|
17
|
|
wolffd@0
|
18 domain = fmarginal.domain;
|
|
wolffd@0
|
19 self = domain(end);
|
|
wolffd@0
|
20 ps = domain(1:end-1);
|
|
wolffd@0
|
21 cnodes = domain(CPD.cpndx);
|
|
wolffd@0
|
22 cps = myintersect(domain, cnodes);
|
|
wolffd@0
|
23 dps = mysetdiff(ps, cps);
|
|
wolffd@0
|
24 dn_use = dps;
|
|
wolffd@0
|
25 if isempty(evidence{self}) dn_use = [dn_use self]; end % if self is hidden we must consider its dimension
|
|
wolffd@0
|
26 dps_as_cps = domain(CPD.dps_as_cps.ndx);
|
|
wolffd@0
|
27 odom = domain(~isemptycell(evidence(domain)));
|
|
wolffd@0
|
28
|
|
wolffd@0
|
29 ns = zeros(1, max(domain));
|
|
wolffd@0
|
30 ns(domain) = CPD.sizes; % CPD.sizes = bnet.node_sizes([ps self]);
|
|
wolffd@0
|
31 ens = ns; % effective node sizes
|
|
wolffd@0
|
32 ens(odom) = 1;
|
|
wolffd@0
|
33 dpsize = prod(ns(dps));
|
|
wolffd@0
|
34
|
|
wolffd@0
|
35 % Extract the params compatible with the observations (if any) on the discrete parents (if any)
|
|
wolffd@0
|
36 dops = myintersect(dps, odom);
|
|
wolffd@0
|
37 dpvals = cat(1, evidence{dops});
|
|
wolffd@0
|
38
|
|
wolffd@0
|
39 subs = ind2subv(ens(dn_use), 1:prod(ens(dn_use)));
|
|
wolffd@0
|
40 dpmap = find_equiv_posns(dops, dn_use);
|
|
wolffd@0
|
41 if ~isempty(dpmap), subs(:,dpmap) = subs(:,dpmap)+repmat(dpvals(:)',[size(subs,1) 1])-1; end
|
|
wolffd@0
|
42 supportedQs = subv2ind(ns(dn_use), subs); subs=subs(1:prod(ens(dps)),1:length(dps));
|
|
wolffd@0
|
43 Qarity = prod(ns(dn_use));
|
|
wolffd@0
|
44 if isempty(dn_use), Qarity = 1; end
|
|
wolffd@0
|
45
|
|
wolffd@0
|
46 fullm.T = zeros(Qarity, 1);
|
|
wolffd@0
|
47 fullm.T(supportedQs) = fmarginal.T(:);
|
|
wolffd@0
|
48 rs_dim = CPD.sizes; rs_dim(CPD.cpndx) = 1; %
|
|
wolffd@0
|
49 if ~isempty(evidence{self}), rs_dim(end)=1; end % reshaping the marginal
|
|
wolffd@0
|
50 fullm.T = reshape(fullm.T, rs_dim); %
|
|
wolffd@0
|
51
|
|
wolffd@0
|
52 % --------------------------------------------------------------------------------UPDATE--
|
|
wolffd@0
|
53
|
|
wolffd@0
|
54 CPD.nsamples = CPD.nsamples + 1;
|
|
wolffd@0
|
55
|
|
wolffd@0
|
56 % 1) observations vector -> CPD.parents_vals ---------------------------------------------
|
|
wolffd@0
|
57 cpvals = cat(1, evidence{cps});
|
|
wolffd@0
|
58
|
|
wolffd@0
|
59 if ~isempty(dps_as_cps), % ...get in the dp_as_cp parents...
|
|
wolffd@0
|
60 separator = CPD.dps_as_cps.separator;
|
|
wolffd@0
|
61 dp_as_cpmap = find_equiv_posns(dps_as_cps, dps);
|
|
wolffd@0
|
62 for i=1:dpsize,
|
|
wolffd@0
|
63 dp_as_cpvals=zeros(1,sum(ns(dps_as_cps)));
|
|
wolffd@0
|
64 possible_vals = ind2subv(ns(dps),i);
|
|
wolffd@0
|
65 ll=find(ismember(subs(:,dp_as_cpmap), possible_vals(dp_as_cpmap), 'rows')==1);
|
|
wolffd@0
|
66 if ~isempty(ll),
|
|
wolffd@0
|
67 where_one = separator + possible_vals(dp_as_cpmap);
|
|
wolffd@0
|
68 dp_as_cpvals(where_one)=1;
|
|
wolffd@0
|
69 end
|
|
wolffd@0
|
70 CPD.parent_vals(CPD.nsamples,:,i) = [dp_as_cpvals(:); cpvals(:)]';
|
|
wolffd@0
|
71 end
|
|
wolffd@0
|
72 else
|
|
wolffd@0
|
73 CPD.parent_vals(CPD.nsamples,:) = cpvals(:)';
|
|
wolffd@0
|
74 end
|
|
wolffd@0
|
75
|
|
wolffd@0
|
76 % 2) weights vector -> CPD.eso_weights ----------------------------------------------------
|
|
wolffd@0
|
77 if isempty(evidence{self}), % self is hidden
|
|
wolffd@0
|
78 pesi=reshape(sum(fullm.T, length(rs_dim)),[dpsize,1]);
|
|
wolffd@0
|
79 else
|
|
wolffd@0
|
80 pesi=reshape(fullm.T,[dpsize,1]);
|
|
wolffd@0
|
81 end
|
|
wolffd@0
|
82 assert(approxeq(sum(pesi),1)); % check
|
|
wolffd@0
|
83
|
|
wolffd@0
|
84 % 3) estimate (if R is hidden) or recover (if R is obs) self'value-------------------------
|
|
wolffd@0
|
85 if isempty(evidence{self}) % P(self|Pa1,...,Pak)=fmarginal/prob(Pa1,...,Pak)
|
|
wolffd@0
|
86 r=reshape(mk_stochastic(fullm.T), [dpsize ns(self)]); % matrix size: prod{j,ns(Paj)} x ns(self)
|
|
wolffd@0
|
87 else
|
|
wolffd@0
|
88 r = zeros(dpsize,ns(self));
|
|
wolffd@0
|
89 for i=1:dpsize, if pesi(i)~=0, r(i,evidence{self}) = 1; end; end
|
|
wolffd@0
|
90 end
|
|
wolffd@0
|
91 for i=1:dpsize, if pesi(i)~=0, assert(approxeq(sum(r(i,:)),1)); end; end % check
|
|
wolffd@0
|
92
|
|
wolffd@0
|
93 % 4) save the previous values --------------------------------------------------------------
|
|
wolffd@0
|
94 for i=1:dpsize
|
|
wolffd@0
|
95 CPD.eso_weights(CPD.nsamples,:,i)=pesi(i);
|
|
wolffd@0
|
96 CPD.self_vals(CPD.nsamples,:,i) = r(i,:);
|
|
wolffd@0
|
97 end
|
