diff options
Diffstat (limited to 'sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD')
21 files changed, 1115 insertions, 0 deletions
diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/CVS/Entries b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/CVS/Entries new file mode 100644 index 00000000..0afc5821 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/CVS/Entries @@ -0,0 +1,7 @@ +/hhmmQ_CPD.m/1.1.1.1/Tue Sep 24 04:19:26 2002// +/log_prior.m/1.1.1.1/Wed May 29 15:59:54 2002// +/maximize_params.m/1.1.1.1/Tue Sep 24 13:10:18 2002// +/reset_ess.m/1.1.1.1/Wed May 29 15:59:54 2002// +/update_CPT.m/1.1.1.1/Tue Sep 24 02:58:18 2002// +/update_ess.m/1.1.1.1/Thu Jul 24 13:41:34 2003// +D/Old//// diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/CVS/Repository b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/CVS/Repository new file mode 100644 index 00000000..f226b7fc --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/CVS/Repository @@ -0,0 +1 @@ +FullBNT/BNT/CPDs/@hhmmQ_CPD diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/CVS/Root b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/CVS/Root new file mode 100644 index 00000000..f3bd14a6 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/CVS/Root @@ -0,0 +1 @@ +:ext:nsaunier@bnt.cvs.sourceforge.net:/cvsroot/bnt diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/CVS/Entries b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/CVS/Entries new file mode 100644 index 00000000..06bd5c8c --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/CVS/Entries @@ -0,0 +1,10 @@ +/hhmmQ_CPD.m/1.1.1.1/Mon Jun 24 18:19:00 2002// +/log_prior.m/1.1.1.1/Mon Jun 24 18:19:00 2002// +/maximize_params.m/1.1.1.1/Mon Jun 24 18:19:00 2002// +/reset_ess.m/1.1.1.1/Mon Jun 24 18:19:00 2002// +/update_CPT.m/1.1.1.1/Tue Sep 24 02:30:32 2002// +/update_ess.m/1.1.1.1/Mon Jun 24 18:19:00 2002// +/update_ess2.m/1.1.1.1/Mon Jun 24 21:20:52 2002// +/update_ess3.m/1.1.1.1/Mon Jun 24 22:08:08 2002// +/update_ess4.m/1.1.1.1/Mon Jun 24 22:23:32 2002// +D diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/CVS/Repository b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/CVS/Repository new file mode 100644 index 00000000..8e7c978a --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/CVS/Repository @@ -0,0 +1 @@ +FullBNT/BNT/CPDs/@hhmmQ_CPD/Old diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/CVS/Root b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/CVS/Root new file mode 100644 index 00000000..f3bd14a6 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/CVS/Root @@ -0,0 +1 @@ +:ext:nsaunier@bnt.cvs.sourceforge.net:/cvsroot/bnt diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/hhmmQ_CPD.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/hhmmQ_CPD.m new file mode 100644 index 00000000..24ef464b --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/hhmmQ_CPD.m @@ -0,0 +1,126 @@ +function CPD = hhmmQ_CPD(bnet, self, Qnodes, d, D, varargin) +% HHMMQ_CPD Make the CPD for a Q node at depth D of a D-level hierarchical HMM +% CPD = hhmmQ_CPD(bnet, self, Qnodes, d, D, ...) +% +% Fd(t-1) \ Q1:d-1(t) +% \ | +% \ v +% Qd(t-1) -> Qd(t) +% / +% / +% Fd+1(t-1) +% +% We assume parents are ordered (numbered) as follows: +% Qd(t-1), Fd+1(t-1), Fd(t-1), Q1(t), ..., Qd(t) +% +% The parents of Qd(t) can either be just Qd-1(t) or the whole stack Q1:d-1(t) (allQ) +% In either case, we will call them Qps. +% If d=1, Qps does not exist. Also, the F1(t-1) -> Q1(t) arc is optional. +% If the arc is missing, startprob does not need to be specified, +% since the toplevel is assumed to never reset (F1 does not exist). +% If d=D, Fd+1(t-1) does not exist (there is no signal from below). +% +% optional args [defaults] +% +% transprob - transprob(i,k,j) = prob transition from i to j given Qps = k ['leftright'] +% selfprob - prob of a transition from i to i given Qps=k [0.1] +% startprob - startprob(k,j) = prob start in j given Qps = k ['leftstart'] +% startargs - other args to be passed to the sub tabular_CPD for learning startprob +% transargs - other args will be passed to the sub tabular_CPD for learning transprob +% allQ - 1 means use all Q nodes above d as parents, 0 means just level d-1 [0] +% F1toQ1 - 1 means add F1(t-1) -> Q1(t) arc, 0 means level 1 never resets [0] +% +% For d=1, startprob(1,j) is only needed if F1toQ1=1 +% Also, transprob(i,j) can be used instead of transprob(i,1,j). +% +% hhmmQ_CPD is a subclass of tabular_CPD so we inherit inference methods like CPD_to_pot, etc. +% +% We create isolated tabular_CPDs with no F parents to learn transprob/startprob +% so we can avail of e.g., entropic or Dirichlet priors. +% In the future, we will be able to represent the transprob using a tree_CPD. +% +% For details, see "Linear-time inference in hierarchical HMMs", Murphy and Paskin, NIPS'01. + + +ss = bnet.nnodes_per_slice; +%assert(self == Qnodes(d)+ss); +ns = bnet.node_sizes(:); +CPD.Qsizes = ns(Qnodes); +CPD.d = d; +CPD.D = D; +allQ = 0; + +% find out which parents to use, to get right size +for i=1:2:length(varargin) + switch varargin{i}, + case 'allQ', allQ = varargin{i+1}; + end +end + +if d==1 + CPD.Qps = []; +else + if allQ + CPD.Qps = Qnodes(1:d-1); + else + CPD.Qps = Qnodes(d-1); + end +end + +Qsz = ns(self); +Qpsz = prod(ns(CPD.Qps)); + +% set default arguments +startprob = 'leftstart'; +transprob = 'leftright'; +startargs = {}; +transargs = {}; +CPD.F1toQ1 = 0; +selfprob = 0.1; + +for i=1:2:length(varargin) + switch varargin{i}, + case 'transprob', transprob = varargin{i+1}; + case 'selfprob', selfprob = varargin{i+1}; + case 'startprob', startprob = varargin{i+1}; + case 'startargs', startargs = varargin{i+1}; + case 'transargs', transargs = varargin{i+1}; + case 'F1toQ1', CPD.F1toQ1 = varargin{i+1}; + end +end + +Qps = CPD.Qps + ss; +old_self = self-ss; + +if strcmp(transprob, 'leftright') + LR = mk_leftright_transmat(Qsz, selfprob); + transprob = repmat(reshape(LR, [1 Qsz Qsz]), [Qpsz 1 1]); % transprob(k,i,j) + transprob = permute(transprob, [2 1 3]); % now transprob(i,k,j) +end +transargs{end+1} = 'CPT'; +transargs{end+1} = transprob; +CPD.sub_CPD_trans = mk_isolated_tabular_CPD([old_self Qps], ns([old_self Qps self]), transargs); +S = struct(CPD.sub_CPD_trans); +CPD.transprob = myreshape(S.CPT, [Qsz Qpsz Qsz]); + + +if strcmp(startprob, 'leftstart') + startprob = zeros(Qpsz, Qsz); + startprob(:,1) = 1; +end + +if (d==1) & ~CPD.F1toQ1 + CPD.sub_CPD_start = []; + CPD.startprob = []; +else + startargs{end+1} = 'CPT'; + startargs{end+1} = startprob; + CPD.sub_CPD_start = mk_isolated_tabular_CPD(Qps, ns([Qps self]), startargs); + S = struct(CPD.sub_CPD_start); + CPD.startprob = myreshape(S.CPT, [Qpsz Qsz]); +end + +CPD = class(CPD, 'hhmmQ_CPD', tabular_CPD(bnet, self)); + +CPD = update_CPT(CPD); + diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/log_prior.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/log_prior.m new file mode 100644 index 00000000..d44bec5e --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/log_prior.m @@ -0,0 +1,8 @@ +function L = log_prior(CPD) +% LOG_PRIOR Return log P(theta) for a hhmm CPD +% L = log_prior(CPD) + +L = log_prior(CPD.sub_CPD_trans); +if ~isempty(CPD.sub_CPD_start) + L = L + log_prior(CPD.sub_CPD_start); +end diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/maximize_params.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/maximize_params.m new file mode 100644 index 00000000..0e4632aa --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/maximize_params.m @@ -0,0 +1,40 @@ +function CPD = maximize_params(CPD, temp) +% MAXIMIZE_PARAMS Set the params of a hhmmQ node to their ML/MAP values. +% CPD = maximize_params(CPD, temperature) + +Qsz = CPD.Qsizes(CPD.d); +Qpsz = prod(CPD.Qsizes(CPD.Qps)); + +if ~isempty(CPD.sub_CPD_start) + CPD.sub_CPD_start = maximize_params(CPD.sub_CPD_start, temp); + S = struct(CPD.sub_CPD_start); + CPD.startprob = myreshape(S.CPT, [Qpsz Qsz]); + %CPD.startprob = S.CPT; +end + +if 1 + % If we are in a state that can only go the end state, + % we will never see a transition to another (non-end) state, + % so counts(i,k,j)=0 (and termprob(k,i)=1). + % We set counts(i,k,i)=1 in this case. + % This will cause remove_hhmm_end_state to return a + % stochastic matrix, but otherwise has no effect on EM. + counts = get_field(CPD.sub_CPD_trans, 'counts'); + counts = reshape(counts, [Qsz Qpsz Qsz]); + for k=1:Qpsz + for i=1:Qsz + if sum(counts(i,k,:))==0 % never witnessed a transition out of i + counts(i,k,i)=1; % add self loop + %fprintf('CPDQ d=%d i=%d k=%d\n', CPD.d, i, k); + end + end + end + CPD.sub_CPD_trans = set_fields(CPD.sub_CPD_trans, 'counts', counts(:)); +end + +CPD.sub_CPD_trans = maximize_params(CPD.sub_CPD_trans, temp); +S = struct(CPD.sub_CPD_trans); +%CPD.transprob = S.CPT; +CPD.transprob = myreshape(S.CPT, [Qsz Qpsz Qsz]); + +CPD = update_CPT(CPD); diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/reset_ess.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/reset_ess.m new file mode 100644 index 00000000..45a70ad7 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/reset_ess.m @@ -0,0 +1,8 @@ +function CPD = reset_ess(CPD) +% RESET_ESS Reset the Expected Sufficient Statistics of a hhmm Q node. +% CPD = reset_ess(CPD) + +if ~isempty(CPD.sub_CPD_start) + CPD.sub_CPD_start = reset_ess(CPD.sub_CPD_start); +end +CPD.sub_CPD_trans = reset_ess(CPD.sub_CPD_trans); diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_CPT.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_CPT.m new file mode 100644 index 00000000..503c225b --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_CPT.m @@ -0,0 +1,74 @@ +function CPD = update_CPT(CPD) +% Compute the big CPT for an HHMM Q node (including F parents) given internal transprob and startprob +% function CPD = update_CPT(CPD) + +Qsz = CPD.Qsz; +Qpsz = CPD.Qpsz; + +if ~isempty(CPD.Fbelow_ndx) + if ~isempty(CPD.Fself_ndx) % general case + % Fb(t-1) Fself(t-1) P(Q(t)=j| Q(t-1)=i, Qps(t)=k) + % ------------------------------------------------------ + % 1 1 delta(i,j) + % 2 1 transprob(i,k,j) + % 1 2 impossible + % 2 2 startprob(k,j) + CPT = zeros(Qsz, 2, 2, Qpsz, Qsz); + I = repmat(eye(Qsz), [1 1 Qpsz]); % i,j,k + I = permute(I, [1 3 2]); % i,k,j + CPT(:, 1, 1, :, :) = I; + CPT(:, 2, 1, :, :) = CPD.transprob; + CPT(:, 1, 2, :, :) = I; + CPT(:, 2, 2, :, :) = repmat(reshape(CPD.startprob, [1 Qpsz Qsz]), [Qsz 1 1]); % replicate over i + else % no F from self, hence no startprob + % Fb(t-1) P(Q(t)=j| Q(t-1)=i, Qps(t)=k) + % ------------------------------------------------------ + % 1 delta(i,j) + % 2 transprob(i,k,j) + + nps = length(CPD.dom_sz)-1; % num parents + CPT = 0*myones(CPD.dom_sz); + %CPT = zeros(Qsz, 2, Qpsz, Qsz); % assumes CPT(Q(t-1), F(t-1), Qps, Q(t)) + % but a member of Qps may preceed Q(t-1) or F(t-1) in the ordering + + I = repmat(eye(Qsz), [1 1 Qpsz]); % i,j,k + I = permute(I, [1 3 2]); % i,k,j + + % the following fails if there is a member of Qps with a lower + % number than F + %CPT(:, 1, :, :) = I; + %CPT(:, 2, :, :) = CPD.transprob; + + ndx = mk_multi_index(nps+1, CPD.Fbelow_ndx, 1); + CPT(ndx{:}) = I; + ndx = mk_multi_index(nps+1, CPD.Fbelow_ndx, 2); + CPT(ndx{:}) = CPD.transprob; + keyboard + end +else % no F signal from below + if ~isempty(CPD.Fself_ndx) + % Q(t-1), Fself(t-1), Qps, Q(t) + + % if condition start on previous concrete state (as in map learning), + % CPT(:, 1, :, :, :) = CPD.transprob(Q(t-1), Qps, Q(t)) + % CPT(:, 2, :, :, :) = CPD.startprob(Q(t-1), Qps, Q(t)) + + % Fself(t-1) P(Q(t-1)=i, Qps(t)=k -> Q(t)=j) + % ------------------------------------------------------ + % 1 transprob(i,k,j) + % 2 startprob(k,j) + CPT = zeros(Qsz, 2, Qpsz, Qsz); + I = repmat(eye(Qsz), [1 1 Qpsz]); % i,j,k + I = permute(I, [1 3 2]); % i,k,j + CPT(:, 1, :, :) = CPD.transprob; + if CPD.fullstartprob + CPT(:, 2, :, :) = CPD.startprob; + else + CPT(:, 2, :, :) = repmat(reshape(CPD.startprob, [1 Qpsz Qsz]), [Qsz 1 1]); % replicate over i + end + else % no F from self + error('An hhmmQ node without any F parents is just a tabular_CPD') + end +end + +CPD = set_fields(CPD, 'CPT', CPT); diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_ess.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_ess.m new file mode 100644 index 00000000..51c2bd1f --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_ess.m @@ -0,0 +1,141 @@ +function CPD = update_ess(CPD, fmarginal, evidence, ns, cnodes, hidden_bitv) +% UPDATE_ESS Update the Expected Sufficient Statistics of a hhmm Q node. +% function CPD = update_ess(CPD, fmarginal, evidence, ns, cnodes, idden_bitv) + +% Figure out the node numbers associated with each parent +% e.g., D=4, d=3, Qps = all Qs above, so dom = [Q3(t-1) F4(t-1) F3(t-1) Q1(t) Q2(t) Q3(t)]. +% so self = Q3(t), old_self = Q3(t-1), CPD.Qps = [1 2], Qps = [Q1(t) Q2(t)] +dom = fmarginal.domain; +self = dom(end); +old_self = dom(1); +Qps = dom(length(dom)-length(CPD.Qps):end-1); + +Qsz = CPD.Qsizes(CPD.d); +Qpsz = prod(CPD.Qsizes(CPD.Qps)); + +% If some of the Q nodes are observed (which happens during supervised training) +% the counts will only be non-zero in positions +% consistent with the evidence. We put the computed marginal responsibilities +% into the appropriate slots of the big counts array. +% (Recall that observed discrete nodes only have a single effective value.) +% (A more general, but much slower, way is to call add_evidence_to_dmarginal.) +% We assume the F nodes are never observed. + +obs_self = ~hidden_bitv(self); +obs_Qps = (~isempty(Qps)) & (~any(hidden_bitv(Qps))); % we assume that all or none of the Q parents are observed + +if obs_self + self_val = evidence{self}; + oldself_val = evidence{old_self}; +end + +if obs_Qps + Qps_val = subv2ind(Qpsz, cat(1, evidence{Qps})); + if Qps_val == 0 + keyboard + end +end + +if CPD.d==1 % no Qps from above + if ~CPD.F1toQ1 % no F from self + % marg(Q1(t-1), F2(t-1), Q1(t)) + % F2(t-1) P(Q1(t)=j | Q1(t-1)=i) + % 1 delta(i,j) + % 2 transprob(i,j) + if obs_self + hor_counts = zeros(Qsz, Qsz); + hor_counts(oldself_val, self_val) = fmarginal.T(2); + else + marg = reshape(fmarginal.T, [Qsz 2 Qsz]); + hor_counts = squeeze(marg(:,2,:)); + end + else + % marg(Q1(t-1), F2(t-1), F1(t-1), Q1(t)) + % F2(t-1) F1(t-1) P(Qd(t)=j| Qd(t-1)=i) + % ------------------------------------------------------ + % 1 1 delta(i,j) + % 2 1 transprob(i,j) + % 1 2 impossible + % 2 2 startprob(j) + if obs_self + marg = myreshape(fmarginal.T, [1 2 2 1]); + hor_counts = zeros(Qsz, Qsz); + hor_counts(oldself_val, self_val) = marg(1,2,1,1); + ver_counts = zeros(Qsz, 1); + %ver_counts(self_val) = marg(1,2,2,1); + ver_counts(self_val) = marg(1,2,2,1) + marg(1,1,2,1); + else + marg = reshape(fmarginal.T, [Qsz 2 2 Qsz]); + hor_counts = squeeze(marg(:,2,1,:)); + %ver_counts = squeeze(sum(marg(:,2,2,:),1)); % sum over i + ver_counts = squeeze(sum(marg(:,2,2,:),1)) + squeeze(sum(marg(:,1,2,:),1)); % sum i,b + end + end % F1toQ1 +else % d ~= 1 + if CPD.d < CPD.D % general case + % marg(Qd(t-1), Fd+1(t-1), Fd(t-1), Qps(t), Qd(t)) + % Fd+1(t-1) Fd(t-1) P(Qd(t)=j| Qd(t-1)=i, Qps(t)=k) + % ------------------------------------------------------ + % 1 1 delta(i,j) + % 2 1 transprob(i,k,j) + % 1 2 impossible + % 2 2 startprob(k,j) + if obs_Qps & obs_self + marg = myreshape(fmarginal.T, [1 2 2 1 1]); + k = 1; + hor_counts = zeros(Qsz, Qpsz, Qsz); + hor_counts(oldself_val, Qps_val, self_val) = marg(1, 2,1, k,1); + ver_counts = zeros(Qpsz, Qsz); + %ver_counts(Qps_val, self_val) = marg(1, 2,2, k,1); + ver_counts(Qps_val, self_val) = marg(1, 2,2, k,1) + marg(1, 1,2, k,1); + elseif obs_Qps & ~obs_self + marg = myreshape(fmarginal.T, [Qsz 2 2 1 Qsz]); + k = 1; + hor_counts = zeros(Qsz, Qpsz, Qsz); + hor_counts(:, Qps_val, :) = marg(:, 2,1, k,:); + ver_counts = zeros(Qpsz, Qsz); + %ver_counts(Qps_val, :) = sum(marg(:, 2,2, k,:), 1); + ver_counts(Qps_val, :) = sum(marg(:, 2,2, k,:), 1) + sum(marg(:, 1,2, k,:), 1); + elseif ~obs_Qps & obs_self + error('not yet implemented') + else % everything is hidden + marg = reshape(fmarginal.T, [Qsz 2 2 Qpsz Qsz]); + hor_counts = squeeze(marg(:,2,1,:,:)); % i,k,j + %ver_counts = squeeze(sum(marg(:,2,2,:,:),1)); % sum over i + ver_counts = squeeze(sum(marg(:,2,2,:,:),1)) + squeeze(sum(marg(:,1,2,:,:),1)); % sum over i,b + end + else % d == D, so no F from below + % marg(QD(t-1), FD(t-1), Qps(t), QD(t)) + % FD(t-1) P(QD(t)=j | QD(t-1)=i, Qps(t)=k) + % 1 transprob(i,k,j) + % 2 startprob(k,j) + if obs_Qps & obs_self + marg = myreshape(fmarginal.T, [1 2 1 1]); + k = 1; + hor_counts = zeros(Qsz, Qpsz, Qsz); + hor_counts(oldself_val, Qps_val, self_val) = marg(1, 1, k,1); + ver_counts = zeros(Qpsz, Qsz); + ver_counts(Qps_val, self_val) = marg(1, 2, k,1); + elseif obs_Qps & ~obs_self + marg = myreshape(fmarginal.T, [Qsz 2 1 Qsz]); + k = 1; + hor_counts = zeros(Qsz, Qpsz, Qsz); + hor_counts(:, Qps_val, :) = marg(:, 1, k,:); + ver_counts = zeros(Qpsz, Qsz); + ver_counts(Qps_val, :) = sum(marg(:, 2, k, :), 1); + elseif ~obs_Qps & obs_self + error('not yet implemented') + else % everything is hidden + marg = reshape(fmarginal.T, [Qsz 2 Qpsz Qsz]); + hor_counts = squeeze(marg(:,1,:,:)); + ver_counts = squeeze(sum(marg(:,2,:,:),1)); % sum over i + end + end +end + +CPD.sub_CPD_trans = update_ess_simple(CPD.sub_CPD_trans, hor_counts); + +if ~isempty(CPD.sub_CPD_start) + CPD.sub_CPD_start = update_ess_simple(CPD.sub_CPD_start, ver_counts); +end + diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_ess2.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_ess2.m new file mode 100644 index 00000000..41fc7380 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_ess2.m @@ -0,0 +1,178 @@ +function CPD = update_ess2(CPD, fmarginal, evidence, ns, cnodes, hidden_bitv) +% UPDATE_ESS Update the Expected Sufficient Statistics of a hhmm Q node. +% function CPD = update_ess(CPD, fmarginal, evidence, ns, cnodes, idden_bitv) + +% Figure out the node numbers associated with each parent +dom = fmarginal.domain; +self = dom(end); % by assumption +old_self = dom(CPD.old_self_ndx); +Fself = dom(CPD.Fself_ndx); +Fbelow = dom(CPD.Fbelow_ndx); +Qps = dom(CPD.Qps_ndx); + +Qsz = CPD.Qsz; +Qpsz = CPD.Qpsz; + + +fmarg = add_ev_to_dmarginal(fmarginal, evidence, ns); + + + +% hor_counts(old_self, Qps, self), +% fmarginal(old_self, Fbelow, Fself, Qps, self) +% hor_counts(i,k,j) = fmarginal(i,2,1,k,j) % below has finished, self has not +% ver_counts(i,k,j) = fmarginal(i,2,2,k,j) % below has finished, and so has self (reset) +% Since any of i,j,k may be observed, we write +% hor_counts(counts_ndx{:}) = fmarginal(fmarg_ndx{:}) +% where e.g., counts_ndx = {1, ':', 2} if Qps is hidden but we observe old_self=1, self=2. +% To create this counts_ndx, we write counts_ndx = mk_multi_ndx(3, obs_dim, obs_val) +% where counts_obs_dim = [1 3], counts_obs_val = [1 2] specifies the values of dimensions 1 and 3. + +counts_obs_dim = []; +fmarg_obs_dim = []; +obs_val = []; +if hidden_bitv(self) + effQsz = Qsz; +else + effQsz = 1; + counts_obs_dim = [counts_obs_dim 3]; + fmarg_obs_dim = [fmarg_obs_dim 5]; + obs_val = [obs_val evidence{self}]; +end + +% e.g., D=4, d=3, Qps = all Qs above, so dom = [Q3(t-1) F4(t-1) F3(t-1) Q1(t) Q2(t) Q3(t)]. +% so self = Q3(t), old_self = Q3(t-1), CPD.Qps = [1 2], Qps = [Q1(t) Q2(t)] +dom = fmarginal.domain; +self = dom(end); +old_self = dom(1); +Qps = dom(length(dom)-length(CPD.Qps):end-1); + +Qsz = CPD.Qsizes(CPD.d); +Qpsz = prod(CPD.Qsizes(CPD.Qps)); + +% If some of the Q nodes are observed (which happens during supervised training) +% the counts will only be non-zero in positions +% consistent with the evidence. We put the computed marginal responsibilities +% into the appropriate slots of the big counts array. +% (Recall that observed discrete nodes only have a single effective value.) +% (A more general, but much slower, way is to call add_evidence_to_dmarginal.) +% We assume the F nodes are never observed. + +obs_self = ~hidden_bitv(self); +obs_Qps = (~isempty(Qps)) & (~any(hidden_bitv(Qps))); % we assume that all or none of the Q parents are observed + +if obs_self + self_val = evidence{self}; + oldself_val = evidence{old_self}; +end + +if obs_Qps + Qps_val = subv2ind(Qpsz, cat(1, evidence{Qps})); + if Qps_val == 0 + keyboard + end +end + +if CPD.d==1 % no Qps from above + if ~CPD.F1toQ1 % no F from self + % marg(Q1(t-1), F2(t-1), Q1(t)) + % F2(t-1) P(Q1(t)=j | Q1(t-1)=i) + % 1 delta(i,j) + % 2 transprob(i,j) + if obs_self + hor_counts = zeros(Qsz, Qsz); + hor_counts(oldself_val, self_val) = fmarginal.T(2); + else + marg = reshape(fmarginal.T, [Qsz 2 Qsz]); + hor_counts = squeeze(marg(:,2,:)); + end + else + % marg(Q1(t-1), F2(t-1), F1(t-1), Q1(t)) + % F2(t-1) F1(t-1) P(Qd(t)=j| Qd(t-1)=i) + % ------------------------------------------------------ + % 1 1 delta(i,j) + % 2 1 transprob(i,j) + % 1 2 impossible + % 2 2 startprob(j) + if obs_self + marg = myreshape(fmarginal.T, [1 2 2 1]); + hor_counts = zeros(Qsz, Qsz); + hor_counts(oldself_val, self_val) = marg(1,2,1,1); + ver_counts = zeros(Qsz, 1); + %ver_counts(self_val) = marg(1,2,2,1); + ver_counts(self_val) = marg(1,2,2,1) + marg(1,1,2,1); + else + marg = reshape(fmarginal.T, [Qsz 2 2 Qsz]); + hor_counts = squeeze(marg(:,2,1,:)); + %ver_counts = squeeze(sum(marg(:,2,2,:),1)); % sum over i + ver_counts = squeeze(sum(marg(:,2,2,:),1)) + squeeze(sum(marg(:,1,2,:),1)); % sum i,b + end + end % F1toQ1 +else % d ~= 1 + if CPD.d < CPD.D % general case + % marg(Qd(t-1), Fd+1(t-1), Fd(t-1), Qps(t), Qd(t)) + % Fd+1(t-1) Fd(t-1) P(Qd(t)=j| Qd(t-1)=i, Qps(t)=k) + % ------------------------------------------------------ + % 1 1 delta(i,j) + % 2 1 transprob(i,k,j) + % 1 2 impossible + % 2 2 startprob(k,j) + if obs_Qps & obs_self + marg = myreshape(fmarginal.T, [1 2 2 1 1]); + k = 1; + hor_counts = zeros(Qsz, Qpsz, Qsz); + hor_counts(oldself_val, Qps_val, self_val) = marg(1, 2,1, k,1); + ver_counts = zeros(Qpsz, Qsz); + %ver_counts(Qps_val, self_val) = marg(1, 2,2, k,1); + ver_counts(Qps_val, self_val) = marg(1, 2,2, k,1) + marg(1, 1,2, k,1); + elseif obs_Qps & ~obs_self + marg = myreshape(fmarginal.T, [Qsz 2 2 1 Qsz]); + k = 1; + hor_counts = zeros(Qsz, Qpsz, Qsz); + hor_counts(:, Qps_val, :) = marg(:, 2,1, k,:); + ver_counts = zeros(Qpsz, Qsz); + %ver_counts(Qps_val, :) = sum(marg(:, 2,2, k,:), 1); + ver_counts(Qps_val, :) = sum(marg(:, 2,2, k,:), 1) + sum(marg(:, 1,2, k,:), 1); + elseif ~obs_Qps & obs_self + error('not yet implemented') + else % everything is hidden + marg = reshape(fmarginal.T, [Qsz 2 2 Qpsz Qsz]); + hor_counts = squeeze(marg(:,2,1,:,:)); % i,k,j + %ver_counts = squeeze(sum(marg(:,2,2,:,:),1)); % sum over i + ver_counts = squeeze(sum(marg(:,2,2,:,:),1)) + squeeze(sum(marg(:,1,2,:,:),1)); % sum over i,b + end + else % d == D, so no F from below + % marg(QD(t-1), FD(t-1), Qps(t), QD(t)) + % FD(t-1) P(QD(t)=j | QD(t-1)=i, Qps(t)=k) + % 1 transprob(i,k,j) + % 2 startprob(k,j) + if obs_Qps & obs_self + marg = myreshape(fmarginal.T, [1 2 1 1]); + k = 1; + hor_counts = zeros(Qsz, Qpsz, Qsz); + hor_counts(oldself_val, Qps_val, self_val) = marg(1, 1, k,1); + ver_counts = zeros(Qpsz, Qsz); + ver_counts(Qps_val, self_val) = marg(1, 2, k,1); + elseif obs_Qps & ~obs_self + marg = myreshape(fmarginal.T, [Qsz 2 1 Qsz]); + k = 1; + hor_counts = zeros(Qsz, Qpsz, Qsz); + hor_counts(:, Qps_val, :) = marg(:, 1, k,:); + ver_counts = zeros(Qpsz, Qsz); + ver_counts(Qps_val, :) = sum(marg(:, 2, k, :), 1); + elseif ~obs_Qps & obs_self + error('not yet implemented') + else % everything is hidden + marg = reshape(fmarginal.T, [Qsz 2 Qpsz Qsz]); + hor_counts = squeeze(marg(:,1,:,:)); + ver_counts = squeeze(sum(marg(:,2,:,:),1)); % sum over i + end + end +end + +CPD.sub_CPD_trans = update_ess_simple(CPD.sub_CPD_trans, hor_counts); + +if ~isempty(CPD.sub_CPD_start) + CPD.sub_CPD_start = update_ess_simple(CPD.sub_CPD_start, ver_counts); +end + diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_ess3.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_ess3.m new file mode 100644 index 00000000..da7ab6bd --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_ess3.m @@ -0,0 +1,80 @@ +function CPD = update_ess(CPD, fmarginal, evidence, ns, cnodes, hidden_bitv) +% UPDATE_ESS Update the Expected Sufficient Statistics of a hhmm Q node. +% function CPD = update_ess(CPD, fmarginal, evidence, ns, cnodes, idden_bitv) +% +% we assume if one of the Qps is observed, all of them are +% We assume the F nodes are already hidden + +% Figure out the node numbers associated with each parent +dom = fmarginal.domain; +self = dom(CPD.self_ndx); +old_self = dom(CPD.old_self_ndx); +%Fself = dom(CPD.Fself_ndx); +%Fbelow = dom(CPD.Fbelow_ndx); +Qps = dom(CPD.Qps_ndx); + +Qsz = CPD.Qsz; +Qpsz = CPD.Qpsz; + + +% hor_counts(old_self, Qps, self), +% fmarginal(old_self, Fbelow, Fself, Qps, self) +% hor_counts(i,k,j) = fmarginal(i,2,1,k,j) % below has finished, self has not +% ver_counts(i,k,j) = fmarginal(i,2,2,k,j) % below has finished, and so has self (reset) +% Since any of i,j,k may be observed, we write +% hor_counts(ndx{:}) = fmarginal(...) +% where e.g., ndx = {1, ':', 2} if Qps is hidden but we observe old_self=1, self=2. + +% ndx{i,k,j} +if hidden_bitv(old_self) + ndx{1} = ':'; +else + ndx{1} = evidence{old_self}; +end +if hidden_bitv(Qps) + ndx{2} = ':'; +else + ndx{2} = subv2ind(Qpsz, cat(1, evidence{Qps})); +end +if hidden_bitv(self) + ndx{3} = ':'; +else + ndx{3} = evidence{self}; +end + +fmarg = add_ev_to_dmarginal(fmarginal, evidence, ns); +% marg(Qold(t-1), Fbelow(t-1), Fself(t-1), Qps(t), Qself(t)) +hor_counts = zeros(Qsz, Qpsz, Qsz); +ver_counts = zeros(Qpsz, Qsz); + +if ~isempty(CPD.Fbelow_ndx) + if ~isempty(CPD.Fself_ndx) % general case + fmarg.T = myreshape(fmarg.T, [Qsz 2 2 Qpsz Qsz]); + marg_ndx = {ndx{1}, 2, 1, ndx{2}, ndx{3}}; + hor_counts(ndx{:}) = fmarg.T(marg_ndx{:}); + ver_counts(ndx{2:3}) = ... % sum over Fbelow and Qold=i + sum(fmarg.T({ndx{1}, 1, 2, ndx{2}, ndx{3}}),1) + .. + sum(fmarg.T({ndx{1}, 2, 2, ndx{2}, ndx{3}}),1); + else % no F from self, hence no startprob + fmarg.T = myreshape(fmarg.T, [Qsz 2 Qpsz Qsz]); + hor_counts(ndx{:}) = fmarg.T({ndx{1}, 2, ndx{2}, ndx{3}}); + end +else % no F signal from below + if ~isempty(CPD.Fself_ndx) % self F + fmarg.T = myreshape(fmarg.T, [Qsz 2 Qpsz Qsz]); + hor_counts(ndx{:}) = fmarg.T({ndx{1}, 1, ndx{2}, ndx{3}}); + ver_counts(ndx{2:3}) = ... % sum over Qold=i + sum(fmarg.T({ndx{1}, 2, ndx{2}, ndx{3}}),1); + else % no F from self + error('An hhmmQ node without any F parents is just a tabular_CPD') + end +end + + +CPD.sub_CPD_trans = update_ess_simple(CPD.sub_CPD_trans, hor_counts); + +if ~isempty(CPD.sub_CPD_start) + CPD.sub_CPD_start = update_ess_simple(CPD.sub_CPD_start, ver_counts); +end + + diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_ess4.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_ess4.m new file mode 100644 index 00000000..c826da1c --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/Old/update_ess4.m @@ -0,0 +1,95 @@ +function CPD = update_ess(CPD, fmarginal, evidence, ns, cnodes, hidden_bitv) +% UPDATE_ESS Update the Expected Sufficient Statistics of a hhmm Q node. +% function CPD = update_ess(CPD, fmarginal, evidence, ns, cnodes, idden_bitv) +% +% we assume if one of the Qps is observed, all of them are +% We assume the F nodes are already hidden + +% Figure out the node numbers associated with each parent +dom = fmarginal.domain; +self = dom(CPD.self_ndx); +old_self = dom(CPD.old_self_ndx); +%Fself = dom(CPD.Fself_ndx); +%Fbelow = dom(CPD.Fbelow_ndx); +Qps = dom(CPD.Qps_ndx); + +Qsz = CPD.Qsz; +Qpsz = CPD.Qpsz; + + +% hor_counts(old_self, Qps, self), +% fmarginal(old_self, Fbelow, Fself, Qps, self) +% hor_counts(i,k,j) = fmarginal(i,2,1,k,j) % below has finished, self has not +% ver_counts(i,k,j) = fmarginal(i,2,2,k,j) % below has finished, and so has self (reset) +% Since any of i,j,k may be observed, we write +% hor_counts(i_counts_ndx, kndx, jndx) = fmarginal(i_fmarg_ndx...) +% where i_fmarg_ndx = 1 and i_counts_ndx = i if old_self is observed to have value i, +% i_fmarg_ndx = 1:Qsz and i_counts_ndx = 1:Qsz if old_self is hidden, etc. + + +if hidden_bitv(old_self) + i_counts_ndx = 1:Qsz; + i_fmarg_ndx = 1:Qsz; + eff_oldQsz = Qsz; +else + i_counts_ndx = evidence{old_self}; + i_fmarg_ndx = 1; + eff_oldQsz = 1; +end + +if all(hidden_bitv(Qps)) % we assume all are hidden or all are observed + k_counts_ndx = 1:Qpsz; + k_fmarg_ndx = 1:Qpsz; + eff_Qpsz = Qpsz; +else + k_counts_ndx = subv2ind(Qpsz, cat(1, evidence{Qps})); + k_fmarg_ndx = 1; + eff_Qpsz = 1; +end + +if hidden_bitv(self) + j_counts_ndx = 1:Qsz; + j_fmarg_ndx = 1:Qsz; + eff_Qsz = Qsz; +else + j_counts_ndx = evidence{self}; + j_fmarg_ndx = 1; + eff_Qsz = 1; +end + +hor_counts = zeros(Qsz, Qpsz, Qsz); +ver_counts = zeros(Qpsz, Qsz); + +if ~isempty(CPD.Fbelow_ndx) + if ~isempty(CPD.Fself_ndx) % general case + fmarg.T = myreshape(fmarg.T, [eff_oldQsz 2 2 eff_Qpsz eff_Qsz]); + hor_counts(i_counts_ndx, k_counts_ndx, j_counts_ndx) = ... + fmarg.T(:, i_fmarg_ndx, 2, 1, k_fmarg_ndx, j_fmarg_ndx); + ver_counts(k_counts_ndx, j_counts_ndx) = ... % sum over Fbelow and Qold + sum(fmarg.T(:, 1, 2, k_fmarg_ndx, j_fmarg_ndx), 1) + ... + sum(fmarg.T(:, 2, 2, k_fmarg_ndx, j_fmarg_ndx), 1); + else % no F from self, hence no startprob + fmarg.T = myreshape(fmarg.T, [eff_oldQsz 2 eff_Qpsz eff_Qsz]); + hor_counts(i_counts_ndx, k_counts_ndx, j_counts_ndx) = ... + fmarg.T(i_fmarg_ndx, 2, k_fmarg_ndx, j_fmarg_ndx); + end +else % no F signal from below + if ~isempty(CPD.Fself_ndx) % self F + fmarg.T = myreshape(fmarg.T, [eff_oldQsz 2 eff_Qpsz eff_Qsz]); + hor_counts(i_counts_ndx, k_counts_ndx, j_counts_ndx) = ... + fmarg.T(i_fmarg_ndx, 1, k_fmarg_ndx, j_fmarg_ndx); + ver_counts(k_counts_ndx, j_counts_ndx) = ... % sum over Qold + sum(fmarg.T(:, 2, k_fmarg_ndx, j_fmarg_ndx), 1); + else % no F from self + error('An hhmmQ node without any F parents is just a tabular_CPD') + end +end + + +CPD.sub_CPD_trans = update_ess_simple(CPD.sub_CPD_trans, hor_counts); + +if ~isempty(CPD.sub_CPD_start) + CPD.sub_CPD_start = update_ess_simple(CPD.sub_CPD_start, ver_counts); +end + + diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/hhmmQ_CPD.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/hhmmQ_CPD.m new file mode 100644 index 00000000..6f289602 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/hhmmQ_CPD.m @@ -0,0 +1,132 @@ +function CPD = hhmmQ_CPD(bnet, self, varargin) +% HHMMQ_CPD Make the CPD for a Q node in a hierarchical HMM +% CPD = hhmmQ_CPD(bnet, self, ...) +% +% Fself(t-1) Qps(t) +% \ | +% \ v +% Qold(t-1) -> Q(t) +% / +% / +% Fbelow(t-1) +% +% Let ss = slice size = num. nodes per slice. +% This node is Q(t), and has mandatory parents Qold(t-1) (assumed to be numbered Q(t)-ss) +% and optional parents Fbelow, Fself, Qps. +% We require parents to be ordered (numbered) as follows: +% Qold, Fbelow, Fself, Qps, Q. +% +% If Fself=2, we use the transition matrix, else we use the prior matrix. +% If Fself node is omitted (eg. top level), we always use the transition matrix. +% If Fbelow=2, we may change state, otherwise we must stay in the same state. +% If Fbelow node is omitted (eg., bottom level), we may change state at every step. +% If Qps (Q parents) are specified, all parameters are conditioned on their joint value. +% We may choose any subset of nodes to condition on, as long as they as numbered lower than self. +% +% optional args [defaults] +% +% Fself - node number <= ss +% Fbelow - node number <= ss +% Qps - node numbers (all <= 2*ss) - uses 2TBN indexing +% transprob - transprob(i,k,j) = prob transition from i to j given Qps = k ['leftright'] +% selfprob - prob of a transition from i to i given Qps=k [0.1] +% startprob - startprob(k,j) = prob start in j given Qps = k ['leftstart'] +% startargs - other args to be passed to the sub tabular_CPD for learning startprob +% transargs - other args will be passed to the sub tabular_CPD for learning transprob +% fullstartprob - 1 means startprob depends on Q(t-1) [0] +% hhmmQ_CPD is a subclass of tabular_CPD so we inherit inference methods like CPD_to_pot, etc. +% +% We create isolated tabular_CPDs with no F parents to learn transprob/startprob +% so we can avail of e.g., entropic or Dirichlet priors. +% In the future, we will be able to represent the transprob using a tree_CPD. +% +% For details, see "Linear-time inference in hierarchical HMMs", Murphy and Paskin, NIPS'01. + + +ss = bnet.nnodes_per_slice; +ns = bnet.node_sizes(:); + +% set default arguments +Fself = []; +Fbelow = []; +Qps = []; +startprob = 'leftstart'; +transprob = 'leftright'; +startargs = {}; +transargs = {}; +selfprob = 0.1; +fullstartprob = 0; + +for i=1:2:length(varargin) + switch varargin{i}, + case 'Fself', Fself = varargin{i+1}; + case 'Fbelow', Fbelow = varargin{i+1}; + case 'Qps', Qps = varargin{i+1}; + case 'transprob', transprob = varargin{i+1}; + case 'selfprob', selfprob = varargin{i+1}; + case 'startprob', startprob = varargin{i+1}; + case 'startargs', startargs = varargin{i+1}; + case 'transargs', transargs = varargin{i+1}; + case 'fullstartprob', fullstartprob = varargin{i+1}; + end +end + +CPD.fullstartprob = fullstartprob; + +ps = parents(bnet.dag, self); +ndsz = ns(:)'; +CPD.dom_sz = [ndsz(ps) ns(self)]; +CPD.Fself_ndx = find_equiv_posns(Fself, ps); +CPD.Fbelow_ndx = find_equiv_posns(Fbelow, ps); +%CPD.Qps_ndx = find_equiv_posns(Qps+ss, ps); +CPD.Qps_ndx = find_equiv_posns(Qps, ps); +old_self = self-ss; +CPD.old_self_ndx = find_equiv_posns(old_self, ps); + +Qps = ps(CPD.Qps_ndx); +CPD.Qsz = ns(self); +CPD.Qpsz = prod(ns(Qps)); +CPD.Qpsizes = ns(Qps); +Qsz = CPD.Qsz; +Qpsz = CPD.Qpsz; + +if strcmp(transprob, 'leftright') + LR = mk_leftright_transmat(Qsz, selfprob); + transprob = repmat(reshape(LR, [1 Qsz Qsz]), [Qpsz 1 1]); % transprob(k,i,j) + transprob = permute(transprob, [2 1 3]); % now transprob(i,k,j) +end +transargs{end+1} = 'CPT'; +transargs{end+1} = transprob; +CPD.sub_CPD_trans = mk_isolated_tabular_CPD(ns([old_self Qps self]), transargs); +S = struct(CPD.sub_CPD_trans); +%CPD.transprob = myreshape(S.CPT, [Qsz Qpsz Qsz]); +CPD.transprob = S.CPT; + + +if strcmp(startprob, 'leftstart') + startprob = zeros(Qpsz, Qsz); + startprob(:,1) = 1; +end +if isempty(CPD.Fself_ndx) + CPD.sub_CPD_start = []; + CPD.startprob = []; +else + startargs{end+1} = 'CPT'; + startargs{end+1} = startprob; + if CPD.fullstartprob + CPD.sub_CPD_start = mk_isolated_tabular_CPD(ns([self Qps self]), startargs); + S = struct(CPD.sub_CPD_start); + %CPD.startprob = myreshape(S.CPT, [Qsz Qpsz Qsz]); + CPD.startprob = S.CPT; + else + CPD.sub_CPD_start = mk_isolated_tabular_CPD(ns([Qps self]), startargs); + S = struct(CPD.sub_CPD_start); + %CPD.startprob = myreshape(S.CPT, [CPD.Qpsizes Qsz]); + CPD.startprob = S.CPT; + end +end + +CPD = class(CPD, 'hhmmQ_CPD', tabular_CPD(bnet, self)); + +CPD = update_CPT(CPD); + diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/log_prior.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/log_prior.m new file mode 100644 index 00000000..d44bec5e --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/log_prior.m @@ -0,0 +1,8 @@ +function L = log_prior(CPD) +% LOG_PRIOR Return log P(theta) for a hhmm CPD +% L = log_prior(CPD) + +L = log_prior(CPD.sub_CPD_trans); +if ~isempty(CPD.sub_CPD_start) + L = L + log_prior(CPD.sub_CPD_start); +end diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/maximize_params.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/maximize_params.m new file mode 100644 index 00000000..541a50be --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/maximize_params.m @@ -0,0 +1,40 @@ +function CPD = maximize_params(CPD, temp) +% MAXIMIZE_PARAMS Set the params of a hhmmQ node to their ML/MAP values. +% CPD = maximize_params(CPD, temperature) + +Qsz = CPD.Qsz; +Qpsz = CPD.Qpsz; + +if ~isempty(CPD.sub_CPD_start) + CPD.sub_CPD_start = maximize_params(CPD.sub_CPD_start, temp); + S = struct(CPD.sub_CPD_start); + CPD.startprob = myreshape(S.CPT, [Qpsz Qsz]); + %CPD.startprob = S.CPT; +end + +if 1 + % If we are in a state that can only go the end state, + % we will never see a transition to another (non-end) state, + % so counts(i,k,j)=0 (and termprob(k,i)=1). + % We set counts(i,k,i)=1 in this case. + % This will cause remove_hhmm_end_state to return a + % stochastic matrix, but otherwise has no effect on EM. + counts = get_field(CPD.sub_CPD_trans, 'counts'); + counts = reshape(counts, [Qsz Qpsz Qsz]); + for k=1:Qpsz + for i=1:Qsz + if sum(counts(i,k,:))==0 % never witnessed a transition out of i + counts(i,k,i)=1; % add self loop + %fprintf('CPDQ d=%d i=%d k=%d\n', CPD.d, i, k); + end + end + end + CPD.sub_CPD_trans = set_fields(CPD.sub_CPD_trans, 'counts', counts(:)); +end + +CPD.sub_CPD_trans = maximize_params(CPD.sub_CPD_trans, temp); +S = struct(CPD.sub_CPD_trans); +%CPD.transprob = S.CPT; +CPD.transprob = myreshape(S.CPT, [Qsz Qpsz Qsz]); + +CPD = update_CPT(CPD); diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/reset_ess.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/reset_ess.m new file mode 100644 index 00000000..45a70ad7 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/reset_ess.m @@ -0,0 +1,8 @@ +function CPD = reset_ess(CPD) +% RESET_ESS Reset the Expected Sufficient Statistics of a hhmm Q node. +% CPD = reset_ess(CPD) + +if ~isempty(CPD.sub_CPD_start) + CPD.sub_CPD_start = reset_ess(CPD.sub_CPD_start); +end +CPD.sub_CPD_trans = reset_ess(CPD.sub_CPD_trans); diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/update_CPT.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/update_CPT.m new file mode 100644 index 00000000..9ed1a352 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/update_CPT.m @@ -0,0 +1,70 @@ +function CPD = update_CPT(CPD) +% Compute the big CPT for an HHMM Q node (including F parents) given internal transprob and startprob +% function CPD = update_CPT(CPD) + +Qsz = CPD.Qsz; +Qpsz = CPD.Qpsz; + +if ~isempty(CPD.Fbelow_ndx) + if ~isempty(CPD.Fself_ndx) % general case + % Fb(t-1) Fself(t-1) P(Q(t)=j| Q(t-1)=i, Qps(t)=k) + % ------------------------------------------------------ + % 1 1 delta(i,j) + % 2 1 transprob(i,k,j) + % 1 2 impossible + % 2 2 startprob(k,j) + CPT = zeros(Qsz, 2, 2, Qpsz, Qsz); + I = repmat(eye(Qsz), [1 1 Qpsz]); % i,j,k + I = permute(I, [1 3 2]); % i,k,j + CPT(:, 1, 1, :, :) = I; + CPT(:, 2, 1, :, :) = CPD.transprob; + CPT(:, 1, 2, :, :) = I; + CPT(:, 2, 2, :, :) = repmat(reshape(CPD.startprob, [1 Qpsz Qsz]), ... + [Qsz 1 1]); % replicate over i + else % no F from self, hence no startprob + % Fb(t-1) P(Q(t)=j| Q(t-1)=i, Qps(t)=k) + % ------------------------------------------------------ + % 1 delta(i,j) + % 2 transprob(i,k,j) + + nps = length(CPD.dom_sz)-1; % num parents + CPT = 0*myones(CPD.dom_sz); + %CPT = zeros(Qsz, 2, Qpsz, Qsz); % assumes CPT(Q(t-1), F(t-1), Qps, Q(t)) + % but a member of Qps may preceed Q(t-1) or F(t-1) in the ordering + + for k=1:CPD.Qpsz + Qps_vals = ind2subv(CPD.Qpsizes, k); + ndx = mk_multi_index(nps+1, [CPD.Fbelow_ndx CPD.Qps_ndx], [1 Qps_vals]); + CPT(ndx{:}) = eye(Qsz); % CPT(:,2,k,:) or CPT(:,k,2,:) etc + end + ndx = mk_multi_index(nps+1, CPD.Fbelow_ndx, 2); + CPT(ndx{:}) = CPD.transprob; % we assume transprob is in topo order + end +else % no F signal from below + if ~isempty(CPD.Fself_ndx) + % Q(t-1), Fself(t-1), Qps, Q(t) + + % Fself(t-1) P(Q(t-1)=i, Qps(t)=k -> Q(t)=j) + % ------------------------------------------------------ + % 1 transprob(i,k,j) + % 2 startprob(k,j) + + nps = length(CPD.dom_sz)-1; % num parents + CPT = 0*myones(CPD.dom_sz); + ndx = mk_multi_index(nps+1, CPD.Fself_ndx, 1); + CPT(ndx{:}) = CPD.transprob; + if CPD.fullstartprob + ndx = mk_multi_index(nps+1, CPD.Fself_ndx, 2); + CPT(ndx{:}) = CPD.startprob; + else + for i=1:CPD.Qsz + ndx = mk_multi_index(nps+1, [CPD.Fself_ndx CPD.old_self_ndx], [2 i]); + CPT(ndx{:}) = CPD.startprob; + end + end + else % no F from self + error('An hhmmQ node without any F parents is just a tabular_CPD') + end +end + +CPD = set_fields(CPD, 'CPT', CPT); diff --git a/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/update_ess.m b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/update_ess.m new file mode 100644 index 00000000..07dfc72e --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@hhmmQ_CPD/update_ess.m @@ -0,0 +1,86 @@ +function CPD = update_ess(CPD, fmarginal, evidence, ns, cnodes, hidden_bitv) +% UPDATE_ESS Update the Expected Sufficient Statistics of a hhmm Q node. +% function CPD = update_ess(CPD, fmarginal, evidence, ns, cnodes, idden_bitv) +% +% we assume if one of the Qps is observed, all of them are +% We assume the F nodes are already hidden + +% Figure out the node numbers associated with each parent +dom = fmarginal.domain; +self = dom(end); +old_self = dom(CPD.old_self_ndx); +%Fself = dom(CPD.Fself_ndx); +%Fbelow = dom(CPD.Fbelow_ndx); +Qps = dom(CPD.Qps_ndx); + +Qsz = CPD.Qsz; +Qpsz = CPD.Qpsz; + + +% hor_counts(old_self, Qps, self), +% fmarginal(old_self, Fbelow, Fself, Qps, self) +% hor_counts(i,k,j) = fmarginal(i,2,1,k,j) % below has finished, self has not +% ver_counts(i,k,j) = fmarginal(i,2,2,k,j) % below has finished, and so has self (reset) +% Since any of i,j,k may be observed, we write +% hor_counts(i_counts_ndx, kndx, jndx) = fmarginal(i_fmarg_ndx...) +% where i_fmarg_ndx = 1 and i_counts_ndx = i if old_self is observed to have value i, +% i_fmarg_ndx = 1:Qsz and i_counts_ndx = 1:Qsz if old_self is hidden, etc. + + +if hidden_bitv(old_self) + i_counts_ndx = 1:Qsz; + eff_oldQsz = Qsz; +else + i_counts_ndx = evidence{old_self}; + eff_oldQsz = 1; +end + +if all(hidden_bitv(Qps)) % we assume all are hidden or all are observed + k_counts_ndx = 1:Qpsz; + eff_Qpsz = Qpsz; +else + k_counts_ndx = subv2ind(Qpsz, cat(1, evidence{Qps})); + eff_Qpsz = 1; +end + +if hidden_bitv(self) + j_counts_ndx = 1:Qsz; + eff_Qsz = Qsz; +else + j_counts_ndx = evidence{self}; + eff_Qsz = 1; +end + +hor_counts = zeros(Qsz, Qpsz, Qsz); +ver_counts = zeros(Qpsz, Qsz); + +if ~isempty(CPD.Fbelow_ndx) + if ~isempty(CPD.Fself_ndx) % general case + fmarg = myreshape(fmarginal.T, [eff_oldQsz 2 2 eff_Qpsz eff_Qsz]); + hor_counts(i_counts_ndx, k_counts_ndx, j_counts_ndx) = fmarg(:, 2, 1, :, :); + ver_counts(k_counts_ndx, j_counts_ndx) = ... % sum over Fbelow and Qold + sumv(fmarg(:, :, 2, :, :), [1 2]); % require Fself=2 + else % no F from self, hence no startprob + fmarg = myreshape(fmarginal.T, [eff_oldQsz 2 eff_Qpsz eff_Qsz]); + hor_counts(i_counts_ndx, k_counts_ndx, j_counts_ndx) = ... + fmarg(:, 2, :, :); % require Fbelow = 2 + end +else % no F signal from below + if ~isempty(CPD.Fself_ndx) % self F + fmarg = myreshape(fmarginal.T, [eff_oldQsz 2 eff_Qpsz eff_Qsz]); + hor_counts(i_counts_ndx, k_counts_ndx, j_counts_ndx) = fmarg(:, 1, :, :); + ver_counts(k_counts_ndx, j_counts_ndx) = ... % sum over Qold + squeeze(sum(fmarg(:, 2, :, :), 1)); % Fself=2 + else % no F from self + error('An hhmmQ node without any F parents is just a tabular_CPD') + end +end + + +CPD.sub_CPD_trans = update_ess_simple(CPD.sub_CPD_trans, hor_counts); + +if ~isempty(CPD.sub_CPD_start) + CPD.sub_CPD_start = update_ess_simple(CPD.sub_CPD_start, ver_counts); +end + + |