include/dai/treeep.h

   1 /*  This file is part of libDAI - http://www.libdai.org/
   2  *
   3  *  libDAI is licensed under the terms of the GNU General Public License version
   4  *  2, or (at your option) any later version. libDAI is distributed without any
   5  *  warranty. See the file COPYING for more details.
   6  *
   7  *  Copyright (C) 2006-2009  Joris Mooij  [joris dot mooij at libdai dot org]
   8  *  Copyright (C) 2006-2007  Radboud University Nijmegen, The Netherlands
   9  */
  10
  11
  12 /// \file
  13 /// \brief Defines class TreeEP, which implements Tree Expectation Propagation
  14 /// \todo Clean up the TreeEP code
  15
  16
  17 #ifndef __defined_libdai_treeep_h
  18 #define __defined_libdai_treeep_h
  19
  20
  21 #include <vector>
  22 #include <string>
  23 #include <dai/daialg.h>
  24 #include <dai/varset.h>
  25 #include <dai/regiongraph.h>
  26 #include <dai/factorgraph.h>
  27 #include <dai/clustergraph.h>
  28 #include <dai/weightedgraph.h>
  29 #include <dai/jtree.h>
  30 #include <dai/properties.h>
  31 #include <dai/enum.h>
  32
  33
  34 namespace dai {
  35
  36
  37 /// Approximate inference algorithm "Tree Expectation Propagation" [\ref MiQ04]
  38 class TreeEP : public JTree {
  39     private:
  40         /// Maximum difference encountered so far
  41         Real                  _maxdiff;
  42         /// Number of iterations needed
  43         size_t                _iters;
  44
  45     public:
  46         /// Parameters for TreeEP
  47         struct Properties {
  48             /// Enumeration of possible choices for the tree
  49             /** The two possibilities are:
  50              *  - \c ORG: take the maximum spanning tree where the weights are crude
  51              *            estimates of the mutual information between the nodes;
  52              *  - \c ALT: take the maximum spanning tree where the weights are upper
  53              *            bounds on the effective interaction strengths between pairs of nodes.
  54              */
  55             DAI_ENUM(TypeType,ORG,ALT);
  56
  57             /// Verbosity (amount of output sent to stderr)
  58             size_t verbose;
  59
  60             /// Maximum number of iterations
  61             size_t maxiter;
  62
  63             /// Tolerance for convergence test
  64             Real tol;
  65
  66             /// How to choose the tree
  67             TypeType type;
  68         } props;
  69
  70         /// Name of this inference method
  71         static const char *Name;
  72
  73     private:
  74         /// Stores the data structures needed to efficiently update the approximation of an off-tree factor
  75         /// \todo Write documentation
  76         class TreeEPSubTree {
  77             private:
  78                 std::vector<Factor>  _Qa;
  79                 std::vector<Factor>  _Qb;
  80                 RootedTree           _RTree;
  81                 std::vector<size_t>  _a;        // _Qa[alpha]  <->  superTree.Qa[_a[alpha]]
  82                 std::vector<size_t>  _b;        // _Qb[beta]   <->  superTree.Qb[_b[beta]]
  83                                                 // _Qb[beta]   <->  _RTree[beta]
  84                 const Factor *       _I;
  85                 VarSet               _ns;
  86                 VarSet               _nsrem;
  87                 Real                 _logZ;
  88
  89
  90             public:
  91                 TreeEPSubTree() : _Qa(), _Qb(), _RTree(), _a(), _b(), _I(NULL), _ns(), _nsrem(), _logZ(0.0) {}
  92                 TreeEPSubTree( const TreeEPSubTree &x) : _Qa(x._Qa), _Qb(x._Qb), _RTree(x._RTree), _a(x._a), _b(x._b), _I(x._I), _ns(x._ns), _nsrem(x._nsrem), _logZ(x._logZ) {}
  93                 TreeEPSubTree & operator=( const TreeEPSubTree& x ) {
  94                     if( this != &x ) {
  95                         _Qa         = x._Qa;
  96                         _Qb         = x._Qb;
  97                         _RTree      = x._RTree;
  98                         _a          = x._a;
  99                         _b          = x._b;
 100                         _I          = x._I;
 101                         _ns         = x._ns;
 102                         _nsrem      = x._nsrem;
 103                         _logZ       = x._logZ;
 104                     }
 105                     return *this;
 106                 }
 107
 108                 TreeEPSubTree( const RootedTree &subRTree, const RootedTree &jt_RTree, const std::vector<Factor> &jt_Qa, const std::vector<Factor> &jt_Qb, const Factor *I );
 109                 void init();
 110                 void InvertAndMultiply( const std::vector<Factor> &Qa, const std::vector<Factor> &Qb );
 111                 void HUGIN_with_I( std::vector<Factor> &Qa, std::vector<Factor> &Qb );
 112                 Real logZ( const std::vector<Factor> &Qa, const std::vector<Factor> &Qb ) const;
 113                 const Factor *& I() { return _I; }
 114         };
 115
 116         /// Stores a TreeEPSubTree object for each off-tree factor
 117         std::map<size_t, TreeEPSubTree>  _Q;
 118
 119     public:
 120         /// Default constructor
 121         TreeEP() : JTree(), _maxdiff(0.0), _iters(0), props(), _Q() {}
 122
 123         /// Copy constructor
 124         TreeEP( const TreeEP &x ) : JTree(x), _maxdiff(x._maxdiff), _iters(x._iters), props(x.props), _Q(x._Q) {
 125             for( size_t I = 0; I < nrFactors(); I++ )
 126                 if( offtree( I ) )
 127                     _Q[I].I() = &factor(I);
 128         }
 129
 130         /// Assignment operator
 131         TreeEP& operator=( const TreeEP &x ) {
 132             if( this != &x ) {
 133                 JTree::operator=( x );
 134                 _maxdiff = x._maxdiff;
 135                 _iters   = x._iters;
 136                 props    = x.props;
 137                 _Q       = x._Q;
 138                 for( size_t I = 0; I < nrFactors(); I++ )
 139                     if( offtree( I ) )
 140                         _Q[I].I() = &factor(I);
 141             }
 142             return *this;
 143         }
 144
 145         /// Construct from FactorGraph \a fg and PropertySet \a opts
 146         /** \param opts Parameters @see Properties
 147          *  \note The factor graph has to be connected.
 148          *  \throw FACTORGRAPH_NOT_CONNECTED if \a fg is not connected
 149          */
 150         TreeEP( const FactorGraph &fg, const PropertySet &opts );
 151
 152
 153     /// \name General InfAlg interface
 154     //@{
 155         virtual TreeEP* clone() const { return new TreeEP(*this); }
 156         virtual std::string identify() const;
 157         virtual Real logZ() const;
 158         virtual void init();
 159         virtual void init( const VarSet &/*ns*/ ) { init(); }
 160         virtual Real run();
 161         virtual Real maxDiff() const { return _maxdiff; }
 162         virtual size_t Iterations() const { return _iters; }
 163         virtual void setProperties( const PropertySet &opts );
 164         virtual PropertySet getProperties() const;
 165         virtual std::string printProperties() const;
 166     //@}
 167
 168
 169     private:
 170         /// Helper function for constructors
 171         void construct( const RootedTree &tree );
 172         /// Returns \c true if factor \a I is not part of the tree
 173         bool offtree( size_t I ) const { return (fac2OR[I] == -1U); }
 174 };
 175
 176
 177 } // end of namespace dai
 178
 179
 180 #endif