diff options
Diffstat (limited to 'sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD')
29 files changed, 652 insertions, 0 deletions
diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CPD_to_lambda_msg.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CPD_to_lambda_msg.m new file mode 100644 index 00000000..d53e9e0f --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CPD_to_lambda_msg.m @@ -0,0 +1,16 @@ +function lam_msg = CPD_to_lambda_msg(CPD, msg_type, n, ps, msg, p, evidence) +% CPD_TO_LAMBDA_MSG Compute lambda message (discrete) +% lam_msg = compute_lambda_msg(CPD, msg_type, n, ps, msg, p, evidence) +% Pearl p183 eq 4.52 + +switch msg_type + case 'd', + T = prod_CPT_and_pi_msgs(CPD, n, ps, msg, p); + mysize = length(msg{n}.lambda); + lambda = dpot(n, mysize, msg{n}.lambda); + T = multiply_by_pot(T, lambda); + lam_msg = pot_to_marginal(marginalize_pot(T, p)); + lam_msg = lam_msg.T; + case 'g', + error('discrete_CPD can''t create Gaussian msgs') +end diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CPD_to_pi.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CPD_to_pi.m new file mode 100644 index 00000000..5962c92e --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CPD_to_pi.m @@ -0,0 +1,13 @@ +function pi = CPD_to_pi(CPD, msg_type, n, ps, msg, evidence) +% COMPUTE_PI Compute pi vector (discrete) +% pi = compute_pi(CPD, msg_type, n, ps, msg, evidence) +% Pearl p183 eq 4.51 + +switch msg_type + case 'd', + T = prod_CPT_and_pi_msgs(CPD, n, ps, msg); + pi = pot_to_marginal(marginalize_pot(T, n)); + pi = pi.T(:); + case 'g', + error('can only convert discrete CPD to Gaussian pi if observed') +end diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CPD_to_scgpot.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CPD_to_scgpot.m new file mode 100644 index 00000000..3d611536 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CPD_to_scgpot.m @@ -0,0 +1,25 @@ +function pot = CPD_to_scgpot(CPD, domain, ns, cnodes, evidence) +% CPD_TO_SCGPOT Convert a CPD to a CG potential, incorporating any evidence (discrete) +% pot = CPD_to_scgpot(CPD, domain, ns, cnodes, evidence) +% +% domain is the domain of CPD. +% node_sizes(i) is the size of node i. +% cnodes +% evidence{i} is the evidence on the i'th node. + +%odom = domain(~isemptycell(evidence(domain))); + +%vals = cat(1, evidence{odom}); +%map = find_equiv_posns(odom, domain); +%index = mk_multi_index(length(domain), map, vals); +CPT = CPD_to_CPT(CPD); +%CPT = CPT(index{:}); +CPT = CPT(:); +%ns(odom) = 1; +potarray = cell(1, length(CPT)); +for i=1:length(CPT) + %p = CPT(i); + potarray{i} = scgcpot(0, 0, CPT(i)); + %scpot{i} = scpot(0, 0); +end +pot = scgpot(domain, [], [], ns, potarray); diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CVS/Entries b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CVS/Entries new file mode 100644 index 00000000..57599441 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CVS/Entries @@ -0,0 +1,15 @@ +/CPD_to_lambda_msg.m/1.1.1.1/Wed May 29 15:59:52 2002// +/CPD_to_pi.m/1.1.1.1/Wed May 29 15:59:52 2002// +/CPD_to_scgpot.m/1.1.1.1/Wed May 29 15:59:52 2002// +/README/1.1.1.1/Wed May 29 15:59:52 2002// +/convert_CPD_to_table_hidden_ps.m/1.1.1.1/Wed May 29 15:59:52 2002// +/convert_obs_CPD_to_table.m/1.1.1.1/Wed May 29 15:59:52 2002// +/convert_to_pot.m/1.1.1.1/Fri Feb 20 22:00:38 2004// +/convert_to_sparse_table.c/1.1.1.1/Wed May 29 15:59:52 2002// +/convert_to_table.m/1.1.1.1/Wed May 29 15:59:52 2002// +/discrete_CPD.m/1.1.1.1/Wed May 29 15:59:52 2002// +/dom_sizes.m/1.1.1.1/Wed May 29 15:59:52 2002// +/log_prob_node.m/1.1.1.1/Wed May 29 15:59:52 2002// +/prob_node.m/1.1.1.1/Wed May 29 15:59:52 2002// +/sample_node.m/1.1.1.1/Wed May 29 15:59:52 2002// +D diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CVS/Entries.Log b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CVS/Entries.Log new file mode 100644 index 00000000..9c6f22e4 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CVS/Entries.Log @@ -0,0 +1,2 @@ +A D/Old//// +A D/private//// diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CVS/Repository b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CVS/Repository new file mode 100644 index 00000000..f3418ec7 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CVS/Repository @@ -0,0 +1 @@ +FullBNT/BNT/CPDs/@discrete_CPD diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CVS/Root b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CVS/Root new file mode 100644 index 00000000..f3bd14a6 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/CVS/Root @@ -0,0 +1 @@ +:ext:nsaunier@bnt.cvs.sourceforge.net:/cvsroot/bnt diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/CVS/Entries b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/CVS/Entries new file mode 100644 index 00000000..15bb91c3 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/CVS/Entries @@ -0,0 +1,5 @@ +/convert_to_pot.m/1.1.1.1/Wed May 29 15:59:52 2002// +/convert_to_table.m/1.1.1.1/Wed May 29 15:59:52 2002// +/prob_CPD.m/1.1.1.1/Wed May 29 15:59:52 2002// +/prob_node.m/1.1.1.1/Wed May 29 15:59:52 2002// +D diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/CVS/Repository b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/CVS/Repository new file mode 100644 index 00000000..df41b4fd --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/CVS/Repository @@ -0,0 +1 @@ +FullBNT/BNT/CPDs/@discrete_CPD/Old diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/CVS/Root b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/CVS/Root new file mode 100644 index 00000000..f3bd14a6 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/CVS/Root @@ -0,0 +1 @@ +:ext:nsaunier@bnt.cvs.sourceforge.net:/cvsroot/bnt diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/convert_to_pot.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/convert_to_pot.m new file mode 100644 index 00000000..3f178e1c --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/convert_to_pot.m @@ -0,0 +1,44 @@ +function pot = convert_to_pot(CPD, pot_type, domain, evidence) +% CONVERT_TO_POT Convert a tabular CPD to one or more potentials +% pots = convert_to_pot(CPD, pot_type, domain, evidence) +% +% pots{i} = CPD evaluated using evidence(domain(:,i)) +% If 'domains' is a single row vector, pots will be an object, not a cell array. + +ncases = size(domain,2); +assert(ncases==1); % not yet vectorized + +sz = dom_sizes(CPD); +ns = zeros(1, max(domain)); +ns(domain) = sz; + +local_ev = evidence(domain); +obs_bitv = ~isemptycell(local_ev); +odom = domain(obs_bitv); +T = convert_to_table(CPD, domain, local_ev, obs_bitv); + +switch pot_type + case 'u', + pot = upot(domain, sz, T, 0*myones(sz)); + case 'd', + ns(odom) = 1; + pot = dpot(domain, ns(domain), T); + case {'c','g'}, + % Since we want the output to be a Gaussian, the whole family must be observed. + % In other words, the potential is really just a constant. + p = T; + %p = prob_node(CPD, evidence(domain(end)), evidence(domain(1:end-1))); + ns(domain) = 0; + pot = cpot(domain, ns(domain), log(p)); + case 'cg', + T = T(:); + ns(odom) = 1; + can = cell(1, length(T)); + for i=1:length(T) + can{i} = cpot([], [], log(T(i))); + end + pot = cgpot(domain, [], ns, can); + otherwise, + error(['unrecognized pot type ' pot_type]) +end + diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/convert_to_table.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/convert_to_table.m new file mode 100644 index 00000000..65121122 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/convert_to_table.m @@ -0,0 +1,23 @@ +function T = convert_to_table(CPD, domain, local_ev, obs_bitv) +% CONVERT_TO_TABLE Convert a discrete CPD to a table +% function T = convert_to_table(CPD, domain, local_ev, obs_bitv) +% +% We convert the CPD to a CPT, and then lookup the evidence on the discrete parents. +% The resulting table can easily be converted to a potential. + + +CPT = CPD_to_CPT(CPD); +obs_child_only = ~any(obs_bitv(1:end-1)) & obs_bitv(end); + +if obs_child_only + sz = size(CPT); + CPT = reshape(CPT, prod(sz(1:end-1)), sz(end)); + o = local_ev{end}; + T = CPT(:, o); +else + odom = domain(obs_bitv); + vals = cat(1, local_ev{find(obs_bitv)}); % undo cell array + map = find_equiv_posns(odom, domain); + index = mk_multi_index(length(domain), map, vals); + T = CPT(index{:}); +end diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/prob_CPD.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/prob_CPD.m new file mode 100644 index 00000000..c0a79bda --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/prob_CPD.m @@ -0,0 +1,25 @@ +function p = prob_CPD(CPD, domain, ns, cnodes, evidence) +% PROB_CPD Compute prob of a node given evidence on the parents (discrete) +% p = prob_CPD(CPD, domain, ns, cnodes, evidence) +% +% domain is the domain of CPD. +% node_sizes(i) is the size of node i. +% cnodes = all the cts nodes +% evidence{i} is the evidence on the i'th node. + +ps = domain(1:end-1); +self = domain(end); +CPT = CPD_to_CPT(CPD); + +if isempty(ps) + T = CPT; +else + assert(~any(isemptycell(evidence(ps)))); + pvals = cat(1, evidence{ps}); + i = subv2ind(ns(ps), pvals(:)'); + T = reshape(CPT, [prod(ns(ps)) ns(self)]); + T = T(i,:); +end +p = T(evidence{self}); + + diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/prob_node.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/prob_node.m new file mode 100644 index 00000000..1a39fc79 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/Old/prob_node.m @@ -0,0 +1,51 @@ +function [P, p] = prob_node(CPD, self_ev, pev) +% PROB_NODE Compute prod_m P(x(i,m)| x(pi_i,m), theta_i) for node i (discrete) +% [P, p] = prob_node(CPD, self_ev, pev) +% +% self_ev(m) is the evidence on this node in case m. +% pev(i,m) is the evidence on the i'th parent in case m (if there are any parents). +% (These may also be cell arrays.) +% +% p(m) = P(x(i,m)| x(pi_i,m), theta_i) +% P = prod p(m) + +if iscell(self_ev), usecell = 1; else usecell = 0; end + +ncases = length(self_ev); +sz = dom_sizes(CPD); + +nparents = length(sz)-1; +if nparents == 0 + assert(isempty(pev)); +else + assert(isequal(size(pev), [nparents ncases])); +end + +n = length(sz); +dom = 1:n; +p = zeros(1, ncases); +if nparents == 0 + for m=1:ncases + if usecell + evidence = {self_ev{m}}; + else + evidence = num2cell(self_ev(m)); + end + T = convert_to_table(CPD, dom, evidence); + p(m) = T; + end +else + for m=1:ncases + if usecell + evidence = cell(1,n); + evidence(1:n-1) = pev(:,m); + evidence(n) = self_ev(m); + else + evidence = num2cell([pev(:,m)', self_ev(m)]); + end + T = convert_to_table(CPD, dom, evidence); + p(m) = T; + end +end +P = prod(p); + diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/README b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/README new file mode 100644 index 00000000..c0c5a3b3 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/README @@ -0,0 +1,5 @@ +Any CPD on a discrete child with discrete parents +can be represented as a table (although this might be quite big). +discrete_CPD uses this tabular representation to implement various +functions. Subtypes are free to implement more efficient versions. + diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_CPD_to_table_hidden_ps.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_CPD_to_table_hidden_ps.m new file mode 100644 index 00000000..c8f44f7e --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_CPD_to_table_hidden_ps.m @@ -0,0 +1,20 @@ +function T = convert_CPD_to_table_hidden_ps(CPD, child_obs) +% CONVERT_CPD_TO_TABLE_HIDDEN_PS Convert a discrete CPD to a table +% T = convert_CPD_to_table_hidden_ps(CPD, child_obs) +% +% This is like convert_to_table, except that we are guaranteed that +% none of the parents have evidence on them. +% child_obs may be an integer (1,2,...) or []. + +CPT = CPD_to_CPT(CPD); +if isempty(child_obs) + T = CPT(:); +else + sz = dom_sizes(CPD); + if length(sz)==1 % no parents + T = CPT(child_obs); + else + CPT = reshape(CPT, prod(sz(1:end-1)), sz(end)); + T = CPT(:, child_obs); + end +end diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_obs_CPD_to_table.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_obs_CPD_to_table.m new file mode 100644 index 00000000..04004088 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_obs_CPD_to_table.m @@ -0,0 +1,13 @@ +function T = convert_to_table(CPD, domain, evidence) +% CONVERT_TO_TABLE Convert a discrete CPD to a table +% T = convert_to_table(CPD, domain, evidence) +% +% We convert the CPD to a CPT, and then lookup the evidence on the discrete parents. +% The resulting table can easily be converted to a potential. + +CPT = CPD_to_CPT(CPD); +odom = domain(~isemptycell(evidence(domain))); +vals = cat(1, evidence{odom}); +map = find_equiv_posns(odom, domain); +index = mk_multi_index(length(domain), map, vals); +T = CPT(index{:}); diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_to_pot.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_to_pot.m new file mode 100644 index 00000000..ecc57d49 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_to_pot.m @@ -0,0 +1,62 @@ +function pot = convert_to_pot(CPD, pot_type, domain, evidence) +% CONVERT_TO_POT Convert a discrete CPD to a potential +% pot = convert_to_pot(CPD, pot_type, domain, evidence) +% +% pots = CPD evaluated using evidence(domain) + +ncases = size(domain,2); +assert(ncases==1); % not yet vectorized + +sz = dom_sizes(CPD); +ns = zeros(1, max(domain)); +ns(domain) = sz; + +CPT1 = CPD_to_CPT(CPD); +spar = issparse(CPT1); +odom = domain(~isemptycell(evidence(domain))); +if spar + T = convert_to_sparse_table(CPD, domain, evidence); +else + T = convert_to_table(CPD, domain, evidence); +end + +switch pot_type + case 'u', + pot = upot(domain, sz, T, 0*myones(sz)); + case 'd', + ns(odom) = 1; + pot = dpot(domain, ns(domain), T); + case {'c','g'}, + % Since we want the output to be a Gaussian, the whole family must be observed. + % In other words, the potential is really just a constant. + p = T; + %p = prob_node(CPD, evidence(domain(end)), evidence(domain(1:end-1))); + ns(domain) = 0; + pot = cpot(domain, ns(domain), log(p)); + + case 'cg', + T = T(:); + ns(odom) = 1; + can = cell(1, length(T)); + for i=1:length(T) + if T(i) == 0 + can{i} = cpot([], [], -Inf); % bug fix by Bob Welch 20/2/04 + else + can{i} = cpot([], [], log(T(i))); + end; + end + pot = cgpot(domain, [], ns, can); + + case 'scg' + T = T(:); + ns(odom) = 1; + pot_array = cell(1, length(T)); + for i=1:length(T) + pot_array{i} = scgcpot([], [], T(i)); + end + pot = scgpot(domain, [], [], ns, pot_array); + + otherwise, + error(['unrecognized pot type ' pot_type]) +end + diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_to_sparse_table.c b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_to_sparse_table.c new file mode 100644 index 00000000..369f5b7e --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_to_sparse_table.c @@ -0,0 +1,154 @@ +/* convert_to_sparse_table.c convert a sparse discrete CPD with evidence into sparse table */ +/* convert_to_pot.m located in ../CPDs/discrete_CPD call it */ +/* 3 input */ +/* CPD prhs[0] with 1D sparse CPT */ +/* domain prhs[1] */ +/* evidence prhs[2] */ +/* 1 output */ +/* T plhs[0] sparse table */ + +#include <math.h> +#include "mex.h" + +void ind_subv(int index, const int *cumprod, const int n, int *bsubv){ + int i; + + for (i = n-1; i >= 0; i--) { + bsubv[i] = ((int)floor(index / cumprod[i])); + index = index % cumprod[i]; + } +} + +int subv_ind(const int n, const int *cumprod, const int *subv){ + int i, index=0; + + for(i=0; i<n; i++){ + index += subv[i] * cumprod[i]; + } + return index; +} + +void reset_nzmax(mxArray *spArray, const int old_nzmax, const int new_nzmax){ + double *ptr; + void *newptr; + int *ir, *jc; + int nbytes; + + if(new_nzmax == old_nzmax) return; + nbytes = new_nzmax * sizeof(*ptr); + ptr = mxGetPr(spArray); + newptr = mxRealloc(ptr, nbytes); + mxSetPr(spArray, newptr); + nbytes = new_nzmax * sizeof(*ir); + ir = mxGetIr(spArray); + newptr = mxRealloc(ir, nbytes); + mxSetIr(spArray, newptr); + jc = mxGetJc(spArray); + jc[0] = 0; + jc[1] = new_nzmax; + mxSetNzmax(spArray, new_nzmax); +} + + +void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]){ + int i, j, NS, NZB, count, bdim, match, domain, bindex, sindex, nzCounts=0; + int *observed, *bsubv, *ssubv, *bir, *sir, *bjc, *sjc, *mask, *ssize, *bcumprod, *scumprod; + double *pDomain, *pSize, *bpr, *spr; + mxArray *pTemp; + + pTemp = mxGetField(prhs[0], 0, "CPT"); + bpr = mxGetPr(pTemp); + bir = mxGetIr(pTemp); + bjc = mxGetJc(pTemp); + NZB = bjc[1]; + pTemp = mxGetField(prhs[0], 0, "sizes"); + pSize = mxGetPr(pTemp); + + pDomain = mxGetPr(prhs[1]); + bdim = mxGetNumberOfElements(prhs[1]); + + mask = malloc(bdim * sizeof(int)); + ssize = malloc(bdim * sizeof(int)); + observed = malloc(bdim * sizeof(int)); + + for(i=0; i<bdim; i++){ + ssize[i] = (int)pSize[i]; + } + + count = 0; + for(i=0; i<bdim; i++){ + domain = (int)pDomain[i] - 1; + pTemp = mxGetCell(prhs[2], domain); + if(pTemp){ + mask[count] = i; + ssize[i] = 1; + observed[count] = (int)mxGetScalar(pTemp) - 1; + count++; + } + } + + if(count == 0){ + pTemp = mxGetField(prhs[0], 0, "CPT"); + plhs[0] = mxDuplicateArray(pTemp); + free(mask); + free(ssize); + free(observed); + return; + } + + bsubv = malloc(bdim * sizeof(int)); + ssubv = malloc(count * sizeof(int)); + bcumprod = malloc(bdim * sizeof(int)); + scumprod = malloc(bdim * sizeof(int)); + + NS = 1; + for(i=0; i<bdim; i++){ + NS *= ssize[i]; + } + + plhs[0] = mxCreateSparse(NS, 1, NS, mxREAL); + spr = mxGetPr(plhs[0]); + sir = mxGetIr(plhs[0]); + sjc = mxGetJc(plhs[0]); + sjc[0] = 0; + sjc[1] = NS; + + bcumprod[0] = 1; + scumprod[0] = 1; + for(i=0; i<bdim-1; i++){ + bcumprod[i+1] = bcumprod[i] * (int)pSize[i]; + scumprod[i+1] = scumprod[i] * ssize[i]; + } + + nzCounts = 0; + for(i=0; i<NZB; i++){ + bindex = bir[i]; + ind_subv(bindex, bcumprod, bdim, bsubv); + for(j=0; j<count; j++){ + ssubv[j] = bsubv[mask[j]]; + } + match = 1; + for(j=0; j<count; j++){ + if((ssubv[j]) != observed[j]){ + match = 0; + break; + } + } + if(match){ + spr[nzCounts] = bpr[i]; + sindex = subv_ind(bdim, scumprod, bsubv); + sir[nzCounts] = sindex; + nzCounts++; + } + } + + reset_nzmax(plhs[0], NS, nzCounts); + free(mask); + free(ssize); + free(observed); + free(bsubv); + free(ssubv); + free(bcumprod); + free(scumprod); +} + diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_to_table.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_to_table.m new file mode 100644 index 00000000..dc5bcd40 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/convert_to_table.m @@ -0,0 +1,15 @@ +function T = convert_to_table(CPD, domain, evidence) +% CONVERT_TO_TABLE Convert a discrete CPD to a table +% T = convert_to_table(CPD, domain, evidence) +% +% We convert the CPD to a CPT, and then lookup the evidence on the discrete parents. +% The resulting table can easily be converted to a potential. + +domain = domain(:); +CPT = CPD_to_CPT(CPD); +odom = domain(~isemptycell(evidence(domain))); +vals = cat(1, evidence{odom}); +map = find_equiv_posns(odom, domain); +index = mk_multi_index(length(domain), map, vals); +T = CPT(index{:}); +T = T(:); diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/discrete_CPD.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/discrete_CPD.m new file mode 100644 index 00000000..b4250831 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/discrete_CPD.m @@ -0,0 +1,6 @@ +function CPD = discrete_CPD(clamped, dom_sizes) +% DISCRETE_CPD Virtual constructor for generic discrete CPD +% CPD = discrete_CPD(clamped, dom_sizes) + +CPD.dom_sizes = dom_sizes; +CPD = class(CPD, 'discrete_CPD', generic_CPD(clamped)); diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/dom_sizes.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/dom_sizes.m new file mode 100644 index 00000000..2ee750de --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/dom_sizes.m @@ -0,0 +1,5 @@ +function sz = dom_sizes(CPD) +% DOM_SIZES Return the size of each node in the domain +% sz = dom_sizes(CPD) + +sz = CPD.dom_sizes; diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/log_prob_node.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/log_prob_node.m new file mode 100644 index 00000000..315464a9 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/log_prob_node.m @@ -0,0 +1,12 @@ +function L = log_prob_node(CPD, self_ev, pev) +% LOG_PROB_NODE Compute sum_m log P(x(i,m)| x(pi_i,m), theta_i) for node i (discrete) +% L = log_prob_node(CPD, self_ev, pev) +% +% self_ev(m) is the evidence on this node in case m. +% pev(i,m) is the evidence on the i'th parent in case m (if there are any parents). +% (These may also be cell arrays.) + +[P, p] = prob_node(CPD, self_ev, pev); % P may underflow, so we use p +tiny = exp(-700); +p = p + (p==0)*tiny; % replace 0s by tiny +L = sum(log(p)); diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/private/CVS/Entries b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/private/CVS/Entries new file mode 100644 index 00000000..da678df7 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/private/CVS/Entries @@ -0,0 +1,2 @@ +/prod_CPT_and_pi_msgs.m/1.1.1.1/Wed May 29 15:59:52 2002// +D diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/private/CVS/Repository b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/private/CVS/Repository new file mode 100644 index 00000000..2b3c1c9d --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/private/CVS/Repository @@ -0,0 +1 @@ +FullBNT/BNT/CPDs/@discrete_CPD/private diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/private/CVS/Root b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/private/CVS/Root new file mode 100644 index 00000000..f3bd14a6 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/private/CVS/Root @@ -0,0 +1 @@ +:ext:nsaunier@bnt.cvs.sourceforge.net:/cvsroot/bnt diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/private/prod_CPT_and_pi_msgs.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/private/prod_CPT_and_pi_msgs.m new file mode 100644 index 00000000..fe8f6a20 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/private/prod_CPT_and_pi_msgs.m @@ -0,0 +1,18 @@ +function T = prod_CPT_and_pi_msgs(CPD, n, ps, msgs, except) +% PROD_CPT_AND_PI_MSGS Multiply the CPD and all the pi messages from parents, perhaps excepting one +% T = prod_CPY_and_pi_msgs(CPD, n, ps, msgs, except) + +if nargin < 5, except = -1; end + +dom = [ps n]; +%ns = sparse(1, max(dom)); +ns = zeros(1, max(dom)); +CPT = CPD_to_CPT(CPD); +ns(dom) = mysize(CPT); +T = dpot(dom, ns(dom), CPT); +for i=1:length(ps) + p = ps(i); + if p ~= except + T = multiply_by_pot(T, dpot(p, ns(p), msgs{n}.pi_from_parent{i})); + end +end diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/prob_node.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/prob_node.m new file mode 100644 index 00000000..275870c8 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/prob_node.m @@ -0,0 +1,81 @@ +function [P, p] = prob_node(CPD, self_ev, pev) +% PROB_NODE Compute prod_m P(x(i,m)| x(pi_i,m), theta_i) for node i (discrete) +% [P, p] = prob_node(CPD, self_ev, pev) +% +% self_ev(m) is the evidence on this node in case m. +% pev(i,m) is the evidence on the i'th parent in case m (if there are any parents). +% (These may also be cell arrays.) +% +% p(m) = P(x(i,m)| x(pi_i,m), theta_i) +% P = prod p(m) + +if iscell(self_ev), usecell = 1; else usecell = 0; end + +ncases = length(self_ev); +sz = dom_sizes(CPD); + +nparents = length(sz)-1; +if nparents == 0 + assert(isempty(pev)); +else + assert(isequal(size(pev), [nparents ncases])); +end + +n = length(sz); +dom = 1:n; +p = zeros(1, ncases); +if isa(CPD, 'tabular_CPD') + % speed up by looking up CPT using index Zhang Yimin 2001-12-31 + if usecell + if nparents == 0 + data = [cell2num(self_ev)]; + else + data = [cell2num(pev); cell2num(self_ev)]; + end + else + if nparents == 0 + data = [self_ev]; + else + data = [pev; self_ev]; + end + end + + indices = subv2ind(sz, data'); % each row of data' is a case + + CPT=CPD_to_CPT(CPD); + p = CPT(indices); + + %get the prob list + %cpt_size = prod(sz); + %prob_list=reshape(CPT, cpt_size, 1); + %for m=1:ncases %here we assume we get evidence for node and all its parents + % idx=indices(m); + % p(m)=prob_list(idx); + %end + +else % eg. softmax + + for m=1:ncases + if usecell + if nparents == 0 + evidence = {self_ev{m}}; + else + evidence = cell(1,n); + evidence(1:n-1) = pev(:,m); + evidence(n) = self_ev(m); + end + else + if nparents == 0 + evidence = num2cell(self_ev(m)); + else + evidence = num2cell([pev(:,m)', self_ev(m)]); + end + end + T = convert_to_table(CPD, dom, evidence); + p(m) = T; + end +end + +P = prod(p); + + diff --git a/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/sample_node.m b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/sample_node.m new file mode 100644 index 00000000..9e0ed994 --- /dev/null +++ b/sourcecodes/bnt-master/BNT/CPDs/@discrete_CPD/sample_node.m @@ -0,0 +1,34 @@ +function y = sample_node(CPD, pvals) +% SAMPLE_NODE Draw a random sample from P(Xi | x(pi_i), theta_i) (discrete) +% y = sample_node(CPD, parent_evidence) +% +% parent_evidence{i} is the value of the i'th parent + +if 0 +n = length(pvals)+1; +dom = 1:n; +evidence = cell(1,n); +evidence(1:n-1) = pvals; +T = convert_to_table(CPD, dom, evidence); +y = sample_discrete(T); +end + + +CPT = CPD_to_CPT(CPD); +sz = mysize(CPT); +nparents = length(sz)-1; +switch nparents + case 0, T = CPT; + case 1, T = CPT(pvals{1}, :); + case 2, T = CPT(pvals{1}, pvals{2}, :); + case 3, T = CPT(pvals{1}, pvals{2}, pvals{3}, :); + case 4, T = CPT(pvals{1}, pvals{2}, pvals{3}, pvals{4}, :); + otherwise, + pvals = cat(1, pvals{:}); + psz = sz(1:end-1); + ssz = sz(end); + i = subv2ind(psz, pvals(:)'); + T = reshape(CPT, [prod(psz) ssz]); + T = T(i,:); +end +y = sample_discrete(T); |