standardmixturepdf.hpp

dysii是一款非常出色的滤波函数库

#ifndef INDII_ML_AUX_STANDARDMIXTUREPDF_HPP
#define INDII_ML_AUX_STANDARDMIXTUREPDF_HPP

#include "MixturePdf.hpp"

namespace indii {
  namespace ml {
    namespace aux {
/**
 * Mixture probability density with standard calculation of covariance.
 *
 * @author Lawrence Murray <[e-mail]>
 * @version $Rev: 404 $
 * @date $Date: 2008-03-05 14:52:55 +0000 (Wed, 05 Mar 2008) $
 *
 * @param T Component type, should be derivative of Pdf.
 *
 * @see MixturePdf for more information regarding the serialization
 * and parallelisation features of this class.
 */
template <class T>
class StandardMixturePdf : public MixturePdf<T> {
public:
  /**
   * Default constructor.
   *
   * Initialises the mixture with zero dimensions. This should
   * generally only be used when the object is to be restored from a
   * serialization. Indeed, there is no other way to resize the
   * mixture to nonzero dimensionality except by subsequently
   * restoring from a serialization.
   */
  StandardMixturePdf();

  /**
   * Constructor. One or more components should be added with
   * addComponent() after construction.
   *
   * @param N Dimensionality of the distribution.
   */
  StandardMixturePdf(const unsigned int N);

  /**
   * Destructor.
   */
  virtual ~StandardMixturePdf();

  virtual void setDimensions(const unsigned int N,
      const bool preserve = false);

  /**
   * Get the covariance of the distribution. The covariance is defined
   * as:
   *
   * \f[ \frac{1}{W}\sum_{i=1}^{K} w^{(i)}(\Sigma^{(i)} +
   * \mathbf{\mu}^{(i)} (\mathbf{\mu}^{(i)})^{T}) -
   * \mathbf{\mu}\mathbf{\mu}^{T} \f]
   *
   * where \f$\mathbf{\mu}^{(i)}\f$, \f$\Sigma^{(i)}\f$ and
   * \f$w^{(i)}\f$ are the mean, covariance and weight of the
   * \f$i\f$th component, respectively, and \f$\mathbf{\mu}\f$ the
   * overall mean.
   *
   * @return \f$\Sigma\f$; covariance of the distribution.
   */
  virtual const symmetric_matrix& getCovariance();

  /**
   * Get the covariance of the full distribution.
   *
   * @return \f$\Sigma\f$; covariance of the full distribution.
   */
  symmetric_matrix getDistributedCovariance();

protected:
  virtual void dirty();

private:
  /**
   * Last calculated covariance.
   */
  symmetric_matrix sigma;

  /**
   * Last calculated unnormalized covariance.
   */
  symmetric_matrix Zsigma;

  /**
   * Is the last calculated covariance up to date?
   */
  bool haveSigma;

  /**
   * Calculate covariance from current components.
   *
   * The covariance is calculated as:
   *
   * \f[\Sigma =
   * \frac{1}{W}\sum_{i=1}^{K}w_i(\Sigma_i+\mathbf{\mu}_i\mathbf{\mu}_i^T)
   * - \bar{\mathbf{\mu}}\bar{\mathbf{\mu}}^T
   * \f]
   *
   * where \f$\mathbf{\mu}_i\f$ and \f$\Sigma_i\f$ are the mean and
   * covariance of the \f$i\f$th component, respectively, and
   * \f$\bar{\mathbf{\mu}}\f$ is the mean of means.
   */
  void calculateCovariance();

  /**
   * Serialize.
   */
  template<class Archive>
  void save(Archive& ar, const unsigned int version) const;

  /**
   * Restore from serialization.
   */
  template<class Archive>
  void load(Archive& ar, const unsigned int version);

  /*
   * Boost.Serialization requirements.
   */
  BOOST_SERIALIZATION_SPLIT_MEMBER()
  friend class boost::serialization::access;

};

    }
  }
}

#include "boost/serialization/base_object.hpp"

namespace ublas = boost::numeric::ublas;

template <class T>
indii::ml::aux::StandardMixturePdf<T>::StandardMixturePdf() :
    indii::ml::aux::MixturePdf<T>(0), sigma(0), Zsigma(0) {
  haveSigma = false;
}

template <class T>
indii::ml::aux::StandardMixturePdf<T>::StandardMixturePdf(
    const unsigned int N) : indii::ml::aux::MixturePdf<T>(N), sigma(N),
    Zsigma(N) {
  haveSigma = false;
}

template <class T>
indii::ml::aux::StandardMixturePdf<T>::~StandardMixturePdf() {
  //
}

template <class T>
void indii::ml::aux::StandardMixturePdf<T>::setDimensions(
    const unsigned int N, const bool preserve) {
  MixturePdf<T>::setDimensions(N, preserve);

  Zsigma.resize(N, false);
  sigma.resize(N, false);
}

template <class T>
const indii::ml::aux::symmetric_matrix& indii::ml::aux::StandardMixturePdf<
    T>::getCovariance() {
  if (!haveSigma) {
    calculateCovariance();
  }
  return sigma;
}

template <class T>
indii::ml::aux::symmetric_matrix indii::ml::aux::StandardMixturePdf<
    T>::getDistributedCovariance() {
  const unsigned int N = this->getDimensions();

  boost::mpi::communicator world;
  matrix Zsigma_d(N,N), sigma_d(N,N);
  vector mu_d(N);

  if (this->getTotalWeight() > 0.0) {
    if (!haveSigma) {
      calculateCovariance();
    }
  } else {
    Zsigma.clear();
  }

  noalias(mu_d) = this->getDistributedExpectation();
  noalias(Zsigma_d) = boost::mpi::all_reduce(world, matrix(Zsigma),
      std::plus<matrix>());
  noalias(sigma_d) = Zsigma_d / this->getDistributedTotalWeight() -
      outer_prod(mu_d, mu_d);

  return ublas::symmetric_adaptor<matrix, ublas::lower>(sigma_d);
}

template <class T>
void indii::ml::aux::StandardMixturePdf<T>::dirty() {
  MixturePdf<T>::dirty();
  haveSigma = false;
}

template <class T>
void indii::ml::aux::StandardMixturePdf<T>::calculateCovariance() {
  /* pre-condition */
  assert (this->getTotalWeight() > 0.0);

  typename MixturePdf<T>::weighted_component_iterator iter, end;
  double w;
  const vector& mu = this->getExpectation();

  Zsigma.clear();
  iter = this->getComponents().begin();
  end = this->getComponents().end();
  while (iter != end) {
    w = iter->w;
    const vector& mu_i = iter->x.getExpectation();
    const symmetric_matrix& sigma_i = iter->x.getCovariance();

    noalias(Zsigma) += w * (sigma_i + outer_prod(mu_i, mu_i));
    iter++;
  }
  noalias(sigma) = Zsigma / this->getTotalWeight() - outer_prod(mu, mu);
  haveSigma = true;
}

template <class T>
template <class Archive>
void indii::ml::aux::StandardMixturePdf<T>::save(Archive& ar,
    const unsigned int version) const {
  ar & boost::serialization::base_object<
      indii::ml::aux::MixturePdf<T> >(*this);
}

template <class T>
template <class Archive>
void indii::ml::aux::StandardMixturePdf<T>::load(Archive& ar,
    const unsigned int version) {
  ar & boost::serialization::base_object<
      indii::ml::aux::MixturePdf<T> >(*this);

  sigma.resize(this->getDimensions(), false);
  Zsigma.resize(this->getDimensions(), false);
  haveSigma = false;
}

#endif