include/dai/emalg.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) 2009  Charles Vaske  [cvaske at soe dot ucsc dot edu]
   8  *  Copyright (C) 2009  University of California, Santa Cruz
   9  */
  10
  11
  12 #ifndef __defined_libdai_emalg_h
  13 #define __defined_libdai_emalg_h
  14
  15
  16 #include <vector>
  17 #include <map>
  18
  19 #include <dai/factor.h>
  20 #include <dai/daialg.h>
  21 #include <dai/evidence.h>
  22 #include <dai/index.h>
  23 #include <dai/properties.h>
  24
  25
  26 /// \file
  27 /// \brief Defines classes related to Expectation Maximization: EMAlg, ParameterEstimation, CondProbEstimation and SharedParameters
  28 /// \todo Describe EM file format
  29
  30
  31 namespace dai {
  32
  33
  34 /// Interface for a parameter estimation method.
  35 /** This parameter estimation interface is based on sufficient statistics.
  36  *  Implementations are responsible for collecting data from a probability
  37  *  vector passed to it from a SharedParameters container object.
  38  *
  39  *  Implementations of this interface should register a factory function
  40  *  via the static ParameterEstimation::registerMethod function.
  41  *  The default registry only contains CondProbEstimation, named
  42  *  "ConditionalProbEstimation".
  43  */
  44 class ParameterEstimation {
  45     public:
  46         /// A pointer to a factory function.
  47         typedef ParameterEstimation* (*ParamEstFactory)( const PropertySet& );
  48
  49         /// Virtual destructor for deleting pointers to derived classes.
  50         virtual ~ParameterEstimation() {}
  51         /// Virtual copy constructor.
  52         virtual ParameterEstimation* clone() const = 0;
  53
  54         /// General factory method for construction of ParameterEstimation subclasses.
  55         static ParameterEstimation* construct( const std::string &method, const PropertySet &p );
  56
  57         /// Register a subclass so that it can be used with ParameterEstimation::construct.
  58         static void registerMethod( const std::string &method, const ParamEstFactory &f ) {
  59             if( _registry == NULL )
  60                 loadDefaultRegistry();
  61             (*_registry)[method] = f;
  62         }
  63
  64         /// Estimate the factor using the accumulated sufficient statistics and reset.
  65         virtual Prob estimate() = 0;
  66
  67         /// Accumulate the sufficient statistics for p.
  68         virtual void addSufficientStatistics( const Prob &p ) = 0;
  69
  70         /// Returns the size of the Prob that should be passed to addSufficientStatistics.
  71         virtual size_t probSize() const = 0;
  72
  73     private:
  74         /// A static registry containing all methods registered so far.
  75         static std::map<std::string, ParamEstFactory> *_registry;
  76
  77         /// Registers default ParameterEstimation subclasses (currently, only CondProbEstimation).
  78         static void loadDefaultRegistry();
  79 };
  80
  81
  82 /// Estimates the parameters of a conditional probability table, using pseudocounts.
  83 class CondProbEstimation : private ParameterEstimation {
  84     private:
  85         /// Number of states of the variable of interest
  86         size_t _target_dim;
  87         /// Current pseudocounts
  88         Prob _stats;
  89         /// Initial pseudocounts
  90         Prob _initial_stats;
  91
  92     public:
  93         /// Constructor
  94         /** For a conditional probability \f$ P( X | Y ) \f$,
  95          *  \param target_dimension should equal \f$ | X | \f$
  96          *  \param pseudocounts has length \f$ |X| \cdot |Y| \f$
  97          */
  98         CondProbEstimation( size_t target_dimension, const Prob &pseudocounts );
  99
 100         /// Virtual constructor, using a PropertySet.
 101         /** Some keys in the PropertySet are required.
 102          *  For a conditional probability \f$ P( X | Y ) \f$,
 103          *     - target_dimension should be equal to \f$ | X | \f$
 104          *     - total_dimension should be equal to \f$ |X| \cdot |Y| \f$
 105          *
 106          *  An optional key is:
 107          *     - pseudo_count, which specifies the initial counts (defaults to 1)
 108          */
 109         static ParameterEstimation* factory( const PropertySet &p );
 110
 111         /// Virtual copy constructor
 112         virtual ParameterEstimation* clone() const { return new CondProbEstimation( _target_dim, _initial_stats ); }
 113
 114         /// Virtual destructor
 115         virtual ~CondProbEstimation() {}
 116
 117         /// Returns an estimate of the conditional probability distribution.
 118         /** The format of the resulting Prob keeps all the values for
 119          *  \f$ P(X | Y=y) \f$ in sequential order in the array.
 120          */
 121         virtual Prob estimate();
 122
 123         /// Accumulate sufficient statistics from the expectations in p.
 124         virtual void addSufficientStatistics( const Prob &p );
 125
 126         /// Returns the required size for arguments to addSufficientStatistics.
 127         virtual size_t probSize() const { return _stats.size(); }
 128 };
 129
 130
 131 /// A single factor or set of factors whose parameters should be estimated.
 132 /** To ensure that parameters can be shared between different factors during
 133  *  EM learning, each factor's values are reordered to match a desired variable
 134  *  ordering. The ordering of the variables in a factor may therefore differ
 135  *  from the canonical ordering used in libDAI. The SharedParameters
 136  *  class couples one or more factors (together with the specified orderings
 137  *  of the variables) with a ParameterEstimation object, taking care of the
 138  *  necessary permutations of the factor entries / parameters.
 139  */
 140 class SharedParameters {
 141     public:
 142         /// Convenience label for an index into a factor in a FactorGraph.
 143         typedef size_t FactorIndex;
 144         /// Convenience label for a grouping of factor orientations.
 145         typedef std::map<FactorIndex, std::vector<Var> > FactorOrientations;
 146
 147     private:
 148         /// Maps factor indices to the corresponding VarSets
 149         std::map<FactorIndex, VarSet> _varsets;
 150         /// Maps factor indices to the corresponding Permute objects that permute the desired ordering into the canonical ordering
 151         std::map<FactorIndex, Permute> _perms;
 152         /// Maps factor indices to the corresponding desired variable orderings
 153         FactorOrientations _varorders;
 154         /// Parameter estimation method to be used
 155         ParameterEstimation *_estimation;
 156         /// Indicates whether the object pointed to by _estimation should be deleted upon destruction
 157         bool _deleteEstimation;
 158
 159         /// Calculates the permutation that permutes the variables in varorder into the canonical ordering
 160         static Permute calculatePermutation( const std::vector<Var> &varorder, VarSet &outVS );
 161
 162         /// Initializes _varsets and _perms from _varorders
 163         void setPermsAndVarSetsFromVarOrders();
 164
 165     public:
 166         /// Copy constructor
 167         SharedParameters( const SharedParameters &sp );
 168
 169         /// Constructor
 170         /** \param varorders  all the factor orientations for this parameter
 171          *  \param estimation a pointer to the parameter estimation method
 172          *  \param deletePE whether the parameter estimation object should be deleted in the destructor
 173          */
 174         SharedParameters( const FactorOrientations &varorders, ParameterEstimation *estimation, bool deletePE=false );
 175
 176         /// Constructor for making an object from a stream and a factor graph
 177         SharedParameters( std::istream &is, const FactorGraph &fg_varlookup );
 178
 179         /// Destructor
 180         ~SharedParameters() {
 181             if( _deleteEstimation )
 182                 delete _estimation;
 183         }
 184
 185         /// Collect the necessary statistics from expected values
 186         void collectSufficientStatistics( InfAlg &alg );
 187
 188         /// Estimate and set the shared parameters
 189         void setParameters( FactorGraph &fg );
 190
 191         /// Returns the parameters
 192         void collectParameters( const FactorGraph &fg, std::vector<Real> &outVals, std::vector<Var> &outVarOrder );
 193 };
 194
 195
 196 /// A MaximizationStep groups together several parameter estimation tasks into a single unit.
 197 class MaximizationStep {
 198     private:
 199         std::vector<SharedParameters> _params;
 200
 201     public:
 202         /// Default constructor
 203         MaximizationStep() : _params() {}
 204
 205         /// Constructor from a vector of SharedParameters objects
 206         MaximizationStep( std::vector<SharedParameters> &maximizations ) : _params(maximizations) {}
 207
 208         /// Constructor from an input stream and a corresponding factor graph
 209         MaximizationStep( std::istream &is, const FactorGraph &fg_varlookup );
 210
 211         /// Collect the beliefs from this InfAlg as expectations for the next Maximization step.
 212         void addExpectations( InfAlg &alg );
 213
 214         /// Using all of the currently added expectations, make new factors with maximized parameters and set them in the FactorGraph.
 215         void maximize( FactorGraph &fg );
 216
 217         /// @name Iterator interface
 218         //@{
 219         typedef std::vector<SharedParameters>::iterator iterator;
 220         typedef std::vector<SharedParameters>::const_iterator const_iterator;
 221         iterator begin() { return _params.begin(); }
 222         const_iterator begin() const { return _params.begin(); }
 223         iterator end() { return _params.end(); }
 224         const_iterator end() const { return _params.end(); }
 225         //@}
 226 };
 227
 228
 229 /// EMAlg performs Expectation Maximization to learn factor parameters.
 230 /** This requires specifying:
 231  *     - Evidence (instances of observations from the graphical model),
 232  *     - InfAlg for performing the E-step, which includes the factor graph,
 233  *     - a vector of MaximizationSteps steps to be performed.
 234  *
 235  *  This implementation can perform incremental EM by using multiple
 236  *  MaximizationSteps.  An expectation step is performed between execution
 237  *  of each MaximizationStep.  A call to iterate() will cycle through all
 238  *  MaximizationSteps.
 239  *
 240  *  Having multiple and separate maximization steps allows for maximizing some
 241  *  parameters, performing another E step, and then maximizing separate
 242  *  parameters, which may result in faster convergence in some cases.
 243  */
 244 class EMAlg {
 245     private:
 246         /// All the data samples used during learning
 247         const Evidence &_evidence;
 248
 249         /// How to do the expectation step
 250         InfAlg &_estep;
 251
 252         /// The maximization steps to take
 253         std::vector<MaximizationStep> _msteps;
 254
 255         /// Number of iterations done
 256         size_t _iters;
 257
 258         /// History of likelihoods
 259         std::vector<Real> _lastLogZ;
 260
 261         /// Maximum number of iterations
 262         size_t _max_iters;
 263
 264         /// Convergence tolerance
 265         Real _log_z_tol;
 266
 267     public:
 268         /// Key for setting maximum iterations @see setTermConditions
 269         static const std::string MAX_ITERS_KEY;
 270         /// Default maximum iterations @see setTermConditions
 271         static const size_t MAX_ITERS_DEFAULT;
 272         /// Key for setting likelihood termination condition @see setTermConditions
 273         static const std::string LOG_Z_TOL_KEY;
 274         /// Default likelihood tolerance @see setTermConditions
 275         static const Real LOG_Z_TOL_DEFAULT;
 276
 277         /// Construct an EMAlg from all these objects
 278         EMAlg( const Evidence &evidence, InfAlg &estep, std::vector<MaximizationStep> &msteps, const PropertySet &termconditions )
 279           : _evidence(evidence), _estep(estep), _msteps(msteps), _iters(0), _lastLogZ(), _max_iters(MAX_ITERS_DEFAULT), _log_z_tol(LOG_Z_TOL_DEFAULT)
 280         {
 281               setTermConditions( termconditions );
 282         }
 283
 284         /// Construct an EMAlg from an Evidence object, an InfAlg object, and an input stream
 285         EMAlg( const Evidence &evidence, InfAlg &estep, std::istream &mstep_file );
 286
 287         /// Change the conditions for termination
 288         /** There are two possible parameters in the PropertySet
 289          *    - max_iters maximum number of iterations
 290          *    - log_z_tol proportion of increase in logZ
 291          *
 292          *  \see hasSatisifiedTermConditions()
 293          */
 294         void setTermConditions( const PropertySet &p );
 295
 296         /// Determine if the termination conditions have been met.
 297         /** There are two sufficient termination conditions:
 298          *    -# the maximum number of iterations has been performed
 299          *    -# the ratio of logZ increase over previous logZ is less than the
 300          *       tolerance, i.e.,
 301          *       \f$ \frac{\log(Z_t) - \log(Z_{t-1})}{| \log(Z_{t-1}) | } < \mathrm{tol} \f$.
 302          */
 303         bool hasSatisfiedTermConditions() const;
 304
 305         /// Return the last calculated log likelihood
 306         Real getLogZ() const { return _lastLogZ.back(); }
 307
 308         /// Returns number of iterations done so far
 309         size_t getCurrentIters() const { return _iters; }
 310
 311         /// Get the iteration method used
 312         const InfAlg& eStep() const { return _estep; }
 313
 314         /// Perform an iteration over all maximization steps
 315         Real iterate();
 316
 317         /// Perform an iteration over a single MaximizationStep
 318         Real iterate( MaximizationStep &mstep );
 319
 320         /// Iterate until termination conditions are satisfied
 321         void run();
 322
 323         /// @name Iterator interface
 324         //@{
 325         typedef std::vector<MaximizationStep>::iterator s_iterator;
 326         typedef std::vector<MaximizationStep>::const_iterator const_s_iterator;
 327         s_iterator s_begin() { return _msteps.begin(); }
 328         const_s_iterator s_begin() const { return _msteps.begin(); }
 329         s_iterator s_end() { return _msteps.end(); }
 330         const_s_iterator s_end() const { return _msteps.end(); }
 331         //@}
 332 };
 333
 334
 335 } // end of namespace dai
 336
 337
 338 #endif