From 32894cd743223b28fa4be7bfc162265bdd4b9464 Mon Sep 17 00:00:00 2001 From: Viktor Gal Date: Thu, 22 Feb 2018 21:12:29 +0100 Subject: [PATCH] fix GMM broke by #4077 and some more cleanup of Distributions API and implementation --- src/shogun/clustering/GMM.cpp | 57 ++++++++----------- .../distributions/DiscreteDistribution.h | 2 +- src/shogun/distributions/Distribution.cpp | 2 +- src/shogun/distributions/Distribution.h | 3 +- src/shogun/distributions/EMMixtureModel.cpp | 5 +- src/shogun/distributions/Gaussian.cpp | 22 +++---- src/shogun/distributions/Gaussian.h | 3 +- 7 files changed, 39 insertions(+), 55 deletions(-) diff --git a/src/shogun/clustering/GMM.cpp b/src/shogun/clustering/GMM.cpp index 7023063e2a3..50d76c67fd1 100644 --- a/src/shogun/clustering/GMM.cpp +++ b/src/shogun/clustering/GMM.cpp @@ -1,12 +1,13 @@ /* * This software is distributed under BSD 3-clause license (see LICENSE file). * - * Authors: Soeren Sonnenburg, Alesis Novik, Weijie Lin, Sergey Lisitsyn, - * Heiko Strathmann, Evgeniy Andreev, Chiyuan Zhang, Evan Shelhamer, + * Authors: Soeren Sonnenburg, Alesis Novik, Weijie Lin, Sergey Lisitsyn, + * Heiko Strathmann, Evgeniy Andreev, Chiyuan Zhang, Evan Shelhamer, * Wuwei Lin, Marcus Edel */ #include +#include #include #include #include @@ -197,8 +198,8 @@ float64_t CGMM::train_smem(int32_t max_iter, int32_t max_cand, float64_t min_cov if (m_components.size()<3) SG_ERROR("Can't run SMEM with less than 3 component mixture model.\n") - CDotFeatures* dotdata=(CDotFeatures *) features; - int32_t num_vectors=dotdata->get_num_vectors(); + CDotFeatures* dotdata = features->as(); + auto num_vectors = dotdata->get_num_vectors(); float64_t cur_likelihood=train_em(min_cov, max_em_iter, min_change); @@ -219,9 +220,9 @@ float64_t CGMM::train_smem(int32_t max_iter, int32_t max_cand, float64_t min_cov while (iter::twonorm(components[0]->get_mean().vector, dim_n)*0.1/ CMath::sqrt((float64_t)dim_n); - auto temp_mean = components[2]->get_mean(); - auto temp_mean_result = components[1]->get_mean(); - linalg::add(temp_mean_result, temp_mean, temp_mean_result, alpha1, alpha2); - components[1]->set_mean(temp_mean_result); + SGVector mean(dim_n); + linalg::add(components[1]->get_mean(), components[2]->get_mean(), mean, alpha1, alpha2); + components[1]->set_mean(mean); for (int32_t i=0; i eigenvalues(dim_n); - SGMatrix eigenvectors(dim_n, dim_n); - linalg::eigen_solver(c1, eigenvalues, eigenvectors); + linalg::eigen_solver_symmetric(c1, eigenvalues, c1); components[1]->set_d(eigenvalues); components[1]->set_u(c1); @@ -616,11 +615,10 @@ void CGMM::max_likelihood(SGMatrix alpha, float64_t min_cov) linalg::scale(cov_sum, cov_sum, 1.0 / alpha_sum); SGVector d0(num_dim); - SGMatrix eigenvectors(num_dim, num_dim); - linalg::eigen_solver(cov_sum, d0, eigenvectors); + linalg::eigen_solver_symmetric(cov_sum, d0, cov_sum); - for (int32_t j = 0; j < num_dim; j++) - d0[j] = CMath::max(min_cov, d0[j]); + for (auto& v: d0) + v = CMath::max(min_cov, v); m_components[i]->set_d(d0); m_components[i]->set_u(cov_sum); @@ -695,7 +693,7 @@ float64_t CGMM::get_likelihood_example(int32_t num_example) ASSERT(features->get_feature_class() == C_DENSE); ASSERT(features->get_feature_type() == F_DREAL); - for (int32_t i=0; i point= ((CDenseFeatures*) features)->get_feature_vector(num_example); result+=CMath::exp(m_components[i]->compute_log_PDF(point)+CMath::log(m_coefficients[i])); @@ -761,27 +759,20 @@ void CGMM::set_comp(vector components) SGMatrix CGMM::alpha_init(SGMatrix init_means) { CDotFeatures* dotdata=(CDotFeatures *) features; - int32_t num_vectors=dotdata->get_num_vectors(); + auto num_vectors=dotdata->get_num_vectors(); SGVector label_num(init_means.num_cols); + linalg::range_fill(label_num); - for (int32_t i=0; i(1, new CEuclideanDistance(), new CMulticlassLabels(label_num)); knn->train(new CDenseFeatures(init_means)); - CMulticlassLabels* init_labels=(CMulticlassLabels*) knn->apply(features); + auto init_labels = knn->apply(features)->as(); - SGMatrix alpha(num_vectors, int32_t(m_components.size())); - linalg::zero(alpha); - - for (int32_t i=0; i alpha(num_vectors, index_t(m_components.size())); + for (auto i: range(num_vectors)) alpha[i * m_components.size() + init_labels->get_int_label(i)] = 1; SG_UNREF(init_labels); - - delete knn; - return alpha; } @@ -791,7 +782,7 @@ SGVector CGMM::sample() "must be positive\n", m_components.size()); float64_t rand_num = CMath::random(0.0, 1.0); float64_t cum_sum=0; - for (int32_t i=0; i=rand_num) @@ -808,7 +799,7 @@ SGVector CGMM::cluster(SGVector point) SGVector answer(m_components.size()+1); answer.vector[m_components.size()]=0; - for (int32_t i=0; icompute_log_PDF(point)+CMath::log(m_coefficients[i]); answer.vector[m_components.size()]+=CMath::exp(answer.vector[i]); diff --git a/src/shogun/distributions/DiscreteDistribution.h b/src/shogun/distributions/DiscreteDistribution.h index e831b7cfd0c..7fc08fcd9a2 100644 --- a/src/shogun/distributions/DiscreteDistribution.h +++ b/src/shogun/distributions/DiscreteDistribution.h @@ -97,7 +97,7 @@ class CDiscreteDistribution : public CDistribution * * @param alpha_k "belongingness" values of various data points */ - virtual void update_params_em(SGVector alpha_k)=0; + virtual void update_params_em(const SGVector alpha_k)=0; }; } #endif /* _DISCRETEDISTRIBUTION_H__ */ \ No newline at end of file diff --git a/src/shogun/distributions/Distribution.cpp b/src/shogun/distributions/Distribution.cpp index 250cbd416e8..6568d9e90a7 100644 --- a/src/shogun/distributions/Distribution.cpp +++ b/src/shogun/distributions/Distribution.cpp @@ -70,7 +70,7 @@ SGVector CDistribution::get_likelihood_for_all_examples() return result; } -float64_t CDistribution::update_params_em(float64_t* alpha_k, int32_t len) +float64_t CDistribution::update_params_em(const SGVector alpha_k) { SG_WARNING("Not implemented in this class. This class cannot be used for Mixture models.\n") SG_NOTIMPLEMENTED diff --git a/src/shogun/distributions/Distribution.h b/src/shogun/distributions/Distribution.h index d1068f3fed4..954fcb034e4 100644 --- a/src/shogun/distributions/Distribution.h +++ b/src/shogun/distributions/Distribution.h @@ -188,10 +188,9 @@ class CDistribution : public CSGObject * abstract base method * * @param alpha_k "belongingness" values of various data points - * @param len length of alpha_k array * @return sum of alpha_k values */ - virtual float64_t update_params_em(float64_t* alpha_k, int32_t len); + virtual float64_t update_params_em(const SGVector alpha_k); /** obtain from generic * diff --git a/src/shogun/distributions/EMMixtureModel.cpp b/src/shogun/distributions/EMMixtureModel.cpp index d602dd1377a..ea0ced3d4b6 100644 --- a/src/shogun/distributions/EMMixtureModel.cpp +++ b/src/shogun/distributions/EMMixtureModel.cpp @@ -70,15 +70,14 @@ float64_t CEMMixtureModel::expectation_step() void CEMMixtureModel::maximization_step() { // for each component - float64_t* alpha_j=NULL; float64_t sum_weights=0; for (int32_t j=0;jget_element(j)->as(); // update mean covariance of components - alpha_j=data.alpha.matrix+j*data.alpha.num_rows; - float64_t weight_j=jth_component->update_params_em(alpha_j,data.alpha.num_rows); + SGVector alpha_j(data.alpha.matrix+j*data.alpha.num_rows, data.alpha.num_rows, false); + float64_t weight_j=jth_component->update_params_em(alpha_j); // update weights sum_weights+=weight_j; diff --git a/src/shogun/distributions/Gaussian.cpp b/src/shogun/distributions/Gaussian.cpp index 6a726f6addd..3f8d98b129a 100644 --- a/src/shogun/distributions/Gaussian.cpp +++ b/src/shogun/distributions/Gaussian.cpp @@ -1,7 +1,7 @@ /* * This software is distributed under BSD 3-clause license (see LICENSE file). * - * Authors: Soeren Sonnenburg, Weijie Lin, Alesis Novik, Heiko Strathmann, + * Authors: Soeren Sonnenburg, Weijie Lin, Alesis Novik, Heiko Strathmann, * Evgeniy Andreev, Viktor Gal, Evan Shelhamer, Björn Esser */ #include @@ -48,8 +48,8 @@ void CGaussian::init() { case FULL: case DIAG: - for (int32_t i=0; i alpha_k) { - CDotFeatures* dotdata=dynamic_cast(features); - REQUIRE( - dotdata, "dynamic cast from CFeatures to CDotFeatures returned NULL\n"); + CDotFeatures* dotdata=features->as(); int32_t num_dim=dotdata->get_dim_feature_space(); // compute mean - float64_t alpha_k_sum=0; SGVector mean(num_dim); - mean.fill_vector(mean.vector, mean.vlen, 0); - for (int32_t i = 0; i < len; i++) + linalg::zero(mean); + + for (auto i: range(alpha_k.vlen)) { alpha_k_sum+=alpha_k[i]; SGVector v=dotdata->get_computed_dot_feature_vector(i); @@ -156,7 +154,7 @@ float64_t CGaussian::update_params_em(float64_t* alpha_k, int32_t len) cov_sum.zero(); } - for (int32_t j=0; j v=dotdata->get_computed_dot_feature_vector(j); linalg::add(v, mean, v, -1.0, 1.0); @@ -299,7 +297,6 @@ void CGaussian::set_d(const SGVector d) SGMatrix CGaussian::get_cov() { SGMatrix cov(m_mean.vlen, m_mean.vlen); - cov.zero(); if (m_cov_type==FULL) { @@ -308,7 +305,6 @@ SGMatrix CGaussian::get_cov() SGMatrix temp_holder(m_mean.vlen, m_mean.vlen); SGMatrix diag_holder(m_mean.vlen, m_mean.vlen); - diag_holder.zero(); for (int32_t i = 0; i < m_d.vlen; i++) diag_holder(i, i) = m_d.vector[i]; #ifdef HAVE_LAPACK diff --git a/src/shogun/distributions/Gaussian.h b/src/shogun/distributions/Gaussian.h index 54b5bdad0d1..9059826f78c 100644 --- a/src/shogun/distributions/Gaussian.h +++ b/src/shogun/distributions/Gaussian.h @@ -100,10 +100,9 @@ class CGaussian : public CDistribution * this distribution is a part * * @param alpha_k "belongingness" values of various data points - * @param len length of alpha_k array * @return sum of values in alpha_k */ - virtual float64_t update_params_em(float64_t* alpha_k, int32_t len); + virtual float64_t update_params_em(const SGVector alpha_k); /** compute PDF *