Added examples example_sprinkler_gibbs and example_sprinkler_em
authorJoris Mooij <joris.mooij@tuebingen.mpg.de>
Tue, 12 Jan 2010 16:39:33 +0000 (17:39 +0100)
committerJoris Mooij <joris.mooij@tuebingen.mpg.de>
Tue, 12 Jan 2010 16:39:33 +0000 (17:39 +0100)
13 files changed:
ChangeLog
Makefile
examples/example_sprinkler.cpp
examples/example_sprinkler_em.cpp [new file with mode: 0644]
examples/example_sprinkler_gibbs.cpp [new file with mode: 0644]
examples/sprinkler.em [new file with mode: 0644]
include/dai/doc.h
include/dai/emalg.h
include/dai/factorgraph.h
include/dai/gibbs.h
src/gibbs.cpp
swig/example_sprinkler.m
swig/example_sprinkler.py

index 38505d7..90076a8 100644 (file)
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,4 @@
+* Added examples example_sprinkler_gibbs and example_sprinkler_em
 * Optimized BipartiteGraph::isConnected() by using Boost Graph Library implementation
 * Strengthened convergence criteria of various algorithms
 * Implemented FBP::logZ()
index 049ab92..d2064f6 100644 (file)
--- a/Makefile
+++ b/Makefile
@@ -13,7 +13,7 @@ include Makefile.conf
 
 # Set version and date
 DAI_VERSION="git HEAD"
-DAI_DATE="January 11, 2010 - or later"
+DAI_DATE="January 12, 2010 - or later"
 
 # Directories of libDAI sources
 # Location libDAI headers
@@ -111,7 +111,7 @@ MEX:=$(MEX) $(CCLIB) $(CCINC) $(MEXFLAGS)
 
 all : $(TARGETS)
 
-examples : examples/example$(EE) examples/example_bipgraph$(EE) examples/example_varset$(EE) examples/example_permute$(EE) examples/example_sprinkler$(EE)
+examples : examples/example$(EE) examples/example_bipgraph$(EE) examples/example_varset$(EE) examples/example_permute$(EE) examples/example_sprinkler$(EE) examples/example_sprinkler_gibbs$(EE) examples/example_sprinkler_em$(EE)
 
 matlabs : matlab/dai$(ME) matlab/dai_readfg$(ME) matlab/dai_writefg$(ME) matlab/dai_potstrength$(ME)
 
@@ -222,6 +222,12 @@ examples/example_permute$(EE) : examples/example_permute.cpp $(HEADERS) $(LIB)/l
 examples/example_sprinkler$(EE) : examples/example_sprinkler.cpp $(HEADERS) $(LIB)/libdai$(LE)
        $(CC) $(CCO)examples/example_sprinkler$(EE) examples/example_sprinkler.cpp $(LIBS)
 
+examples/example_sprinkler_gibbs$(EE) : examples/example_sprinkler_gibbs.cpp $(HEADERS) $(LIB)/libdai$(LE)
+       $(CC) $(CCO)examples/example_sprinkler_gibbs$(EE) examples/example_sprinkler_gibbs.cpp $(LIBS)
+
+examples/example_sprinkler_em$(EE) : examples/example_sprinkler_em.cpp $(HEADERS) $(LIB)/libdai$(LE)
+       $(CC) $(CCO)examples/example_sprinkler_em$(EE) examples/example_sprinkler_em.cpp $(LIBS)
+
 
 # TESTS
 ########
@@ -321,7 +327,7 @@ ifneq ($(OS),WINDOWS)
 clean :
        -rm *$(OE)
        -rm matlab/*$(ME)
-       -rm examples/example$(EE) examples/example_bipgraph$(EE) examples/example_varset$(EE) examples/example_permute$(EE) examples/example_sprinkler$(EE)
+       -rm examples/example$(EE) examples/example_bipgraph$(EE) examples/example_varset$(EE) examples/example_permute$(EE) examples/example_sprinkler$(EE) examples/example_sprinkler_gibbs$(EE) examples/example_sprinkler_em$(EE)
        -rm tests/testdai$(EE) tests/testem/testem$(EE) tests/testbbp$(EE)
        -rm utils/fg2dot$(EE) utils/createfg$(EE) utils/fginfo$(EE)
        -rm -R doc
index 0d1e570..74a7441 100644 (file)
@@ -4,12 +4,13 @@
  *  2, or (at your option) any later version. libDAI is distributed without any
  *  warranty. See the file COPYING for more details.
  *
- *  Copyright (C) 2008-2009  Joris Mooij  [joris dot mooij at libdai dot org]
+ *  Copyright (C) 2008-2010  Joris Mooij  [joris dot mooij at libdai dot org]
  */
 
 
 #include <dai/factorgraph.h>
 #include <iostream>
