source: trunk/yat/classifier/NCC.h @ 1115

#ifndef _theplu_yat_classifier_ncc_ 
#define _theplu_yat_classifier_ncc_ 

// $Id$

/*
  Copyright (C) 2005 Markus Ringnér, Peter Johansson
  Copyright (C) 2006 Jari Häkkinen, Markus Ringnér, Peter Johansson
  Copyright (C) 2007 Peter Johansson

  This file is part of the yat library, http://trac.thep.lu.se/yat

  The yat library is free software; you can redistribute it and/or
  modify it under the terms of the GNU General Public License as
  published by the Free Software Foundation; either version 2 of the
  License, or (at your option) any later version.

  The yat library is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
  General Public License for more details.

  You should have received a copy of the GNU General Public License
  along with this program; if not, write to the Free Software
  Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA
  02111-1307, USA.
*/

#include "DataLookup1D.h"
#include "DataLookup2D.h"
#include "DataLookupWeighted1D.h"
#include "MatrixLookup.h"
#include "MatrixLookupWeighted.h"
#include "SupervisedClassifier.h"
#include "Target.h"

#include "yat/statistics/Averager.h"
#include "yat/statistics/AveragerWeighted.h"
#include "yat/utility/matrix.h"
#include "yat/utility/vector.h"
#include "yat/utility/stl_utility.h"
#include "yat/utility/yat_assert.h"

#include<iostream>
#include<iterator>
#include <map>
#include <cmath>
#include <stdexcept>

namespace theplu {
namespace yat {
namespace classifier {  


  ///
  /// @brief Class for Nearest Centroid Classification.
  ///
  /// The template argument Distance should be a class modelling
  /// the concept \ref concept_distance.
  ///
  template <typename Distance>
  class NCC : public SupervisedClassifier
  {
  
  public:
    ///
    /// Constructor taking the training data and the target vector as
    /// input
    ///
    NCC(const MatrixLookup&, const Target&);
    
    ///
    /// Constructor taking the training data with weights and the
    /// target vector as input.
    ///
    NCC(const MatrixLookupWeighted&, const Target&);

    virtual ~NCC();

    ///
    /// @return the centroids for each class as columns in a matrix.
    ///
    const utility::matrix& centroids(void) const;

    const DataLookup2D& data(void) const;

    SupervisedClassifier* make_classifier(const DataLookup2D&, 
                                          const Target&) const;
    
    ///
    /// Train the classifier using the training data. Centroids are
    /// calculated for each class.
    ///
    /// @return true if training succedeed.
    ///
    void train();

    
    ///
    /// Calculate the distance to each centroid for test samples
    ///
    void predict(const DataLookup2D&, utility::matrix&) const;
    
    
  private:

    void predict_unweighted(const MatrixLookup&, utility::matrix&) const;
    void predict_weighted(const MatrixLookupWeighted&, utility::matrix&) const;   

    utility::matrix* centroids_;
    bool centroids_nan_;
    Distance distance_;

    // data_ has to be of type DataLookup2D to accomodate both
    // MatrixLookup and MatrixLookupWeighted
    const DataLookup2D& data_;
  };

  ///
  /// The output operator for the NCC class.
  ///
  //  std::ostream& operator<< (std::ostream&, const NCC&);
  

  // templates

  template <typename Distance>
  NCC<Distance>::NCC(const MatrixLookup& data, const Target& target) 
    : SupervisedClassifier(target), centroids_(0), centroids_nan_(false), data_(data) 
  {
  }

  template <typename Distance>
  NCC<Distance>::NCC(const MatrixLookupWeighted& data, const Target& target)
    : SupervisedClassifier(target), centroids_(0), centroids_nan_(false), data_(data)
  {
  }

  template <typename Distance>
  NCC<Distance>::~NCC()   
  {
    if(centroids_)
      delete centroids_;
  }

  template <typename Distance>
  const utility::matrix& NCC<Distance>::centroids(void) const
  {
    return *centroids_;
  }
  

  template <typename Distance>
  const DataLookup2D& NCC<Distance>::data(void) const
  {
    return data_;
  }
  
  template <typename Distance>
  SupervisedClassifier* 
  NCC<Distance>::make_classifier(const DataLookup2D& data, const Target& target) const 
  {     
    NCC* ncc=0;
    try {
      if(data.weighted()) {
        ncc=new NCC<Distance>(dynamic_cast<const MatrixLookupWeighted&>(data),
                              target);
      }
      else {
        ncc=new NCC<Distance>(dynamic_cast<const MatrixLookup&>(data),
                              target);
      }
    }
    catch (std::bad_cast) {
      std::string str = "Error in NCC<Distance>::make_classifier: DataLookup2D of unexpected class.";
      throw std::runtime_error(str);
    }
    return ncc;
  }


