examples/example_sprinkler_gibbs.cpp

   1 /*  This file is part of libDAI - http://www.libdai.org/
   2  *
   3  *  Copyright (c) 2006-2011, The libDAI authors. All rights reserved.
   4  *
   5  *  Use of this source code is governed by a BSD-style license that can be found in the LICENSE file.
   6  */
   7
   8
   9 #include <dai/factorgraph.h>
  10 #include <dai/gibbs.h>
  11 #include <dai/properties.h>
  12 #include <iostream>
  13 #include <fstream>
  14
  15 using namespace std;
  16 using namespace dai;
  17
  18 int main() {
  19     // This example program illustrates how to use Gibbs sampling
  20     // to sample from a joint probability distribution described
  21     // by a factor graph, using the sprinkler network example discussed at
  22     // http://www.cs.ubc.ca/~murphyk/Bayes/bnintro.html
  23     //
  24     // The file sprinkler.fg has to be generated by first running
  25     // example_sprinkler
  26
  27     // Read the factorgraph from the file
  28     FactorGraph SprinklerNetwork;
  29     SprinklerNetwork.ReadFromFile( "sprinkler.fg" );
  30     cout << "Sprinkler network read from sprinkler.fg" << endl;
  31
  32     // Output some information about the factorgraph
  33     cout << SprinklerNetwork.nrVars() << " variables" << endl;
  34     cout << SprinklerNetwork.nrFactors() << " factors" << endl;
  35
  36     // Prepare a Gibbs sampler
  37     PropertySet gibbsProps;
  38     gibbsProps.set("maxiter", size_t(100));   // number of Gibbs sampler iterations
  39     gibbsProps.set("burnin", size_t(0));
  40     gibbsProps.set("verbose", size_t(0));
  41     Gibbs gibbsSampler( SprinklerNetwork, gibbsProps );
  42
  43     // Open a .tab file for writing
  44     ofstream outfile;
  45     outfile.open( "sprinkler.tab" );
  46     if( !outfile.is_open() )
  47         throw "Cannot write to file!";
  48
  49     // Write header (consisting of variable labels)
  50     for( size_t i = 0; i < SprinklerNetwork.nrVars(); i++ )
  51         outfile << (i == 0 ? "" : "\t") << SprinklerNetwork.var(i).label();
  52     outfile << endl << endl;
  53
  54     // Draw samples from joint distribution using Gibbs sampling
  55     // and write them to the .tab file
  56     size_t nrSamples = 1000;
  57     std::vector<size_t> state;
  58     for( size_t t = 0; t < nrSamples; t++ ) {
  59         gibbsSampler.init();
  60         gibbsSampler.run();
  61         state = gibbsSampler.state();
  62         for( size_t i = 0; i < state.size(); i++ )
  63             outfile << (i == 0 ? "" : "\t") << state[i];
  64         outfile << endl;
  65     }
  66     cout << nrSamples << " samples written to sprinkler.tab" << endl;
  67
  68     // Close the .tab file
  69     outfile.close();
  70
  71     return 0;
  72 }