diff --git a/fastlib2/contrib/nvasil/dual_manifold/dual_manifold_engine.h b/fastlib2/contrib/nvasil/dual_manifold/dual_manifold_engine.h index d0ca7b119c..db645358ce 100644 --- a/fastlib2/contrib/nvasil/dual_manifold/dual_manifold_engine.h +++ b/fastlib2/contrib/nvasil/dual_manifold/dual_manifold_engine.h @@ -20,16 +20,47 @@ #define DUAL_MANIFOLD_ENGINE_ #include "dual_manifold_objective.h" #include "../l_bfgs/l_bfgs.h" - +#include "../l_bfgs/optimization_utils.h" +/** + * DualManifoldEngine treats the problem of non-linear, non negative matrix factorization + * The classic problem of NMF is (\f& D \simeq WH\f&) + * For large scale problems D is sparse and it has a lot of zeros that represent + * missing or don't care data. This algorithm solves the following set of optimization + * problems: + * (\f$ \max trHH^T subject to A_i \bullet D = A_i \bullet (WH)) + * (\f$ \max trWW^T subject to A_i \bullet D = A_i \bullet (WH)) + * where (\f$\bullet \f$) is the matrix dot product + * (\f$ A_i \f$) is a selection matrix that selects the non-zero elements of (\f$D\f$) + * The matrix (\f& D \f&) can be real or non-negative. + * We can also restrict (\f$ W, H \f$) to be non-negative or sparse + * This depends on the definition of the OptimizedFunction + */ template class DualManifoldEngine { public: - void Init(datanode *module); + /** + * pairs_to_consider are (row,column) indices from a given sparse D + * matrix, These are the elements that we care about in our factorization and + * dot_prod_values are the values that we are trying to match + * + */ + void Init(datanode *module, + // index pairs to consider from the matrix (row,column) pairs + ArrayList > &pairs_to_consider, + // The values of the (row, column) values, also known as the dot products + ArrayList &dot_prod_values); void Destruct(); + void ComputeLocalOptimum(); + + private: + double feasibility_tolerance_; + double desired_feasibility_; + double norm_grad_tolerance_; + index_t iterations_; + index_t max_iterations_; + datanode *module; - private: - datanode *module; LBfgs lbfgs1_; LBfgs lbfgs2_; OptimizedFunction optimized_function1_; diff --git a/fastlib2/contrib/nvasil/dual_manifold/dual_manifold_engine_impl.h b/fastlib2/contrib/nvasil/dual_manifold/dual_manifold_engine_impl.h index c9e9c8fb7d..c3f6a7ca03 100644 --- a/fastlib2/contrib/nvasil/dual_manifold/dual_manifold_engine_impl.h +++ b/fastlib2/contrib/nvasil/dual_manifold/dual_manifold_engine_impl.h @@ -18,24 +18,39 @@ template DualManifoldEngine::Init(datanode *module, - ArrayList > &pairs_to_consider, + // index pairs to consider from the matrix (row,column) pairs + ArrayList > &pairs_to_consider, + // The values of the (row, column) values, also known as the dot products ArrayList &dot_prod_values) { module_=module; - lbfg1_.Init(&optimized_function1_, fx_param_node(module_, "lbfgs")); - lbfg2_.Init(&optimized_function2_, fx_param_node(module_, "lbfgs")); + l_bfg1_.Init(&optimized_function1_, fx_param_node(module_, "l_bfgs")); + l_bfg2_.Init(&optimized_function2_, fx_param_node(module_, "l_bfgs")); optimized_function1_.Init(fx_param_node(module_, "opt1"), - lbfgs1_.coordinates(), + l_bfgs1_.coordinates(), pairs_to_consider, dot_prod_values); for(index_t i=0; i +void DualManifoldEngine::ComputeLocalOptimum() { + for(index_t i=0; i > &pairs_to_consider, + // The values of the (row, column) values, also known as the dot products + ArrayList &dot_prod_values); + void ComputeGradient(Matrix &coordinates, Matrix *gradient); + void ComputeObjective(Matrix &coordinates, double *objective); + void ComputeFeasibilityError(Matrix &coordinates, double *error); + double ComputeLagrangian(Matrix &coordinates); + void UpdateLagrangeMult(Matrix &coordinates); + void Project(Matrix *coordinates); + void set_sigma(double sigma); + bool IsDiverging(double objective); + + private: + datanode *module_; + Matrix *auxiliary_mat_; + ArrayList > pairs_to_consider_; + ArrayList dot_prod_values_; + Vector eq_lagrange_mult_; + double sigma_; + index_t num_of_constraints_; +}; + + +#include "mvu_dot_prod_objective_impl.h" +#endif // MVU_DOT_PROD_OBJECTIVE_H_ diff --git a/fastlib2/contrib/nvasil/dual_manifold/mvu_dot_prod_objective_impl.h b/fastlib2/contrib/nvasil/dual_manifold/mvu_dot_prod_objective_impl.h new file mode 100644 index 0000000000..532fe8e25e --- /dev/null +++ b/fastlib2/contrib/nvasil/dual_manifold/mvu_dot_prod_objective_impl.h @@ -0,0 +1,115 @@ +/* + * ===================================================================================== + * + * Filename: mvu_dot_prod_objective_impl.h + * + * Description: + * + * Version: 1.0 + * Created: 04/09/2008 06:24:01 PM EDT + * Revision: none + * Compiler: gcc + * + * Author: Nikolaos Vasiloglou (NV), nvasil@ieee.org + * Company: Georgia Tech Fastlab-ESP Lab + * + * ===================================================================================== + */ + +void MVUDotProdObjective::Init(datanode *module, + Matrix *coordinates, + ArrayList > &pairs_to_consider, + // The values of the (row, column) values, also known as the dot products + ArrayList &dot_prod_values) { + + module_=module; + auxiliary_mat_=coordinates; + pairs_to_consider_.Copy(pairs_to_consider); + dot_prod_values_.Copy(dot_prod_values); + eq_lagrange_mult_.Init(dot_prod_values.size()); + eq_lagrange_mult_.SetAll(0.0); + num_of_constraints_=dot_prod_values_.size(); +} + +void MVUDotProdObjective::ComputeGradient(Matrix &coordinates, Matrix *gradient) { + gradient->CopyValues(coordinates); + // we need to use -CRR^T because we want to maximize CRR^T + la::Scale(-1.0, gradient); + index_t dimension=auxiliary_mat_->n_rows(); + Vector constant; + for (index_t i=0; iGetColumnPtr(ind1); + double *p2=coordinates.GetColumnPtr(ind2); + double dot_prod =la::Dot(dimension, p1, p2); + double diff=dot_prod-dot_prod_values_[i]; + la::AddExpert(dimension, + -eq_lagrange_mult_[i]+sigma_*diff, + p1, + gradient->GetColumnPtr(ind1)); + } +} + +void MVUDotProdObjective::ComputeObjective(Matrix &coordinates, double *objective) { + *objective=0; + index_t dimension = coordinates.n_rows(); + for(index_t i=0; i< coordinates.n_cols(); i++) { + *objective-=la::Dot(dimension, + coordinates.GetColumnPtr(i), + coordinates.GetColumnPtr(i)); + } +} + +void MVUDotProdObjective::ComputeFeasibilityError(Matrix &coordinates, double *error) { + *error=0; + for(index_t i=0; iGetColumnPtr(ind1); + double *p2=coordinates.GetColumnPtr(ind2); + double dot_prod =la::Dot(dimension, p1, p2); + double diff=dot_prod-dot_prod_values_[i]; + error +=diff*diff; + } +} + +double MVUDotProdObjective::ComputeLagrangian(Matrix &coordinates) { + double lagrangian=0; + ComputeObjective(coordinates, &lagrangian); + for(index_t i=0; iGetColumnPtr(ind1); + double *p2=coordinates.GetColumnPtr(ind2); + double dot_prod =la::Dot(dimension, p1, p2); + double diff=dot_prod-dot_prod_values_[i]; + lagrangian+= -eq_lagrange_mult_[i]*diff + sigma_*diff*diff; + } +} + +void MVUDotProdObjective::UpdateLagrangeMult(Matrix &coordinates) { + index_t dimension=coordinates.n_rows(); + for(index_t i=0; iGetColumnPtr(ind1); + double *p2=coordinates.GetColumnPtr(ind2); + double dot_prod =la::Dot(dimension, p1, p2); + double diff=dot_prod-dot_prod_values_[i]; + eq_lagrange_mult_[i]-=sigma_*diff; + } +} + +void MVUDotProdObjective::Project(Matrix *coordinates) { + OptUtils::RemoveMean(coordinates); +} + +void MVUDotProdObjective::set_sigma(double sigma) { + sigma_=sigma; +} + +bool MVUDotProdObjective::IsDiverging(double objective) { + return false; +} + diff --git a/fastlib2/contrib/nvasil/l_bfgs/l_bfgs.h b/fastlib2/contrib/nvasil/l_bfgs/l_bfgs.h index c566e7d961..935e634eb9 100644 --- a/fastlib2/contrib/nvasil/l_bfgs/l_bfgs.h +++ b/fastlib2/contrib/nvasil/l_bfgs/l_bfgs.h @@ -33,6 +33,9 @@ class LBfgs { void GetResults(Matrix *result); void Reset(); void set_coordinates(Matrix &coordinates); + void set_desired_feasibility(double desired_feasibility); + void set_feasibility_tolerance(double feasibility_tolerance); + void set_norm_grad_tolerance(double norm_grad_tolerance); Matrix *coordinates(); private: