BNW using Octave instead of Matlab.

This version of BNW should perform the same as the original version. The only difference is that it uses Octave instead of Matlab when running BayesNet Toolbox during parameter learning.

I am calling this BNW_1.02. It can be accessed at:
compbio.uthsc.edu/BNW_1.02

10 files changed, 462 insertions, 0 deletions

diff --git a/sourcecodes/bnt-master/BNT/examples/dynamic/Old/CVS/Entries b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/CVS/Entries
new file mode 100644
index 00000000..b5b0dc43
--- /dev/null
+++ b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/CVS/Entries
@@ -0,0 +1,8 @@
+/chmm1.m/1.1.1.1/Wed May 29 15:59:54 2002//
+/cmp_inference.m/1.1.1.1/Wed May 29 15:59:54 2002//
+/kalman1.m/1.1.1.1/Wed May 29 15:59:54 2002//
+/old.water1.m/1.1.1.1/Wed May 29 15:59:54 2002//
+/online1.m/1.1.1.1/Wed May 29 15:59:54 2002//
+/online2.m/1.1.1.1/Wed May 29 15:59:54 2002//
+/scg_dbn.m/1.1.1.1/Wed May 29 15:59:54 2002//
+D
diff --git a/sourcecodes/bnt-master/BNT/examples/dynamic/Old/CVS/Repository b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/CVS/Repository
new file mode 100644
index 00000000..b23f8f6d
--- /dev/null
+++ b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/CVS/Repository
@@ -0,0 +1 @@
+FullBNT/BNT/examples/dynamic/Old
diff --git a/sourcecodes/bnt-master/BNT/examples/dynamic/Old/CVS/Root b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/CVS/Root
new file mode 100644
index 00000000..f3bd14a6
--- /dev/null
+++ b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/CVS/Root
@@ -0,0 +1 @@
+:ext:nsaunier@bnt.cvs.sourceforge.net:/cvsroot/bnt
diff --git a/sourcecodes/bnt-master/BNT/examples/dynamic/Old/chmm1.m b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/chmm1.m
new file mode 100644
index 00000000..d4195c97
--- /dev/null
+++ b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/chmm1.m
@@ -0,0 +1,40 @@
+% Compare the speeds of various inference engines on a coupled HMM
+
+N = 2;
+Q = 2;
+rand('state', 0);
+randn('state', 0);
+discrete = 1;
+if discrete
+  Y = 2; % size of output alphabet
+else
+  Y = 1;
+end
+coupled = 1;
+[bnet, onodes] = mk_chmm(N, Q, Y, discrete, coupled);
+ss = N*2;
+
+T = 3;
+
+
+engine = {};
+tic; engine{end+1} = jtree_dbn_inf_engine(bnet, 'observed', onodes);  toc
+%tic; engine{end+1} = jtree_ndxSD_dbn_inf_engine(bnet, onodes);  toc
+%tic; engine{end+1} = jtree_ndxB_dbn_inf_engine(bnet, onodes);  toc
+engine{end+1} = hmm_inf_engine(bnet, onodes);
+%engine{end+1} = dhmm_inf_engine(bnet, onodes);
+tic; engine{end+1} = jtree_unrolled_dbn_inf_engine(bnet, T, onodes);  toc
+
+%engine{end+1} = bk_inf_engine(bnet, 'ff', onodes);
+%engine{end+1} = loopy_dbn_inf_engine(bnet, onodes);
+
+exact = [1 2 3];
+
+filter = 0;
+single = 0;
+maximize = 0;
+
+[err, time, engine] = cmp_inference(bnet, onodes, engine, exact, T, filter, single, maximize);
+%err = cmp_learning(bnet, onodes, engine, exact, T);
+
+
diff --git a/sourcecodes/bnt-master/BNT/examples/dynamic/Old/cmp_inference.m b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/cmp_inference.m
new file mode 100644
index 00000000..b5c936f0
--- /dev/null
+++ b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/cmp_inference.m
@@ -0,0 +1,75 @@
+function [err, time, engine] = cmp_inference(bnet, engine, exact, T, filter, singletons, maximize)
+% CMP_INFERENCE Compare several inference engines on a DBN
+% [err, time, engine] = cmp_inference(bnet, engine, exact, T, filter, singletons, maximize)
+%
+% engine{i} is the i'th inference engine.
+% 'exact' specifies which engines do exact inference - 
+%   we check that these all give the same results.
+% 'T' is the length of the random sequence we generate.
+% If filter=1, we do filtering, else smoothing (default: smoothing)
+% If singletons=1, we compare marginal_nodes, else marginal_family (default: family)
+%
+% err(e,n,t) = sum_i | Pr_exact(X(n,t)=i) - Pr_e(X(n,t)=i) |
+%   where Pr_e = prob. according to engine e
+% time(e) = elapsed time for doing inference with engine e
+
+err = [];
+
+if nargin < 5, filter = 0; end
+if nargin < 6, singletons = 0; end
+if nargin < 7, maximize = 0; end
+
+check_ll = 1;
+
+assert(~maximize);
+
+E = length(engine);
+ref = exact(1); % reference
+
+ss = length(bnet.intra);
+ev = sample_dbn(bnet, 'length', T);
+evidence = cell(ss,T);
+onodes = bnet.observed;
+evidence(onodes,:) = ev(onodes, :);
+
+assert(~filter);
+for i=1:E
+  tic;
+  %[engine{i}, ll(i)] = enter_evidence(engine{i}, evidence, 'maximize', maximize);
+  [engine{i}, ll(i)] = enter_evidence(engine{i}, evidence);
+  time(i)=toc;
+  fprintf('engine %d took %6.4f seconds\n', i, time(i));
+end
+
+cmp = mysetdiff(exact, ref);
+if check_ll
+for i=cmp(:)'
+  if ~approxeq(ll(ref), ll(i))
+    error(['engine ' num2str(i) ' has wrong ll'])
+  end
+end
+end
+ll
+
+hnodes = mysetdiff(1:ss, onodes);
+m = cell(1,E);
+for t=1:T
+  for n=hnodes(:)'
+    for e=1:E
+      if singletons
+	m{e} = marginal_nodes(engine{e}, n, t);
+      else
+	m{e} = marginal_family(engine{e}, n, t);
+      end
+    end
+    for e=1:E
+      assert(isequal(m{e}.domain, m{ref}.domain));
+    end
+    for e=cmp(:)'
+      if ~approxeq(m{ref}.T(:), m{e}.T(:))
+	str= sprintf('engine %d is wrong; n=%d, t=%d', e, n, t);
+	error(str)
+      end
+    end
+  end
+end
diff --git a/sourcecodes/bnt-master/BNT/examples/dynamic/Old/kalman1.m b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/kalman1.m
new file mode 100644
index 00000000..c068ab3e
--- /dev/null
+++ b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/kalman1.m
@@ -0,0 +1,127 @@
+% Make a linear dynamical system
+%   X1 -> X2
+%   |     | 
+%   v     v
+%   Y1    Y2 
+
+intra = zeros(2);
+intra(1,2) = 1;
+inter = zeros(2);
+inter(1,1) = 1;
+n = 2;
+
+X = 2; % size of hidden state
+Y = 2; % size of observable state
+
+ns = [X Y];
+dnodes = [];
+onodes = [2];
+eclass1 = [1 2];
+eclass2 = [3 2];
+bnet = mk_dbn(intra, inter, ns, 'discrete', dnodes, 'eclass1', eclass1, 'eclass2', eclass2, ...
+	      'observed', onodes);
+
+x0 = rand(X,1);
+V0 = eye(X);
+C0 = rand(Y,X);
+R0 = eye(Y);
+A0 = rand(X,X);
+Q0 = eye(X);
+
+bnet.CPD{1} = gaussian_CPD(bnet, 1, 'mean', x0, 'cov', V0, 'cov_prior_weight', 0);
+bnet.CPD{2} = gaussian_CPD(bnet, 2, 'mean', zeros(Y,1), 'cov', R0, 'weights', C0, ...
+			   'clamp_mean', 1, 'cov_prior_weight', 0);
+bnet.CPD{3} = gaussian_CPD(bnet, 3, 'mean', zeros(X,1), 'cov', Q0, 'weights', A0, ...
+			   'clamp_mean', 1, 'cov_prior_weight', 0);
+
+
+T = 5; % fixed length sequences
+
+clear engine;
+engine{1} = kalman_inf_engine(bnet);
+engine{2} = jtree_unrolled_dbn_inf_engine(bnet, T);
+engine{3} = jtree_dbn_inf_engine(bnet);
+N = length(engine);
+
+% inference
+
+ev = sample_dbn(bnet, T);
+evidence = cell(n,T);
+evidence(onodes,:) = ev(onodes, :);
+
+t = 1;
+query = [1 3];
+m = cell(1, N);
+ll = zeros(1, N);
+for i=1:N
+  [engine{i}, ll(i)] = enter_evidence(engine{i}, evidence);
+  m{i} = marginal_nodes(engine{i}, query, t);
+end
+
+% compare all engines to engine{1}
+for i=2:N
+  assert(approxeq(m{1}.mu, m{i}.mu));
+  assert(approxeq(m{1}.Sigma, m{i}.Sigma));
+  assert(approxeq(ll(1), ll(i)));
+end
+
+if 0
+for i=2:N
+  approxeq(m{1}.mu, m{i}.mu)
+  approxeq(m{1}.Sigma, m{i}.Sigma)
+  approxeq(ll(1), ll(i))
+end
+end
+
+% learning
+
+ncases = 5;
+cases = cell(1, ncases);
+for i=1:ncases
+  ev = sample_dbn(bnet, T);
+  cases{i} = cell(n,T);
+  cases{i}(onodes,:) = ev(onodes, :);
+end
+
+max_iter = 2;
+bnet2 = cell(1,N);
+LLtrace = cell(1,N);
+for i=1:N
+  [bnet2{i}, LLtrace{i}] = learn_params_dbn_em(engine{i}, cases, 'max_iter', max_iter);
+end
+
+for i=1:N
+  temp = bnet2{i};
+  for e=1:3
+    CPD{i,e} = struct(temp.CPD{e});
+  end
+end
+
+for i=2:N
+  assert(approxeq(LLtrace{i}, LLtrace{1}));
+  for e=1:3
+    assert(approxeq(CPD{i,e}.mean, CPD{1,e}.mean));
+    assert(approxeq(CPD{i,e}.cov, CPD{1,e}.cov));
+    assert(approxeq(CPD{i,e}.weights, CPD{1,e}.weights));
+  end
+end
+
+
+% Compare to KF toolbox
+
+data = zeros(Y, T, ncases);
+for i=1:ncases
+  data(:,:,i) = cell2num(cases{i}(onodes, :));
+end   
+[A2, C2, Q2, R2, x2, V2, LL2trace] =  learn_kalman(data, A0, C0, Q0, R0, x0, V0, max_iter);
+
+
+e = 1;
+assert(approxeq(x2, CPD{e,1}.mean))
+assert(approxeq(V2, CPD{e,1}.cov))
+assert(approxeq(C2, CPD{e,2}.weights))
+assert(approxeq(R2, CPD{e,2}.cov));
+assert(approxeq(A2, CPD{e,3}.weights))
+assert(approxeq(Q2, CPD{e,3}.cov));
+assert(approxeq(LL2trace, LLtrace{1}))
+
diff --git a/sourcecodes/bnt-master/BNT/examples/dynamic/Old/old.water1.m b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/old.water1.m
new file mode 100644
index 00000000..0a356ef8
--- /dev/null
+++ b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/old.water1.m
@@ -0,0 +1,48 @@
+% Compare the speeds of various inference engines on the water DBN
+
+[bnet, onodes] = mk_water_dbn;
+
+T = 3;
+
+engine = {};
+engine{end+1} = jtree_unrolled_dbn_inf_engine(bnet, T, onodes);
+engine{end+1} = hmm_inf_engine(bnet, onodes);
+engine{end+1} = frontier_inf_engine(bnet, onodes);
+engine{end+1} = jtree_dbn_inf_engine(bnet, onodes);
+engine{end+1} = bk_inf_engine(bnet, 'exact', onodes);
+
+engine{end+1} = bk_inf_engine(bnet, 'ff', onodes);
+engine{end+1} = bk_inf_engine(bnet, { [1 2], [3 4 5 6], [7 8] }, onodes);
+
+N = length(engine);
+exact = 1:5;
+
+
+filter = 0;
+err = cmp_inference(bnet, onodes, engine, exact, T, filter);
+
+% elapsed times for enter_evidence  (matlab 5.3 on PIII with 256MB running Redhat linux)  
+
+% T = 5, 4/20/00
+%     0.6266 unrolled *
+%     0.3490 hmm *
+%     1.1743 frontier
+%     1.4621 old frontier
+%     0.3270 fast frontier *
+%     1.3926 jtree 
+%     1.3790 bk 
+%     0.4916 fast bk
+%     0.4190 fast bk compiled
+%     0.3574 fast jtree *
+
+
+err = cmp_learning(bnet, onodes, engine, exact, T);
+
+% elapsed times for learn_params_dbn_em (matlab 5.3 on PIII with 256MB running Redhat linux)  
+
+% T = 5, 2cases, 2 iter, 4/20/00
+% 3.5750 unrolled
+% 3.7475 hmm
+% 2.1452 fast frontier
+% 2.5724 fast bk compiled
+% 2.3387 fast jtree
diff --git a/sourcecodes/bnt-master/BNT/examples/dynamic/Old/online1.m b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/online1.m
new file mode 100644
index 00000000..05708e10
--- /dev/null
+++ b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/online1.m
@@ -0,0 +1,59 @@
+% Check that online inference gives same results as filtering for various algorithms
+
+N = 3;
+Q = 2;
+ss = N*2;
+
+rand('state', 0);
+randn('state', 0);
+
+
+obs_size = 1;
+discrete_obs = 0;
+bnet = mk_chmm(N, Q, obs_size, discrete_obs);
+ns = bnet.node_sizes_slice;
+
+engine = {};
+engine{end+1} = hmm_inf_engine(bnet);
+E = length(engine);
+
+onodes = (1:N)+N;
+
+T = 4;
+ev = cell(ss,T);
+ev(onodes,:) = num2cell(randn(N, T));
+
+
+filter = 1;
+loglik2 = zeros(1,E);
+for e=1:E
+  [engine2{e}, loglik2(e)] = enter_evidence(engine{e}, ev, 'filter', filter);
+end
+
+loglik = zeros(1,E);
+marg1 = cell(E,N,T);
+for e=1:E
+  ll = zeros(1,T);
+  engine{e} = dbn_init_bel(engine{e});
+  for t=1:T
+    [engine{e}, ll(t)] = dbn_update_bel(engine{e}, ev(:,t), t);
+    for i=1:N
+      marg1{e,i,t} = dbn_marginal_from_bel(engine{e}, i);
+    end
+  end
+  loglik1(e) = sum(ll);
+end
+
+assert(approxeq(loglik1, loglik2))
+
+a = zeros(E,N,T);
+for e=1:E
+  for t=1:T
+    for i=1:N
+      marg2{e,i,t} = marginal_nodes(engine2{e}, i, t);
+      a(e,i,t) = (approxeq(marg2{e,i,t}.T(:), marg1{e,i,t}.T(:)));
+    end
+  end
+end
+
+assert(all(a(:)==1))
diff --git a/sourcecodes/bnt-master/BNT/examples/dynamic/Old/online2.m b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/online2.m
new file mode 100644
index 00000000..6b141f19
--- /dev/null
+++ b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/online2.m
@@ -0,0 +1,33 @@
+N = 1; % regular HMM
+Q = 2;
+ss = 2;
+hnodes = 1;
+onodes = 2;
+
+rand('state', 0);
+randn('state', 0);
+O = 2;
+discrete_obs = 1;
+bnet = mk_chmm(N, Q, O, discrete_obs);
+ns = bnet.node_sizes_slice;
+
+engine = hmm_inf_engine(bnet, onodes);
+
+T = 4;
+ev = cell(ss,T);
+ev(onodes,:) = num2cell(sample_discrete([0.5 0.5], N, T));
+
+
+engine = dbn_init_bel(engine);
+for t=1:T
+  if t==1
+    [engine, ll(t)] = dbn_update_bel1(engine, ev(:,t));
+  else
+    [engine, ll(t)] = dbn_update_bel(engine, ev(:,t-1:t));
+  end
+  % one-step ahead prediction
+  lag = 1;
+  engine2 = dbn_predict_bel(engine, lag);  
+  marg = dbn_marginal_from_bel(engine2, 1)
+  marg = dbn_marginal_from_bel(engine2, 2)
+end
diff --git a/sourcecodes/bnt-master/BNT/examples/dynamic/Old/scg_dbn.m b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/scg_dbn.m
new file mode 100644
index 00000000..0ddabb34
--- /dev/null
+++ b/sourcecodes/bnt-master/BNT/examples/dynamic/Old/scg_dbn.m
@@ -0,0 +1,70 @@
+% to test whether scg inference engine can handl dynameic BN
+% Make a linear dynamical system
+%   X1 -> X2
+%   |     | 
+%   v     v
+%   Y1    Y2 
+
+intra = zeros(2);
+intra(1,2) = 1;
+inter = zeros(2);
+inter(1,1) = 1;
+n = 2;
+
+X = 2; % size of hidden state
+Y = 2; % size of observable state
+
+ns = [X Y];
+dnodes = [];
+onodes = [2];
+eclass1 = [1 2];
+eclass2 = [3 2];
+bnet = mk_dbn(intra, inter, ns, dnodes, eclass1, eclass2);
+
+x0 = rand(X,1);
+V0 = eye(X);
+C0 = rand(Y,X);
+R0 = eye(Y);
+A0 = rand(X,X);
+Q0 = eye(X);
+
+bnet.CPD{1} = gaussian_CPD(bnet, 1, 'mean', x0, 'cov', V0);
+%bnet.CPD{2} = gaussian_CPD(bnet, 2, 'mean', zeros(Y,1), 'cov', R0, 'weights', C0, 'full', 'untied', 'clamped_mean');
+%bnet.CPD{3} = gaussian_CPD(bnet, 3, 'mean', zeros(X,1), 'cov', Q0, 'weights', A0, 'full', 'untied', 'clamped_mean');
+bnet.CPD{2} = gaussian_CPD(bnet, 2, 'mean', zeros(Y,1), 'cov', R0, 'weights', C0);
+bnet.CPD{3} = gaussian_CPD(bnet, 3, 'mean', zeros(X,1), 'cov', Q0, 'weights', A0);
+
+
+T = 5; % fixed length sequences
+
+clear engine;
+%engine{1} = kalman_inf_engine(bnet, onodes);
+engine{1} = scg_unrolled_dbn_inf_engine(bnet, T, onodes);
+engine{2} = jtree_unrolled_dbn_inf_engine(bnet, T);
+
+N = length(engine);
+
+% inference
+
+ev = sample_dbn(bnet, T);
+evidence = cell(n,T);
+evidence(onodes,:) = ev(onodes, :);
+
+t = 2;
+query = [1 3];
+m = cell(1, N);
+ll = zeros(1, N);
+
+engine{1} = enter_evidence(engine{1}, evidence);
+[engine{2}, ll(2)] = enter_evidence(engine{2}, evidence);
+m{1} = marginal_nodes(engine{1}, query);
+m{2} = marginal_nodes(engine{2}, query, t);
+
+
+% compare all engines to engine{1}
+for i=2:N
+  assert(approxeq(m{1}.mu, m{i}.mu));
+  assert(approxeq(m{1}.Sigma, m{i}.Sigma));
+%  assert(approxeq(ll(1), ll(i)));
+end
+