CRNkMedoids.h

Go to the documentation of this file.

/* Copyright 2014-2016 INSA-Lyon, ENS-Lyon
 * 
 * This file is part of libcrn.
 * 
 * libcrn is free software: you can redistribute it and/or modify
 * it under the terms of the GNU Lesser General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 * 
 * libcrn 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 Lesser General Public License for more details.
 * 
 * You should have received a copy of the GNU Lesser General Public License
 * along with libcrn.  If not, see <http://www.gnu.org/licenses/>.
 * 
 * file: CRNkMedoids.h
 * \author Yann LEYDIER
 */

#ifndef CRNkMedoids_HEADER
#define CRNkMedoids_HEADER

#include <CRNException.h>
#include <vector>
#include <map>
#include <tuple>
#include <limits>

namespace crn
{
    namespace kmedoids
    {
        namespace init
        {
            struct Central
            {
                Central(size_t n_classes):k(n_classes) {}
                std::vector<size_t> operator()(const std::vector<std::vector<double>> &distmat) const;
                private:
                    size_t k;
            };

            struct PAM
            {
                PAM(size_t n_classes):k(n_classes) {}
                std::vector<size_t> operator()(const std::vector<std::vector<double>> &distmat) const;
                private:
                    size_t k;
            };
        }
        namespace update
        {
            struct Local
            {
                void operator()(std::vector<size_t> &medoids, const std::vector<std::multimap<double, size_t>> &clusters, const std::vector<std::vector<double>> &distmat) const;
            };

            struct PAM
            {
                void operator()(std::vector<size_t> &medoids, const std::vector<std::multimap<double, size_t>> &clusters, const std::vector<std::vector<double>> &distmat) const;
            };
        }

        template<typename Init, typename Update> std::tuple<std::vector<size_t>, std::vector<std::multimap<double, size_t>>, std::vector<size_t>> Run(Init init, Update update, const std::vector<std::vector<double>> &distmat, size_t maxiter = std::numeric_limits<size_t>::max())
        {
            const auto nelem = distmat.size();
            // check
            for (const auto &row : distmat)
                if (row.size() != nelem)
                    throw ExceptionDimension("kmedoids::Run(): The distance matrix is not square.");
            // init
            auto medoids = init(distmat);
            const auto k = medoids.size();
            /*
            std::cout << "init: ";
            for (auto m : medoids)
                std::cout << m << " ";
            std::cout << std::endl;
            */
            // run
            auto clusters = std::vector<std::multimap<double, size_t>>(k);
            auto distsum = 0.0, precsum = std::numeric_limits<double>::max();
            auto iter = size_t{0};
            while (precsum != distsum)
            {
                precsum = distsum;
                distsum = 0.0;
                // classification
                std::vector<std::multimap<double, size_t>>(k).swap(clusters);
                for (auto o = size_t(0); o < nelem; ++o)
                {
                    auto ndist = std::numeric_limits<double>::max();
                    auto nclass = size_t{0};
                    for (auto m = size_t(0); m < k; ++m)
                        if (distmat[o][medoids[m]] < ndist)
                        {
                            ndist = distmat[o][medoids[m]];
                            nclass = m;
                        }
                    clusters[nclass].insert(std::make_pair(ndist, o));
                    distsum += ndist;
                }
                // update
                update(medoids, clusters, distmat);

                if (++iter >= maxiter)
                    break;
            }
            auto classnum = std::vector<size_t>(nelem, 0);
            for (auto c = size_t(0); c < k; ++c)
                for (const auto &o : clusters[c])
                    classnum[o.second] = c;
            return std::make_tuple(std::move(classnum), std::move(clusters), std::move(medoids));
        }
    }
}

#endif