  template <typename Distance>
  void NCC<Distance>::train()
  {   
    if(centroids_) 
      delete centroids_;
    centroids_= new utility::matrix(data_.rows(), target_.nof_classes());
    // data_ is a MatrixLookup or a MatrixLookupWeighted
    if(data_.weighted()) {
      const MatrixLookupWeighted* weighted_data = 
        dynamic_cast<const MatrixLookupWeighted*>(&data_);      
      for(size_t i=0; i<data_.rows(); i++) {
        std::vector<statistics::AveragerWeighted> class_averager;
        class_averager.resize(target_.nof_classes());
        for(size_t j=0; j<data_.columns(); j++) {
          class_averager[target_(j)].add(weighted_data->data(i,j),
                                         weighted_data->weight(i,j));
        }
        for(size_t c=0;c<target_.nof_classes();c++) {
          if(class_averager[c].sum_w()==0) {
            centroids_nan_=true;
            (*centroids_)(i,c) = std::numeric_limits<double>::quiet_NaN();
          }
          else {
            (*centroids_)(i,c) = class_averager[c].mean();
          }
        }
      }
    }
    else {
      const MatrixLookup* unweighted_data = 
        dynamic_cast<const MatrixLookup*>(&data_);      
      for(size_t i=0; i<data_.rows(); i++) {
        std::vector<statistics::Averager> class_averager;
        class_averager.resize(target_.nof_classes());
        for(size_t j=0; j<data_.columns(); j++) {
          class_averager[target_(j)].add((*unweighted_data)(i,j));
        }
        for(size_t c=0;c<target_.nof_classes();c++) {
          (*centroids_)(i,c) = class_averager[c].mean();
        }
      }
    }
  }

  template <typename Distance>
  void NCC<Distance>::predict(const DataLookup2D& test,                     
                              utility::matrix& prediction) const
  {   
    utility::yat_assert<std::runtime_error>
      (centroids_,"NCC::predict called for untrained classifier");
    utility::yat_assert<std::runtime_error>
      (data_.rows()==test.rows(),
       "NCC::predict test data with incorrect number of rows");
    
    prediction.resize(centroids_->columns(), test.columns());

    // unweighted test data
    if (const MatrixLookup* test_unweighted =
        dynamic_cast<const MatrixLookup*>(&test)) {
      // If weighted training data has resulted in NaN in centroids: weighted calculations
      if(centroids_nan_) { 
        predict_weighted(MatrixLookupWeighted(*test_unweighted),prediction);
      }
      // If unweighted training data: unweighted calculations
      else {
        predict_unweighted(*test_unweighted,prediction);
      }
    }
    // weighted test data: weighted calculations
    else if (const MatrixLookupWeighted* test_weighted =
             dynamic_cast<const MatrixLookupWeighted*>(&test)) { 
      predict_weighted(*test_weighted,prediction);
    }
    else {
      std::string str = 
        "Error in NCC<Distance>::predict: DataLookup2D of unexpected class.";
      throw std::runtime_error(str);
    }
  }
  
  template <typename Distance>
  void NCC<Distance>::predict_unweighted(const MatrixLookup& test, 
                                         utility::matrix& prediction) const
  {
    MatrixLookup unweighted_centroids(*centroids_);
    for(size_t j=0; j<test.columns();j++) {       
      DataLookup1D in(test,j,false);
      for(size_t k=0; k<centroids_->columns();k++) {
        DataLookup1D centroid(unweighted_centroids,k,false);            
        utility::yat_assert<std::runtime_error>(in.size()==centroid.size());
        prediction(k,j) = distance_(in.begin(), in.end(), centroid.begin());
      }
    }
  }

  template <typename Distance>
  void NCC<Distance>::predict_weighted(const MatrixLookupWeighted& test, 
                                          utility::matrix& prediction) const
  {
    MatrixLookupWeighted weighted_centroids(*centroids_);
    for(size_t j=0; j<test.columns();j++) {       
      DataLookupWeighted1D in(test,j,false);
      for(size_t k=0; k<centroids_->columns();k++) {
        DataLookupWeighted1D centroid(weighted_centroids,k,false);
        utility::yat_assert<std::runtime_error>(in.size()==centroid.size());
        prediction(k,j) = distance_(in.begin(), in.end(), centroid.begin());
      }
    }
  }

      
}}} // of namespace classifier, yat, and theplu

#endif