annotate toolboxes/FullBNT-1.0.7/bnt/general/mk_dbn.m @ 0:cc4b1211e677 tip

initial commit to HG from Changeset: 646 (e263d8a21543) added further path and more save "camirversion.m"
author Daniel Wolff
date Fri, 19 Aug 2016 13:07:06 +0200
parents
children
rev   line source
Daniel@0 1 function bnet = mk_dbn(intra, inter, node_sizes, varargin)
Daniel@0 2 % MK_DBN Make a Dynamic Bayesian Network.
Daniel@0 3 %
Daniel@0 4 % BNET = MK_DBN(INTRA, INTER, NODE_SIZES, ...) makes a DBN with arcs
Daniel@0 5 % from i in slice t to j in slice t iff intra(i,j) = 1, and
Daniel@0 6 % from i in slice t to j in slice t+1 iff inter(i,j) = 1,
Daniel@0 7 % for i,j in {1, 2, ..., n}, where n = num. nodes per slice, and t >= 1.
Daniel@0 8 % node_sizes(i) is the number of values node i can take on.
Daniel@0 9 % The nodes are assumed to be in topological order. Use TOPOLOGICAL_SORT if necessary.
Daniel@0 10 % See also mk_bnet.
Daniel@0 11 %
Daniel@0 12 % Optional arguments [default in brackets]
Daniel@0 13 % 'discrete' - list of discrete nodes [1:n]
Daniel@0 14 % 'observed' - the list of nodes which will definitely be observed in every slice of every case [ [] ]
Daniel@0 15 % 'eclass1' - equiv class for slice 1 [1:n]
Daniel@0 16 % 'eclass2' - equiv class for slice 2 [tie nodes with equivalent parents to slice 1]
Daniel@0 17 % equiv_class1(i) = j means node i in slice 1 gets its parameters from bnet.CPD{j},
Daniel@0 18 % i.e., nodes i and j have tied parameters.
Daniel@0 19 % 'intra1' - topology of first slice, if different from others
Daniel@0 20 % 'names' - a cell array of strings to be associated with nodes 1:n [{}]
Daniel@0 21 % This creates an associative array, so you write e.g.
Daniel@0 22 % 'evidence(bnet.names{'bar'}) = 42' instead of 'evidence(2} = 42'
Daniel@0 23 % assuming names = { 'foo', 'bar', ...}.
Daniel@0 24 %
Daniel@0 25 % For backwards compatibility with BNT2, arguments can also be specified as follows
Daniel@0 26 % bnet = mk_dbn(intra, inter, node_sizes, dnodes, eclass1, eclass2, intra1)
Daniel@0 27 %
Daniel@0 28 % After calling this function, you must specify the parameters (conditional probability
Daniel@0 29 % distributions) using bnet.CPD{i} = gaussian_CPD(...) or tabular_CPD(...) etc.
Daniel@0 30
Daniel@0 31
Daniel@0 32 n = length(intra);
Daniel@0 33 ss = n;
Daniel@0 34 bnet.nnodes_per_slice = ss;
Daniel@0 35 bnet.intra = intra;
Daniel@0 36 bnet.inter = inter;
Daniel@0 37 bnet.intra1 = intra;
Daniel@0 38 dag = zeros(2*n);
Daniel@0 39 dag(1:n,1:n) = bnet.intra1;
Daniel@0 40 dag(1:n,(1:n)+n) = bnet.inter;
Daniel@0 41 dag((1:n)+n,(1:n)+n) = bnet.intra;
Daniel@0 42 bnet.dag = dag;
Daniel@0 43 bnet.names = {};
Daniel@0 44
Daniel@0 45 directed = 1;
Daniel@0 46 if ~acyclic(dag,directed)
Daniel@0 47 error('graph must be acyclic')
Daniel@0 48 end
Daniel@0 49
Daniel@0 50
Daniel@0 51 bnet.eclass1 = 1:n;
Daniel@0 52 %bnet.eclass2 = (1:n)+n;
Daniel@0 53 bnet.eclass2 = bnet.eclass1;
Daniel@0 54 for i=1:ss
Daniel@0 55 if isequal(parents(dag, i+ss), parents(dag, i)+ss)
Daniel@0 56 %fprintf('%d has isomorphic parents, eclass %d\n', i, bnet.eclass2(i))
Daniel@0 57 else
Daniel@0 58 bnet.eclass2(i) = max(bnet.eclass2) + 1;
Daniel@0 59 %fprintf('%d has non isomorphic parents, eclass %d\n', i, bnet.eclass2(i))
Daniel@0 60 end
Daniel@0 61 end
Daniel@0 62
Daniel@0 63 dnodes = 1:n;
Daniel@0 64 bnet.observed = [];
Daniel@0 65
Daniel@0 66 if nargin >= 4
Daniel@0 67 args = varargin;
Daniel@0 68 nargs = length(args);
Daniel@0 69 if ~isstr(args{1})
Daniel@0 70 if nargs >= 1, dnodes = args{1}; end
Daniel@0 71 if nargs >= 2, bnet.eclass1 = args{2}; end
Daniel@0 72 if nargs >= 3, bnet.eclass2 = args{3}; end
Daniel@0 73 if nargs >= 4, bnet.intra1 = args{4}; end
Daniel@0 74 else
Daniel@0 75 for i=1:2:nargs
Daniel@0 76 switch args{i},
Daniel@0 77 case 'discrete', dnodes = args{i+1};
Daniel@0 78 case 'observed', bnet.observed = args{i+1};
Daniel@0 79 case 'eclass1', bnet.eclass1 = args{i+1};
Daniel@0 80 case 'eclass2', bnet.eclass2 = args{i+1};
Daniel@0 81 case 'intra1', bnet.intra1 = args{i+1};
Daniel@0 82 %case 'ar_hmm', bnet.ar_hmm = args{i+1}; % should check topology
Daniel@0 83 case 'names', bnet.names = assocarray(args{i+1}, num2cell(1:n));
Daniel@0 84 otherwise,
Daniel@0 85 error(['invalid argument name ' args{i}]);
Daniel@0 86 end
Daniel@0 87 end
Daniel@0 88 end
Daniel@0 89 end
Daniel@0 90
Daniel@0 91
Daniel@0 92 bnet.observed = sort(bnet.observed); % for comparing sets
Daniel@0 93 ns = node_sizes;
Daniel@0 94 bnet.node_sizes_slice = ns(:)';
Daniel@0 95 bnet.node_sizes = [ns(:) ns(:)];
Daniel@0 96
Daniel@0 97 cnodes = mysetdiff(1:n, dnodes);
Daniel@0 98 bnet.dnodes_slice = dnodes;
Daniel@0 99 bnet.cnodes_slice = cnodes;
Daniel@0 100 bnet.dnodes = [dnodes dnodes+n];
Daniel@0 101 bnet.cnodes = [cnodes cnodes+n];
Daniel@0 102
Daniel@0 103 bnet.equiv_class = [bnet.eclass1(:) bnet.eclass2(:)];
Daniel@0 104 bnet.CPD = cell(1,max(bnet.equiv_class(:)));
Daniel@0 105 eclass = bnet.equiv_class(:);
Daniel@0 106 E = max(eclass);
Daniel@0 107 bnet.rep_of_eclass = zeros(1,E);
Daniel@0 108 for e=1:E
Daniel@0 109 mems = find(eclass==e);
Daniel@0 110 bnet.rep_of_eclass(e) = mems(1);
Daniel@0 111 end
Daniel@0 112
Daniel@0 113 ss = n;
Daniel@0 114 onodes = bnet.observed;
Daniel@0 115 hnodes = mysetdiff(1:ss, onodes);
Daniel@0 116 bnet.hidden_bitv = zeros(1,2*ss);
Daniel@0 117 bnet.hidden_bitv(hnodes) = 1;
Daniel@0 118 bnet.hidden_bitv(hnodes+ss) = 1;
Daniel@0 119
Daniel@0 120 bnet.parents = cell(1, 2*ss);
Daniel@0 121 for i=1:ss
Daniel@0 122 bnet.parents{i} = parents(bnet.dag, i);
Daniel@0 123 bnet.parents{i+ss} = parents(bnet.dag, i+ss);
Daniel@0 124 end
Daniel@0 125
Daniel@0 126 bnet.auto_regressive = zeros(1,ss);
Daniel@0 127 % ar(i)=1 means (observed) node i depends on i in the previous slice
Daniel@0 128 for o=bnet.observed(:)'
Daniel@0 129 if any(bnet.parents{o+ss} <= ss)
Daniel@0 130 bnet.auto_regressive(o) = 1;
Daniel@0 131 end
Daniel@0 132 end
Daniel@0 133