Improved documentation of include/dai/cbp.h
authorJoris Mooij <joris.mooij@tuebingen.mpg.de>
Tue, 10 Nov 2009 17:43:47 +0000 (18:43 +0100)
committerJoris Mooij <joris.mooij@tuebingen.mpg.de>
Tue, 10 Nov 2009 17:43:47 +0000 (18:43 +0100)
TODO
include/dai/bbp.h
include/dai/cbp.h
include/dai/doc.h
src/bbp.cpp
src/cbp.cpp

diff --git a/TODO b/TODO
index be66603..223f083 100644 (file)
--- a/TODO
+++ b/TODO
@@ -2,7 +2,6 @@ To do for the next release (0.2.3):
 
 Improve documentation:
 
-       cbp.h
        evidence.h
        emalg.h
        doc.h
index c361afe..741ff2c 100644 (file)
@@ -258,7 +258,7 @@ class BBP {
     public:
     /// \name Constructors/destructors
     //@{
-        /// Construct from a InfAlg \a ia and a PropertySet \a opts
+        /// Construct BBP object from a InfAlg \a ia and a PropertySet \a opts
         BBP( const InfAlg *ia, const PropertySet &opts ) : _bp_dual(ia), _fg(&(ia->fg())), _ia(ia) {
             props.set(opts);
         }
index 370228d..e996506 100644 (file)
@@ -9,9 +9,7 @@
 
 
 /// \file
-/// \brief Defines class CBP [\ref EaG09]
-/// \author Frederik Eaton
-/// \todo Improve documentation
+/// \brief Defines class CBP, which implements Clamped Belief Propagation
 
 
 #ifndef __defined_libdai_cbp_h
 namespace dai {
 
 
-/// Find a variable to clamp using BBP (goes with maximum adjoint)
-/// \see BBP
-std::pair<size_t, size_t> bbpFindClampVar( const InfAlg &in_bp, bool clampingVar, const PropertySet &bbp_props, const BBPCostFunction &cfn, Real *maxVarOut );
-
-
-/// Class for CBP (Clamped Belief Propagation)
-/** This algorithm uses configurable heuristics to choose a variable
- *  x_i and a state x_i*. Inference is done with x_i "clamped" to x_i*
- *  (i.e., conditional on x_i == x_i*), and also with the negation of this
- *  condition. Clamping is done recursively up to a fixed number of
- *  levels (other stopping criteria are also implemented, see
- *  \a recursion property). The resulting approximate marginals are
- *  combined using logZ estimates.
+/// Class for CBP (Clamped Belief Propagation) [\ref EaG09]
+/** This approximate inference algorithm uses configurable heuristics to choose a variable
+ *  \f$ x_i \f$ and a state \f$ x_i^* \f$. Inference is done with \f$ x_i \f$ "clamped" to \f$ x_i^* \f$
+ *  (i.e., conditional on \f$ x_i = x_i^* \f$), and also with the negation of this condition. 
+ *  Clamping is done recursively up to a fixed number of levels (other stopping criteria are 
+ *  also implemented, see the CBP::Properties::RecurseType property). The resulting approximate 
+ *  marginals are combined using estimates of the partition sum.
  *
  *  \author Frederik Eaton
  */
@@ -50,61 +42,25 @@ class CBP : public DAIAlgFG {
         std::vector<Factor> _beliefsV;
         /// Factor beliefs
         std::vector<Factor> _beliefsF;
-        /// Log-partition sum
+        /// Logarithm of partition sum
         Real _logZ;
 
-        /// Counts number of clampings at each leaf node
-        Real _sum_level;
+        /// Numer of iterations needed
+        size_t _iters;
+        /// Maximum difference encountered so far
+        Real _maxdiff;
 
+        /// Number of clampings at each leaf node
+        Real _sum_level;
         /// Number of leaves of recursion tree
         size_t _num_leaves;
 
         /// Output stream where information about the clampings is written
         boost::shared_ptr<std::ofstream> _clamp_ofstream;
 
-        /// Returns BBP cost function used
-        BBPCostFunction BBP_cost_function() { return props.bbp_cfn; }
-
-        /// Prints beliefs, variables and partition sum, in case of a debugging build
-        void printDebugInfo();
-
-        /// Called by 'run', and by itself. Implements the main algorithm.
-        /** Chooses a variable to clamp, recurses, combines the logZ and
-         *  beliefs estimates of the children, and returns the improved
-         *  estimates in \a lz_out and \a beliefs_out to its parent
-         */
-        void runRecurse( InfAlg *bp, Real orig_logZ, std::vector<size_t> clamped_vars_list, size_t &num_leaves,
-                         size_t &choose_count, Real &sum_level, Real &lz_out, std::vector<Factor> &beliefs_out );
-
-        /// Choose the next variable to clamp
-        /** Choose the next variable to clamp, given a converged InfAlg (\a bp),
-         *  and a vector of variables that are already clamped (\a
-         *  clamped_vars_list). Returns the chosen variable in \a i, and
-         *  the set of states in \a xis. If \a maxVarOut is non-NULL and
-         *  props.choose==CHOOSE_BBP then it is used to store the
-         *  adjoint of the chosen variable
-         */
-        virtual bool chooseNextClampVar( InfAlg* bp, std::vector<size_t> &clamped_vars_list, size_t &i, std::vector<size_t> &xis, Real *maxVarOut );
-
-        /// Return the InfAlg to use at each step of the recursion.
-        /// \todo At present, only returns a BP instance
-        InfAlg* getInfAlg();
-
-        /// Numer of iterations needed
-        size_t _iters;
-        /// Maximum difference encountered so far
-        Real _maxdiff;
-
-        /// Sets variable beliefs, factor beliefs and logZ
-        /** \param bs should be a concatenation of the variable beliefs followed by the factor beliefs
-         */
-        void setBeliefs( const std::vector<Factor> &bs, Real logZ );
-
-        /// Constructor helper function
-        void construct();
 
     public:
-        /// Construct CBP object from FactorGraph fg and PropertySet opts
+        /// Construct CBP object from FactorGraph \a fg and PropertySet \a opts
         CBP( const FactorGraph &fg, const PropertySet &opts ) : DAIAlgFG(fg) {
             props.set( opts );
             construct();
@@ -229,32 +185,55 @@ class CBP : public DAIAlgFG {
         } props;
 /* }}} END OF GENERATED CODE */
 
-        /// Returns heuristic used for clamping variable
-        Properties::ChooseMethodType ChooseMethod() { return props.choose; }
-        /// Returns method used for deciding when to stop recursing
-        Properties::RecurseType Recursion() { return props.recursion; }
-        /// Returns clamping type used
-        Properties::ClampType Clamping() { return props.clamp; }
-        /// Returns maximum number of levels of recursion
-        size_t maxClampLevel() { return props.max_levels; }
-        /// Returns props.min_max_adj @see CBP::Properties::min_max_adj
-        Real minMaxAdj() { return props.min_max_adj; }
-        /// Returns tolerance used for controlling recursion depth
-        Real recTol() { return props.rec_tol; }
-};
+    private:
+        /// Prints beliefs, variables and partition sum, in case of a debugging build
+        void printDebugInfo();
 
+        /// Called by run(), and by itself. Implements the main algorithm.
+        /** Chooses a variable to clamp, recurses, combines the partition sum 
+         *  and belief estimates of the children, and returns the improved
+         *  estimates in \a lz_out and \a beliefs_out to its parent.
+         */
+        void runRecurse( InfAlg *bp, Real orig_logZ, std::vector<size_t> clamped_vars_list, size_t &num_leaves,
+                         size_t &choose_count, Real &sum_level, Real &lz_out, std::vector<Factor> &beliefs_out );
 
-/// Given a sorted vector of states \a xis and total state count \a n_states, return a vector of states not in \a xis
-std::vector<size_t> complement( std::vector<size_t>& xis, size_t n_states );
+        /// Choose the next variable to clamp.
+        /** Choose the next variable to clamp, given a converged InfAlg \a bp,
+         *  and a vector of variables that are already clamped (\a
+         *  clamped_vars_list). Returns the chosen variable in \a i, and
+         *  the set of states in \a xis. If \a maxVarOut is non-NULL and
+         *  \a props.choose == \c CHOOSE_BBP then it is used to store the
+         *  adjoint of the chosen variable.
+         */
+        virtual bool chooseNextClampVar( InfAlg* bp, std::vector<size_t> &clamped_vars_list, size_t &i, std::vector<size_t> &xis, Real *maxVarOut );
+
+        /// Return the InfAlg to use at each step of the recursion.
+        /// \todo At present, CBP::getInfAlg() only returns a BP instance
+        InfAlg* getInfAlg();
+
+        /// Sets variable beliefs, factor beliefs and log partition sum
+        /** \param bs should be a concatenation of the variable beliefs followed by the factor beliefs
+         */
+        void setBeliefs( const std::vector<Factor> &bs, Real logZ );
 
-/// Computes \f$\frac{\exp(a)}{\exp(a)+\exp(b)}\f$
-Real unSoftMax( Real a, Real b );
+        /// Constructor helper function
+        void construct();
+};
 
-/// Computes log of sum of exponents, i.e., \f$\log\left(\exp(a) + \exp(b)\right)\f$
-Real logSumExp( Real a, Real b );
 
-/// Compute sum of pairwise L-infinity distances of the first \a nv factors in each vector
-Real dist( const std::vector<Factor>& b1, const std::vector<Factor>& b2, size_t nv );
+/// Find the best variable/factor to clamp using BBP.
+/** Takes a converged inference algorithm as input, runs Gibbs and BP_dual, creates
+ *  and runs a BBP object, finds the best variable/factor (the one with the maximum
+ *  factor adjoint), and returns the corresponding (index,state) pair.
+ *  \param in_bp inference algorithm (compatible with BP) that should have converged;
+ *  \param clampingVar if \c true, finds best variable, otherwise, finds best factor;
+ *  \param bbp_props BBP parameters to use;
+ *  \param cfn BBP cost function to use;
+ *  \param maxVarOut maximum adjoint value (only set if not NULL).
+ *  \see BBP
+ *  \relates CBP
+ */
+std::pair<size_t, size_t> BBPFindClampVar( const InfAlg &in_bp, bool clampingVar, const PropertySet &bbp_props, const BBPCostFunction &cfn, Real *maxVarOut );
 
 
 } // end of namespace dai
index 9b5fb34..78a4f30 100644 (file)
@@ -13,7 +13,6 @@
  *
  *  \todo Improve documentation
  *
- *  \todo Merge COPYING into doxygen documentation
  *  \todo Merge README into doxygen documentation
  *  \todo Document examples, tests and utils
  *
index 1055c35..2d329dc 100644 (file)
@@ -986,7 +986,7 @@ void BBP::run() {
 }
 
 
-Real numericBBPTest( const InfAlg &bp, const vector<size_t> *state, const PropertySet &bbp_props, const BBPCostFunction &cfn, Real h ) {
+Real numericBBPTest( const InfAlg &bp, const std::vector<size_t> *state, const PropertySet &bbp_props, const BBPCostFunction &cfn, Real h ) {
     BBP bbp( &bp, bbp_props );
     // calculate the value of the unperturbed cost function
     Real cf0 = cfn.evaluate( bp, state );
index d77a90f..39bd6ab 100644 (file)
@@ -32,6 +32,48 @@ using boost::shared_ptr;
 const char *CBP::Name = "CBP";
 
 
+/// Given a sorted vector of states \a xis and total state count \a n_states, return a vector of states not in \a xis
+vector<size_t> complement( vector<size_t> &xis, size_t n_states ) {
+    vector<size_t> cmp_xis( 0 );
+    size_t j = 0;
+    for( size_t xi = 0; xi < n_states; xi++ ) {
+        while( j < xis.size() && xis[j] < xi )
+            j++;
+        if( j >= xis.size() || xis[j] > xi )
+            cmp_xis.push_back(xi);
+    }
+    DAI_ASSERT( xis.size()+cmp_xis.size() == n_states );
+    return cmp_xis;
+}
+
+
+/// Computes \f$\frac{\exp(a)}{\exp(a)+\exp(b)}\f$
+Real unSoftMax( Real a, Real b ) {
+    if( a > b )
+        return 1.0 / (1.0 + exp(b-a));
+    else
+        return exp(a-b) / (exp(a-b) + 1.0);
+}
+
+
+/// Computes log of sum of exponents, i.e., \f$\log\left(\exp(a) + \exp(b)\right)\f$
+Real logSumExp( Real a, Real b ) {
+    if( a > b )
+        return a + log1p( exp( b-a ) );
+    else
+        return b + log1p( exp( a-b ) );
+}
+
+
+/// Compute sum of pairwise L-infinity distances of the first \a nv factors in each vector
+Real dist( const vector<Factor> &b1, const vector<Factor> &b2, size_t nv ) {
+    Real d = 0.0;
+    for( size_t k = 0; k < nv; k++ )
+        d += dist( b1[k], b2[k], Prob::DISTLINF );
+    return d;
+}
+
+
 void CBP::setBeliefs( const std::vector<Factor> &bs, Real logZ ) {
     size_t i = 0;
     _beliefsV.clear();
@@ -131,9 +173,9 @@ void CBP::runRecurse( InfAlg *bp, Real orig_logZ, vector<size_t> clamped_vars_li
     vector<size_t> xis;
     Real maxVar = 0.0;
     bool found;
-    bool clampingVar = (Clamping() == Properties::ClampType::CLAMP_VAR);
+    bool clampingVar = (props.clamp == Properties::ClampType::CLAMP_VAR);
 
-    if( Recursion() == Properties::RecurseType::REC_LOGZ && recTol() > 0 && exp( bp->logZ() - orig_logZ ) < recTol() )
+    if( props.recursion == Properties::RecurseType::REC_LOGZ && props.rec_tol > 0 && exp( bp->logZ() - orig_logZ ) < props.rec_tol )
         found = false;
     else
         found = chooseNextClampVar( bp, clamped_vars_list, i, xis, &maxVar );
@@ -198,11 +240,11 @@ void CBP::runRecurse( InfAlg *bp, Real orig_logZ, vector<size_t> clamped_vars_li
     Real p = unSoftMax( lz, cmp_lz );
     Real bp__d = 0.0;
 
-    if( Recursion() == Properties::RecurseType::REC_BDIFF && recTol() > 0 ) {
+    if( props.recursion == Properties::RecurseType::REC_BDIFF && props.rec_tol > 0 ) {
         vector<Factor> combined_b( mixBeliefs( p, b, cmp_b ) );
         Real new_lz = logSumExp( lz,cmp_lz );
         bp__d = dist( bp->beliefs(), combined_b, nrVars() );
-        if( exp( new_lz - orig_logZ) * bp__d < recTol() ) {
+        if( exp( new_lz - orig_logZ) * bp__d < props.rec_tol ) {
             num_leaves++;
             sum_level += clamped_vars_list.size();
             beliefs_out = combined_b;
@@ -224,7 +266,7 @@ void CBP::runRecurse( InfAlg *bp, Real orig_logZ, vector<size_t> clamped_vars_li
     if( props.verbose >= 2 ) {
         Real d = dist( bp->beliefs(), beliefs_out, nrVars() );
         cerr << "Distance (clamping " << i << "): " << d;
-        if( Recursion() == Properties::RecurseType::REC_BDIFF )
+        if( props.recursion == Properties::RecurseType::REC_BDIFF )
             cerr << "; bp_dual predicted " << bp__d;
         cerr << "; max_adjoint = " << maxVar << "; logZ = " << lz_out << " (in " << bp->logZ() << ") (orig " << orig_logZ << "); p = " << p << "; level = " << clamped_vars_list.size() << endl;
     }
@@ -239,10 +281,10 @@ bool CBP::chooseNextClampVar( InfAlg *bp, vector<size_t> &clamped_vars_list, siz
     Real tiny = 1.0e-14;
     if( props.verbose >= 3 )
         cerr << "clamped_vars_list" << clamped_vars_list << endl;
-    if( clamped_vars_list.size() >= maxClampLevel() )
+    if( clamped_vars_list.size() >= props.max_levels )
         return false;
-    if( ChooseMethod() == Properties::ChooseMethodType::CHOOSE_RANDOM ) {
-        if( Clamping() == Properties::ClampType::CLAMP_VAR ) {
+    if( props.choose == Properties::ChooseMethodType::CHOOSE_RANDOM ) {
+        if( props.clamp == Properties::ClampType::CLAMP_VAR ) {
             int t = 0, t1 = 100;
             do {
                 i = rnd( nrVars() );
@@ -282,8 +324,8 @@ bool CBP::chooseNextClampVar( InfAlg *bp, vector<size_t> &clamped_vars_list, siz
             // DAI_ASSERT(!_clamped_vars.count(i)); // not true for >2-ary variables
             DAI_IFVERB(2, endl<<"CHOOSE_RANDOM chose factor "<<i<<" state "<<xis[0]<<endl);
         }
-    } else if( ChooseMethod() == Properties::ChooseMethodType::CHOOSE_MAXENT ) {
-        if( Clamping() == Properties::ClampType::CLAMP_VAR ) {
+    } else if( props.choose == Properties::ChooseMethodType::CHOOSE_MAXENT ) {
+        if( props.clamp == Properties::ClampType::CLAMP_VAR ) {
             Real max_ent = -1.0;
             int win_k = -1, win_xk = -1;
             for( size_t k = 0; k < nrVars(); k++ ) {
@@ -324,14 +366,14 @@ bool CBP::chooseNextClampVar( InfAlg *bp, vector<size_t> &clamped_vars_list, siz
                 return false;
             }
         }
-    } else if( ChooseMethod()==Properties::ChooseMethodType::CHOOSE_BP_L1 ||
-               ChooseMethod()==Properties::ChooseMethodType::CHOOSE_BP_CFN ) {
-        bool doL1 = (ChooseMethod() == Properties::ChooseMethodType::CHOOSE_BP_L1);
+    } else if( props.choose==Properties::ChooseMethodType::CHOOSE_BP_L1 ||
+               props.choose==Properties::ChooseMethodType::CHOOSE_BP_CFN ) {
+        bool doL1 = (props.choose == Properties::ChooseMethodType::CHOOSE_BP_L1);
         vector<size_t> state;
-        if( !doL1 && BBP_cost_function().needGibbsState() )
+        if( !doL1 && props.bbp_cfn.needGibbsState() )
             state = getGibbsState( bp->fg(), 2*bp->Iterations() );
         // try clamping each variable manually
-        DAI_ASSERT( Clamping() == Properties::ClampType::CLAMP_VAR );
+        DAI_ASSERT( props.clamp == Properties::ClampType::CLAMP_VAR );
         Real max_cost = 0.0;
         int win_k = -1, win_xk = -1;
         for( size_t k = 0; k < nrVars(); k++ ) {
@@ -347,7 +389,7 @@ bool CBP::chooseNextClampVar( InfAlg *bp, vector<size_t> &clamped_vars_list, siz
                     for( size_t j = 0; j < nrVars(); j++ )
                         cost += dist( bp->beliefV(j), bp1->beliefV(j), Prob::DISTL1 );
                 else
-                    cost = BBP_cost_function().evaluate( *bp1, &state );
+                    cost = props.bbp_cfn.evaluate( *bp1, &state );
                 if( cost > max_cost || win_k == -1 ) {
                     max_cost = cost;
                     win_k = k;
@@ -360,15 +402,15 @@ bool CBP::chooseNextClampVar( InfAlg *bp, vector<size_t> &clamped_vars_list, siz
         DAI_ASSERT( win_xk >= 0 );
         i = win_k;
         xis.resize( 1, win_xk );
-    } else if( ChooseMethod() == Properties::ChooseMethodType::CHOOSE_BBP ) {
+    } else if( props.choose == Properties::ChooseMethodType::CHOOSE_BBP ) {
         Real mvo;
         if( !maxVarOut )
             maxVarOut = &mvo;
-        bool clampingVar = (Clamping() == Properties::ClampType::CLAMP_VAR);
-        pair<size_t, size_t> cv = bbpFindClampVar( *bp, clampingVar, props.bbp_props, BBP_cost_function(), &mvo );
+        bool clampingVar = (props.clamp == Properties::ClampType::CLAMP_VAR);
+        pair<size_t, size_t> cv = BBPFindClampVar( *bp, clampingVar, props.bbp_props, props.bbp_cfn, &mvo );
 
         // if slope isn't big enough then don't clamp
-        if( mvo < minMaxAdj() )
+        if( mvo < props.min_max_adj )
             return false;
 
         size_t xi = cv.second;
@@ -380,11 +422,11 @@ bool CBP::chooseNextClampVar( InfAlg *bp, vector<size_t> &clamped_vars_list, siz
             << ", maxVar = "<< mvo << ")"
         Prob b = ( clampingVar ? bp->beliefV(i).p() : bp->beliefF(i).p());
         if( b[xi] < tiny ) {
-            cerr << "Warning, at level " << clamped_vars_list.size() << ", bbpFindClampVar found unlikely " << VAR_INFO << endl;
+            cerr << "Warning, at level " << clamped_vars_list.size() << ", BBPFindClampVar found unlikely " << VAR_INFO << endl;
             return false;
         }
         if( abs(b[xi] - 1) < tiny ) {
-            cerr << "Warning at level " << clamped_vars_list.size() << ", bbpFindClampVar found overly likely " << VAR_INFO << endl;
+            cerr << "Warning at level " << clamped_vars_list.size() << ", BBPFindClampVar found overly likely " << VAR_INFO << endl;
             return false;
         }
 
@@ -408,14 +450,7 @@ void CBP::printDebugInfo() {
 }
 
 
-//----------------------------------------------------------------
-
-
-/** Function which takes a factor graph as input, runs Gibbs and BP_dual,
- *  creates and runs a BBP object, finds best variable, returns
- *  (variable,state) pair for clamping
- */
-pair<size_t, size_t> bbpFindClampVar( const InfAlg &in_bp, bool clampingVar, const PropertySet &bbp_props, const BBPCostFunction &cfn, Real *maxVarOut ) {
+pair<size_t, size_t> BBPFindClampVar( const InfAlg &in_bp, bool clampingVar, const PropertySet &bbp_props, const BBPCostFunction &cfn, Real *maxVarOut ) {
     BBP bbp( &in_bp, bbp_props );
     bbp.initCostFnAdj( cfn, NULL );
     bbp.run();
@@ -475,44 +510,6 @@ pair<size_t, size_t> bbpFindClampVar( const InfAlg &in_bp, bool clampingVar, con
 }
 
 
-vector<size_t> complement( vector<size_t> &xis, size_t n_states ) {
-    vector<size_t> cmp_xis( 0 );
-    size_t j = 0;
-    for( size_t xi = 0; xi < n_states; xi++ ) {
-        while( j < xis.size() && xis[j] < xi )
-            j++;
-        if( j >= xis.size() || xis[j] > xi )
-            cmp_xis.push_back(xi);
-    }
-    DAI_ASSERT( xis.size()+cmp_xis.size() == n_states );
-    return cmp_xis;
-}
-
-
-Real unSoftMax( Real a, Real b ) {
-    if( a > b )
-        return 1.0 / (1.0 + exp(b-a));
-    else
-        return exp(a-b) / (exp(a-b) + 1.0);
-}
-
-
-Real logSumExp( Real a, Real b ) {
-    if( a > b )
-        return a + log1p( exp( b-a ) );
-    else
-        return b + log1p( exp( a-b ) );
-}
-
-
-Real dist( const vector<Factor> &b1, const vector<Factor> &b2, size_t nv ) {
-    Real d = 0.0;
-    for( size_t k = 0; k < nv; k++ )
-        d += dist( b1[k], b2[k], Prob::DISTLINF );
-    return d;
-}
-
-
 } // end of namespace dai