git HEAD
--------
+* Added DECMAP algorithm.
+* Added createFactorDelta( const VarSet& vs, size_t state )
* [Peter Rockett] Improved Makefiles for image segmentation example.
* Removed deprecated interfaces.
* Added FactorGraph::isMaximal(size_t) and FactorGraph::maximalFactor(size_t).
WITHFLAGS:=$(WITHFLAGS) -DDAI_WITH_CBP
NAMES:=$(NAMES) bbp cbp bp_dual
endif
+ifdef WITH_DECMAP
+ WITHFLAGS:=$(WITHFLAGS) -DDAI_WITH_DECMAP
+ NAMES:=$(NAMES) decmap
+endif
# Define standard libDAI header dependencies, source file names and object file names
HEADERS=$(foreach name,graph dag bipgraph index var factor varset smallset prob daialg properties alldai enum exceptions util,$(INC)/$(name).h)
emalg$(OE) : $(SRC)/emalg.cpp $(INC)/emalg.h $(INC)/evidence.h $(HEADERS)
$(CC) -c $<
+decmap$(OE) : $(SRC)/decmap.cpp $(INC)/decmap.h $(HEADERS)
+ $(CC) -c $<
+
# EXAMPLES
###########
WITH_MR=true
WITH_GIBBS=true
WITH_CBP=true
+WITH_DECMAP=true
# Build doxygen documentation?
WITH_DOC=
* Double-loop GBP [HAK03];
* Various variants of Loop Corrected Belief Propagation [MoK07, MoR05];
* Gibbs sampler;
- * Conditioned Belief Propagation [EaG09].
+ * Conditioned Belief Propagation [EaG09];
+ * Decimation algorithm.
These inference methods can be used to calculate partition sums, marginals over
subsets of variables, and MAP states (the joint state of variables that has
// using the parameters specified by opts and two additional properties,
// specifying the type of updates the BP algorithm should perform and
// whether they should be done in the real or in the logdomain
- BP bp(fg, opts("updates",string("SEQFIX"))("logdomain",false));
+ BP bp(fg, opts("updates",string("SEQRND"))("logdomain",false));
// Initialize belief propagation algorithm
bp.init();
// Run belief propagation algorithm
//
// Note that inference is set to MAXPROD, which means that the object
// will perform the max-product algorithm instead of the sum-product algorithm
- BP mp(fg, opts("updates",string("SEQFIX"))("logdomain",false)("inference",string("MAXPROD")));
+ BP mp(fg, opts("updates",string("SEQRND"))("logdomain",false)("inference",string("MAXPROD"))("damping",string("0.1")));
// Initialize max-product algorithm
mp.init();
// Run max-product algorithm
// based on the max-product result
vector<size_t> mpstate = mp.findMaximum();
+ // Construct a decimation algorithm object from the FactorGraph fg
+ // using the parameters specified by opts and three additional properties,
+ // specifying that the decimation algorithm should use the max-product
+ // algorithm and should completely reinitalize its state at every step
+ DecMAP decmap(fg, opts("reinit",true)("ianame",string("BP"))("iaopts",string("[damping=0.1,inference=MAXPROD,logdomain=0,maxiter=1000,tol=1e-9,updates=SEQRND,verbose=1]")) );
+ decmap.init();
+ decmap.run();
+ vector<size_t> decmapstate = decmap.findMaximum();
+
// Report variable marginals for fg, calculated by the junction tree algorithm
cout << "Exact variable marginals:" << endl;
for( size_t i = 0; i < fg.nrVars(); i++ ) // iterate over all variables in fg
cout << "Approximate (max-product) MAP state (log probability = " << evalState( fg, mpstate ) << "):" << endl;
for( size_t i = 0; i < mpstate.size(); i++ )
cout << fg.var(i) << ": " << mpstate[i] << endl;
+
+ // Report DecMAP joint state
+ cout << "Approximate DecMAP state (log probability = " << evalState( fg, decmapstate ) << "):" << endl;
+ for( size_t i = 0; i < decmapstate.size(); i++ )
+ cout << fg.var(i) << ": " << decmapstate[i] << endl;
}
return 0;
#ifdef DAI_WITH_CBP
#include <dai/cbp.h>
#endif
+#ifdef DAI_WITH_DECMAP
+ #include <dai/decmap.h>
+#endif
/// Namespace for libDAI
#ifdef DAI_WITH_CBP
CBP::Name,
#endif
+#ifdef DAI_WITH_DECMAP
+ DecMAP::Name,
+#endif
""
};
*/
virtual Real logZ() const = 0;
+ /// Calculates the joint state of all variables that has maximum probability
+ /** \note Before this method is called, run() should have been called.
+ * \throw NOT_IMPLEMENTED if not implemented/supported
+ */
+ virtual std::vector<std::size_t> findMaximum() const { DAI_THROW(NOT_IMPLEMENTED); }
+
/// Returns maximum difference between single variable beliefs in the last iteration.
/** \throw NOT_IMPLEMENTED if not implemented/supported
*/
- virtual Real maxDiff() const = 0;
+ virtual Real maxDiff() const { DAI_THROW(NOT_IMPLEMENTED); };
/// Returns number of iterations done (one iteration passes over the complete factorgraph).
/** \throw NOT_IMPLEMENTED if not implemented/supported
*/
- virtual size_t Iterations() const = 0;
+ virtual size_t Iterations() const { DAI_THROW(NOT_IMPLEMENTED); };
//@}
/// \name Changing the factor graph
--- /dev/null
+/* 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]
+ */
+
+
+/// \file
+/// \brief Defines class DecMAP, which constructs a MAP state by decimation
+
+
+#ifndef __defined_libdai_decmap_h
+#define __defined_libdai_decmap_h
+
+
+#include <dai/daialg.h>
+
+
+namespace dai {
+
+
+/// Approximate inference algorithm DecMAP, which constructs a MAP state by decimation
+/** Decimation involves repeating the following two steps until no free variables remain:
+ * - run an approximate inference algorithm,
+ * - clamp the factor with the lowest entropy to its most probable state
+ */
+class DecMAP : public DAIAlgFG {
+ private:
+ /// Stores the final MAP state
+ std::vector<size_t> _state;
+ /// Stores the log probability of the MAP state
+ Real _logp;
+ /// Maximum difference encountered so far
+ Real _maxdiff;
+ /// Number of iterations needed
+ size_t _iters;
+
+ public:
+ /// Parameters for DecMAP
+ struct Properties {
+ /// Verbosity (amount of output sent to stderr)
+ size_t verbose;
+
+ /// Complete or partial reinitialization of clamped subgraphs?
+ bool reinit;
+
+ /// Name of the algorithm used to calculate the beliefs on clamped subgraphs
+ std::string ianame;
+
+ /// Parameters for the algorithm used to calculate the beliefs on clamped subgraphs
+ PropertySet iaopts;
+ } props;
+
+ /// Name of this inference algorithm
+ static const char *Name;
+
+ public:
+ /// Default constructor
+ DecMAP() : DAIAlgFG(), _state(), _logp(), _maxdiff(), _iters(), props() {}
+
+ /// Construct from FactorGraph \a fg and PropertySet \a opts
+ /** \param fg Factor graph.
+ * \param opts Parameters @see Properties
+ */
+ DecMAP( const FactorGraph &fg, const PropertySet &opts );
+
+
+ /// \name General InfAlg interface
+ //@{
+ virtual DecMAP* clone() const { return new DecMAP(*this); }
+ virtual std::string identify() const;
+ virtual Factor belief( const Var &v ) const { return beliefV( findVar( v ) ); }
+ virtual Factor belief( const VarSet &/*vs*/ ) const;
+ virtual Factor beliefV( size_t i ) const;
+ virtual Factor beliefF( size_t I ) const { return belief( factor(I).vars() ); }
+ virtual std::vector<Factor> beliefs() const;
+ virtual Real logZ() const { return _logp; }
+ virtual void init() { _maxdiff = 0.0; _iters = 0; }
+ virtual void init( const VarSet &/*ns*/ ) { init(); }
+ virtual Real run();
+ virtual std::vector<size_t> findMaximum() const { return _state; }
+ virtual Real maxDiff() const { return _maxdiff; }
+ virtual size_t Iterations() const { return _iters; }
+ virtual void setProperties( const PropertySet &opts );
+ virtual PropertySet getProperties() const;
+ virtual std::string printProperties() const;
+ //@}
+};
+
+
+} // end of namespace dai
+
+
+#endif
* - Various variants of Loop Corrected Belief Propagation
* [\ref MoK07, \ref MoR05];
* - Gibbs sampler;
- * - Conditioned Belief Propagation [\ref EaG09].
+ * - Conditioned Belief Propagation [\ref EaG09];
+ * - Decimation algorithm.
*
* These inference methods can be used to calculate partition sums, marginals
* over subsets of variables, and MAP states (the joint state of variables that
* - Loop Corrected Belief Propagation: dai::MR [\ref MoR05] and dai::LC [\ref MoK07]
* - Gibbs sampling: dai::Gibbs
* - Conditioned Belief Propagation: dai::CBP [\ref EaG09]
+ * - Decimation algorithm: dai::DECMAP
*
* Not all inference tasks are implemented by each method: calculating MAP states
- * is only possible with dai::JTree and dai::BP, calculating partition sums is
+ * is only possible with dai::JTree, dai::BP and dai::DECMAP; calculating partition sums is
* not possible with dai::MR, dai::LC and dai::Gibbs.
*
* \section terminology-learning Parameter learning
Factor createFactorDelta( const Var &v, size_t state );
+/// Returns a Kronecker delta point mass
+/** \param vs Set of variables
+ * \param state The state of \a vs that should get value 1
+ */
+Factor createFactorDelta( const VarSet& vs, size_t state );
+
+
} // end of namespace dai
if( name == CBP::Name )
return new CBP (fg, opts);
#endif
+#ifdef DAI_WITH_DECMAP
+ if( name == DecMAP::Name )
+ return new DecMAP (fg, opts);
+#endif
DAI_THROWE(UNKNOWN_DAI_ALGORITHM,"Unknown libDAI algorithm: " + name);
}
--- /dev/null
+/* 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>
+
+
+namespace dai {
+
+
+using namespace std;
+
+
+const char *DecMAP::Name = "DECMAP";
+
+
+void DecMAP::setProperties( const PropertySet &opts ) {
+ DAI_ASSERT( opts.hasKey("ianame") );
+ DAI_ASSERT( opts.hasKey("iaopts") );
+
+ props.ianame = opts.getStringAs<string>("ianame");
+ props.iaopts = opts.getStringAs<PropertySet>("iaopts");
+ if( opts.hasKey("verbose") )
+ props.verbose = opts.getStringAs<size_t>("verbose");
+ else
+ props.verbose = 0;
+ if( opts.hasKey("reinit") )
+ props.reinit = opts.getStringAs<bool>("reinit");
+ else
+ props.reinit = true;
+}
+
+
+PropertySet DecMAP::getProperties() const {
+ PropertySet opts;
+ opts.set( "verbose", props.verbose );
+ opts.set( "reinit", props.reinit );
+ opts.set( "ianame", props.ianame );
+ opts.set( "iaopts", props.iaopts );
+ return opts;
+}
+
+
+string DecMAP::printProperties() const {
+ stringstream s( stringstream::out );
+ s << "[";
+ s << "verbose=" << props.verbose << ",";
+ s << "reinit=" << props.reinit << ",";
+ s << "ianame=" << props.ianame << ",";
+ s << "iaopts=" << props.iaopts << "]";
+ return s.str();
+}
+
+
+DecMAP::DecMAP( const FactorGraph& fg, const PropertySet& opts ) : DAIAlgFG(fg), _state(), _logp(), _maxdiff(), _iters(), props() {
+ setProperties( opts );
+
+ _state = vector<size_t>( nrVars(), 0 );
+ _logp = -INFINITY;
+}
+
+
+string DecMAP::identify() const {
+ return string(Name) + printProperties();
+}
+
+
+Factor DecMAP::belief( const VarSet& vs ) const {
+ if( vs.size() == 0 )
+ return Factor();
+ else {
+ map<Var, size_t> state;
+ for( VarSet::const_iterator v = vs.begin(); v != vs.end(); v++ )
+ state[*v] = _state[findVar(*v)];
+ return createFactorDelta( vs, calcLinearState( vs, state ) );
+ }
+}
+
+
+Factor DecMAP::beliefV( size_t i ) const {
+ return createFactorDelta( var(i), _state[i] );
+}
+
+
+vector<Factor> DecMAP::beliefs() const {
+ vector<Factor> result;
+ for( size_t i = 0; i < nrVars(); ++i )
+ result.push_back( beliefV(i) );
+ for( size_t I = 0; I < nrFactors(); ++I )
+ result.push_back( beliefF(I) );
+ return result;
+}
+
+
+Real DecMAP::run() {
+ if( props.verbose >= 1 )
+ cerr << "Starting " << identify() << "...";
+ if( props.verbose >= 2 )
+ cerr << endl;
+
+ // the variables which have not been clamped yet
+ SmallSet<size_t> freeVars;
+ for( size_t i = 0; i < nrVars(); i++ )
+ freeVars |= i;
+
+ // prepare the inference algorithm object
+ InfAlg *clamped = newInfAlg( props.ianame, fg(), props.iaopts );
+
+ // decimate until no free variables remain
+ while( freeVars.size() ) {
+ Real md = clamped->run();
+ if( md > _maxdiff )
+ _maxdiff = md;
+ _iters += clamped->Iterations();
+
+ // store the variables that need initialization
+ VarSet varsToInit;
+ SmallSet<size_t> varsToClamp;
+
+ // schedule clamping for the free variables with zero entropy
+ for( SmallSet<size_t>::const_iterator it = freeVars.begin(); it != freeVars.end(); ) {
+ if( clamped->beliefV( *it ).entropy() == 0.0 ) {
+ // this variable should be clamped
+ varsToInit |= var( *it );
+ varsToClamp |= *it;
+ _state[*it] = clamped->beliefV( *it ).p().argmax().first;
+ freeVars.erase( *it );
+ } else
+ it++;
+ }
+
+ // find the free factor with lowest entropy
+ size_t bestI = 0;
+ Real bestEnt = INFINITY;
+ for( size_t I = 0; I < nrFactors(); I++ ) {
+ // check if the factor is still free
+ if( freeVars.intersects( bipGraph().nb2Set(I) ) ) {
+ Real EntI = clamped->beliefF(I).entropy();
+ if( EntI < bestEnt ) {
+ bestI = I;
+ bestEnt = EntI;
+ }
+ }
+ }
+
+ // schedule clamping for the factor with lowest entropy
+ vector<size_t> Istate(1,0);
+ Istate[0] = clamped->beliefF(bestI).p().argmax().first;
+ map<Var, size_t> Istatemap = calcState( factor(bestI).vars(), Istate[0] );
+ foreach( size_t i, bipGraph().nb2Set(bestI) & freeVars ) {
+ varsToInit |= var(i);
+ varsToClamp |= i;
+ _state[i] = Istatemap[var(i)];
+ freeVars.erase(i);
+ }
+
+ // clamp all variables scheduled for clamping
+ foreach( size_t i, varsToClamp )
+ clamped->clamp( i, _state[i], false );
+
+ // initialize clamped for the next run
+ if( props.reinit )
+ clamped->init();
+ else
+ clamped->init( varsToInit );
+ }
+
+ // calculate MAP state
+ map<Var, size_t> state;
+ for( size_t i = 0; i < nrVars(); i++ )
+ state[var(i)] = _state[i];
+ _logp = 0.0;
+ for( size_t I = 0; I < nrFactors(); I++ )
+ _logp += dai::log( factor(I)[calcLinearState( factor(I).vars(), state )] );
+
+ // clean up
+ delete clamped;
+
+ return _maxdiff;
+}
+
+
+} // end of namespace dai
}
+Factor createFactorDelta( const VarSet& vs, size_t state ) {
+ Factor fac( vs, 0.0 );
+ DAI_ASSERT( state < vs.nrStates() );
+ fac.set( state, 1.0 );
+ return fac;
+}
+
+
} // end of namespace dai
projects / libdai.git / commitdiff