comparison toolboxes/FullBNT-1.0.7/bnt/CPDs/@softmax_CPD/update_ess.m @ 0:e9a9cd732c1e tip

first hg version after svn

+function CPD = update_ess(CPD, fmarginal, evidence, ns, cnodes, hidden_bitv)
+% UPDATE_ESS Update the Expected Sufficient Statistics of a softmax node
+% function CPD = update_ess(CPD, fmarginal, evidence, ns, cnodes, hidden_bitv)
+%
+% fmarginal = overall posterior distribution of self and its parents
+% fmarginal(i1,i2...,ik,s)=prob(Pa1=i1,...,Pak=ik, self=s| X)
+% 
+% => 1) prob(self|Pa1,...,Pak)=fmarginal/prob(Pa1,...,Pak) with prob(Pa1,...,Pak)=sum{s,fmarginal}
+%       [self estimation -> CPD.self_vals]
+% 	  2) prob(Pa1,...,Pak) [WIRLS weights -> CPD.eso_weights]
+%
+% Hidden_bitv is ignored
+
+% Written by Pierpaolo Brutti
+
+if ~adjustable_CPD(CPD), return; end
+
+domain     = fmarginal.domain;                              
+self       = domain(end);          
+ps         = domain(1:end-1);                                     
+cnodes     = domain(CPD.cpndx);
+cps        = myintersect(domain, cnodes);                     
+dps        = mysetdiff(ps, cps);                            
+dn_use     = dps;
+if isempty(evidence{self}) dn_use = [dn_use self]; end % if self is hidden we must consider its dimension  
+dps_as_cps = domain(CPD.dps_as_cps.ndx);
+odom       = domain(~isemptycell(evidence(domain))); 
+
+ns = zeros(1, max(domain));
+ns(domain) = CPD.sizes;     % CPD.sizes = bnet.node_sizes([ps self]);
+ens = ns;                   % effective node sizes
+ens(odom) = 1;              
+dpsize = prod(ns(dps));
+
+% Extract the params compatible with the observations (if any) on the discrete parents (if any)
+dops = myintersect(dps, odom);
+dpvals = cat(1, evidence{dops});
+
+subs = ind2subv(ens(dn_use), 1:prod(ens(dn_use)));
+dpmap = find_equiv_posns(dops, dn_use);
+if ~isempty(dpmap), subs(:,dpmap) = subs(:,dpmap)+repmat(dpvals(:)',[size(subs,1) 1])-1; end
+supportedQs = subv2ind(ns(dn_use), subs); subs=subs(1:prod(ens(dps)),1:length(dps));
+Qarity = prod(ns(dn_use));
+if isempty(dn_use), Qarity = 1; end   
+
+fullm.T              = zeros(Qarity, 1);
+fullm.T(supportedQs) = fmarginal.T(:);
+rs_dim = CPD.sizes;    rs_dim(CPD.cpndx) = 1;           %
+if ~isempty(evidence{self}), rs_dim(end)=1; end         % reshaping the marginal
+fullm.T              = reshape(fullm.T, rs_dim);        %
+
+% --------------------------------------------------------------------------------UPDATE--
+
+CPD.nsamples = CPD.nsamples + 1;
+
+% 1) observations vector -> CPD.parents_vals ---------------------------------------------
+cpvals = cat(1, evidence{cps});
+
+if ~isempty(dps_as_cps),   % ...get in the dp_as_cp parents... 
+    separator          = CPD.dps_as_cps.separator;
+    dp_as_cpmap        = find_equiv_posns(dps_as_cps, dps);       
+    for i=1:dpsize,
+        dp_as_cpvals=zeros(1,sum(ns(dps_as_cps)));
+        possible_vals = ind2subv(ns(dps),i);
+        ll=find(ismember(subs(:,dp_as_cpmap), possible_vals(dp_as_cpmap), 'rows')==1);   
+        if ~isempty(ll),
+            where_one = separator + possible_vals(dp_as_cpmap);
+            dp_as_cpvals(where_one)=1;                            
+        end
+        CPD.parent_vals(CPD.nsamples,:,i) = [dp_as_cpvals(:); cpvals(:)]';
+    end
+else
+    CPD.parent_vals(CPD.nsamples,:) = cpvals(:)';
+end
+
+% 2) weights vector -> CPD.eso_weights ----------------------------------------------------
+if isempty(evidence{self}),             % self is hidden
+    pesi=reshape(sum(fullm.T, length(rs_dim)),[dpsize,1]);
+else
+    pesi=reshape(fullm.T,[dpsize,1]);
+end
+assert(approxeq(sum(pesi),1));          % check
+
+% 3) estimate (if R is hidden) or recover (if R is obs) self'value-------------------------
+if isempty(evidence{self})                                  % P(self|Pa1,...,Pak)=fmarginal/prob(Pa1,...,Pak)
+    r=reshape(mk_stochastic(fullm.T), [dpsize ns(self)]);   % matrix size: prod{j,ns(Paj)} x ns(self)      
+else
+    r = zeros(dpsize,ns(self));
+    for i=1:dpsize, if pesi(i)~=0, r(i,evidence{self}) = 1; end; end
+end
+for i=1:dpsize, if pesi(i)~=0, assert(approxeq(sum(r(i,:)),1)); end; end     % check
+
+% 4) save the previous values --------------------------------------------------------------
+for i=1:dpsize
+    CPD.eso_weights(CPD.nsamples,:,i)=pesi(i);
+    CPD.self_vals(CPD.nsamples,:,i) = r(i,:); 
+end