+#include <fstream>
 
 using namespace std;
 using namespace dai;
@@ -64,9 +65,7 @@ int main() {
 
     // Write factorgraph to a file
     SprinklerNetwork.WriteToFile( "sprinkler.fg" );
-
-    // Reread the factorgraph from the file
-    SprinklerNetwork.ReadFromFile( "sprinkler.fg" );
+    cout << "Sprinkler network written to sprinkler.fg" << endl;
 
     // Output some information about the factorgraph
     cout << SprinklerNetwork.nrVars() << " variables" << endl;
diff --git a/examples/example_sprinkler_em.cpp b/examples/example_sprinkler_em.cpp
new file mode 100644 (file)
index 0000000..fee9f19
--- /dev/null
@@ -0,0 +1,69 @@
+/*  This file is part of libDAI - http://www.libdai.org/
+ *
+ *  libDAI is licensed under the terms of the GNU General Public License version
+ *  2, or (at your option) any later version. libDAI is distributed without any
+ *  warranty. See the file COPYING for more details.
+ *
+ *  Copyright (C) 2010  Joris Mooij  [joris dot mooij at libdai dot org]
+ */
+
+#include <dai/alldai.h>
+#include <iostream>
+#include <fstream>
+#include <string>
+
+using namespace std;
+using namespace dai;
+
+int main() {
+    // This example program illustrates how to learn the
+    // parameters of a Bayesian network from a sample of
+    // the sprinkler network discussed at
+    // http://www.cs.ubc.ca/~murphyk/Bayes/bnintro.html
+    //
+    // The factor graph file (sprinkler.fg) has to be generated first
+    // by running example_sprinkler, and the data sample file 
+    // (sprinkler.tab) by running example_sprinkler_gibbs
+
+    // Read the factorgraph from the file
+    FactorGraph SprinklerNetwork;
+    SprinklerNetwork.ReadFromFile( "sprinkler.fg" );
+
+    // Prepare junction-tree object for doing exact inference for E-step
+    PropertySet infprops;
+    infprops.Set( "verbose", (size_t)1 );
+    infprops.Set( "updates", string("HUGIN") );
+    InfAlg* inf = newInfAlg( "JTREE", SprinklerNetwork, infprops );
+    inf->init();
+
+    // Read sample from file
+    Evidence e;
+    ifstream estream( "sprinkler.tab" );
+    e.addEvidenceTabFile( estream, SprinklerNetwork );
+    cout << "Number of samples: " << e.nrSamples() << endl;
+
+    // Read EM specification
+    ifstream emstream( "sprinkler.em" );
+    EMAlg em(e, *inf, emstream);
+
+    // Iterate EM until convergence
+    while( !em.hasSatisfiedTermConditions() ) {
+        Real l = em.iterate();
+        cout << "Iteration " << em.Iterations() << " likelihood: " << l <<endl;
+    }
+
+    // Output true factor graph
+    cout << endl << "True factor graph:" << endl << "##################" << endl;
+    cout.precision(12);
+    cout << SprinklerNetwork;
+
+    // Output learned factor graph
+    cout << endl << "Learned factor graph:" << endl << "#####################" << endl;
+    cout.precision(12);
+    cout << inf->fg();
+
+    // Clean up
+    delete inf;
+
+    return 0;
+}
diff --git a/examples/example_sprinkler_gibbs.cpp b/examples/example_sprinkler_gibbs.cpp
new file mode 100644 (file)
index 0000000..8d90b4e
--- /dev/null
@@ -0,0 +1,73 @@
+/*  This file is part of libDAI - http://www.libdai.org/
+ *
+ *  libDAI is licensed under the terms of the GNU General Public License version
+ *  2, or (at your option) any later version. libDAI is distributed without any
+ *  warranty. See the file COPYING for more details.
+ *
+ *  Copyright (C) 2010  Joris Mooij  [joris dot mooij at libdai dot org]
+ */
+
+
+#include <dai/factorgraph.h>
+#include <dai/gibbs.h>
+#include <dai/properties.h>
+#include <iostream>
+#include <fstream>
+
+using namespace std;
+using namespace dai;
+
+int main() {
+    // This example program illustrates how to use Gibbs sampling
+    // to sample from a joint probability distribution described
+    // by a factor graph, using the sprinkler network example discussed at
+    // http://www.cs.ubc.ca/~murphyk/Bayes/bnintro.html
+    //
+    // The file sprinkler.fg has to be generated by first running
+    // example_sprinkler
+
+    // Read the factorgraph from the file
+    FactorGraph SprinklerNetwork;
+    SprinklerNetwork.ReadFromFile( "sprinkler.fg" );
+    cout << "Sprinkler network read from sprinkler.fg" << endl;
+
+    // Output some information about the factorgraph
+    cout << SprinklerNetwork.nrVars() << " variables" << endl;
+    cout << SprinklerNetwork.nrFactors() << " factors" << endl;
+
+    // Prepare a Gibbs sampler
+    PropertySet gibbsProps;
+    gibbsProps.Set("iters", size_t(100));   // number of Gibbs sampler iterations
+    gibbsProps.Set("burnin", size_t(0));
+    gibbsProps.Set("verbose", size_t(0));
+    Gibbs gibbsSampler( SprinklerNetwork, gibbsProps );
+
+    // Open a .tab file for writing
+    ofstream outfile;
+    outfile.open( "sprinkler.tab" );
+    if( !outfile.is_open() )
+        throw "Cannot write to file!";
+
+    // Write header (consisting of variable labels)
+    for( size_t i = 0; i < SprinklerNetwork.nrVars(); i++ )
+        outfile << (i == 0 ? "" : "\t") << SprinklerNetwork.var(i).label();
+    outfile << endl << endl;
+
+    // Draw samples from joint distribution using Gibbs sampling
+    // and write them to the .tab file
+    size_t nrSamples = 1000;
+    std::vector<size_t> state;
+    for( size_t t = 0; t < nrSamples; t++ ) {
+        gibbsSampler.run();
+        state = gibbsSampler.state();
+        for( size_t i = 0; i < state.size(); i++ )
+            outfile << (i == 0 ? "" : "\t") << state[i];
+        outfile << endl;
+    }
+    cout << nrSamples << " samples written to sprinkler.tab" << endl;
+
+    // Close the .tab file
+    outfile.close();
+
+    return 0;
+}
diff --git a/examples/sprinkler.em b/examples/sprinkler.em
new file mode 100644 (file)
index 0000000..c1cdd20
--- /dev/null
@@ -0,0 +1,15 @@
+1
+
+4
+CondProbEstimation [target_dim=2,total_dim=2,pseudo_count=1]
+1
+0 0
+CondProbEstimation [target_dim=2,total_dim=4,pseudo_count=1]
+1
+1 2 0
+CondProbEstimation [target_dim=2,total_dim=4,pseudo_count=1]
+1
+2 1 0
+CondProbEstimation [target_dim=2,total_dim=8,pseudo_count=1]
+1
+3 3 1 2
index 57b115d..b15aa7c 100644 (file)
  *  In libDAI, both types of graphical models are represented by a slightly more 
  *  general type of graphical model: a factor graph [\ref KFL01].
  *
+ *  An example of a Bayesian network is: 
  *  \dot
  *  digraph bayesnet {
  *    size="1,1";
  *    x2 -> x4;
  *  }
  *  \enddot
- *
+ *  The probability distribution of a Bayesian network factorizes as:
  *  \f[ P(\mathbf{x}) = \prod_{i\in\mathcal{V}} P(x_i \,|\, x_{\mathrm{pa}(i)}) \f]
  *  where \f$\mathrm{pa}(i)\f$ are the parents of node \a i in a DAG.
  *
+ *  The same probability distribution can be represented as a Markov random field:
  *  \dot
  *  graph mrf {
  *    size="1.5,1.5";
  *  }
  *  \enddot
  *
+ *  The probability distribution of a Markov random field factorizes as:
  *  \f[ P(\mathbf{x}) = \frac{1}{Z} \prod_{C\in\mathcal{C}} \psi_C(x_C) \f]
  *  where \f$ \mathcal{C} \f$ are the cliques of an undirected graph, 
  *  \f$ \psi_C(x_C) \f$ are "potentials" or "compatibility functions", and
  *  \f$ Z \f$ is the partition sum which properly normalizes the probability
  *  distribution.
  *
+ *  Finally, the same probability distribution can be represented as a factor graph:
  *  \dot
  *  graph factorgraph {
  *    size="1.8,1";
  *  }
  *  \enddot
  *
+ *  The probability distribution of a factor graph factorizes as:
  *  \f[ P(\mathbf{x}) = \frac{1}{Z} \prod_{I\in \mathcal{F}} f_I(x_I) \f]
  *  where \f$ \mathcal{F} \f$ are the factor nodes of a factor graph (a 
  *  bipartite graph consisting of variable nodes and factor nodes), 
  *  contain the variable labels of the variables on which that factor depends, 
  *  in a specific ordering. This ordering can be different from the canonical 
  *  ordering of the variables used internally in libDAI (which would be sorted 
- *  ascendingly according to the variable labels). The odering of the variables
+ *  ascendingly according to the variable labels). The ordering of the variables
  *  specifies the implicit ordering of the shared parameters: when iterating
  *  over all shared parameters, the corresponding index of the first variable
  *  changes fastest (in the inner loop), and the corresponding index of the
index 9093c78..d9379da 100644 (file)
@@ -293,9 +293,6 @@ class MaximizationStep {
  *  parameters, which may result in faster convergence in some cases.
  *
  *  \author Charles Vaske
- *
- *  \todo Expand the sprinkler example such that it generates a sample of the probability distribution,
- *  and add another example that estimates the factor parameters of the sprinkler network using EM from this sample.
  */
 class EMAlg {
     private:
@@ -407,4 +404,9 @@ class EMAlg {
 } // end of namespace dai
 
 
+/** \example example_sprinkler_em.cpp
+ *  This example shows how to use the EMAlg class.
+ */
+
+
 #endif
index 36bede3..a8b07c5 100644 (file)
@@ -370,6 +370,12 @@ FactorGraph::FactorGraph(FactorInputIterator fact_begin, FactorInputIterator fac
  */
 
 
+/** \example example_sprinkler.cpp
+ *  This example illustrates how to manually construct a factor graph and
+ *  write it to a file.
+ */
+
+
 } // end of namespace dai
 
 
index ceb4f63..0135397 100644 (file)
@@ -121,10 +121,15 @@ class Gibbs : public DAIAlgFG {
 /// Runs Gibbs sampling for \a iters iterations (of which \a burnin for burn-in) on FactorGraph \a fg, and returns the resulting state
 /** \relates Gibbs
  */
-std::vector<size_t> getGibbsState( const FactorGraph &fg, size_t iters, size_t burnin=0 );
+std::vector<size_t> getGibbsState( const FactorGraph &fg, size_t iters );
 
 
 } // end of namespace dai
 
 
+/** \example example_sprinkler_gibbs.cpp
+ *  This example shows how to use the Gibbs class.
+ */
+
+
 #endif
index a8c8645..b62cdf3 100644 (file)
@@ -232,10 +232,10 @@ Factor Gibbs::belief( const VarSet &ns ) const {
 }
 
 
-std::vector<size_t> getGibbsState( const FactorGraph &fg, size_t iters, size_t burnin ) {
+std::vector<size_t> getGibbsState( const FactorGraph &fg, size_t iters ) {
     PropertySet gibbsProps;
     gibbsProps.Set("iters", iters);
-    gibbsProps.Set("burnin", burnin);
+    gibbsProps.Set("burnin", size_t(0));
     gibbsProps.Set("verbose", size_t(0));
     Gibbs gibbs( fg, gibbsProps );
     gibbs.run();
index b714b89..c7b83de 100644 (file)
@@ -61,9 +61,7 @@ SprinklerNetwork = dai.FactorGraph(SprinklerFactors);
 
 % Write factorgraph to a file
 SprinklerNetwork.WriteToFile('sprinkler.fg');
-
-% Reread the factorgraph from the file
-SprinklerNetwork.ReadFromFile('sprinkler.fg');
+fprintf('Sprinkler network written to sprinkler.fg');
 
 % Output some information about the factorgraph
 fprintf('%d variables\n', SprinklerNetwork.nrVars());
index a94a795..b1d2d5d 100644 (file)
@@ -62,9 +62,7 @@ SprinklerNetwork = dai.FactorGraph(SprinklerFactors)
 
 # Write factorgraph to a file
 SprinklerNetwork.WriteToFile('sprinkler.fg')
-
-# Reread the factorgraph from the file
-SprinklerNetwork.ReadFromFile('sprinkler.fg')
+print 'Sprinkler network written to sprinkler.fg'
 
 # Output some information about the factorgraph
 print SprinklerNetwork.nrVars(), 'variables'