00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011 #include <shogun/converter/NeighborhoodPreservingEmbedding.h>
00012 #include <shogun/lib/config.h>
00013 #ifdef HAVE_LAPACK
00014 #include <shogun/mathematics/arpack.h>
00015 #include <shogun/mathematics/lapack.h>
00016 #include <shogun/lib/FibonacciHeap.h>
00017 #include <shogun/mathematics/Math.h>
00018 #include <shogun/io/SGIO.h>
00019 #include <shogun/lib/Time.h>
00020 #include <shogun/distance/Distance.h>
00021 #include <shogun/lib/Signal.h>
00022
00023 using namespace shogun;
00024
00025 CNeighborhoodPreservingEmbedding::CNeighborhoodPreservingEmbedding() :
00026 CLocallyLinearEmbedding()
00027 {
00028 }
00029
00030 CNeighborhoodPreservingEmbedding::~CNeighborhoodPreservingEmbedding()
00031 {
00032 }
00033
00034 const char* CNeighborhoodPreservingEmbedding::get_name() const
00035 {
00036 return "NeighborhoodPreservingEmbedding";
00037 }
00038
00039 SGMatrix<float64_t> CNeighborhoodPreservingEmbedding::construct_embedding(CFeatures* features, SGMatrix<float64_t> matrix, int dimension)
00040 {
00041 CSimpleFeatures<float64_t>* simple_features = (CSimpleFeatures<float64_t>*)features;
00042 ASSERT(simple_features);
00043 int i,j;
00044 int N = simple_features->get_num_vectors();
00045 int dim = simple_features->get_num_features();
00046 ASSERT(dimension<=dim);
00047
00048 SGMatrix<float64_t> feature_matrix = simple_features->get_feature_matrix();
00049 float64_t* XTM = SG_MALLOC(float64_t, dim*N);
00050 float64_t* lhs_M = SG_MALLOC(float64_t, dim*dim);
00051 float64_t* rhs_M = SG_MALLOC(float64_t, dim*dim);
00052
00053 cblas_dgemm(CblasColMajor,CblasNoTrans,CblasNoTrans,dim,N,N,1.0,feature_matrix.matrix,dim,matrix.matrix,N,0.0,XTM,dim);
00054 cblas_dgemm(CblasColMajor,CblasNoTrans,CblasTrans,dim,dim,N,1.0,XTM,dim,feature_matrix.matrix,dim,0.0,lhs_M,dim);
00055 cblas_dgemm(CblasColMajor,CblasNoTrans,CblasTrans,dim,dim,N,1.0,feature_matrix.matrix,dim,feature_matrix.matrix,dim,0.0,rhs_M,dim);
00056
00057 float64_t* evals = SG_MALLOC(float64_t, dim);
00058 float64_t* evectors = SG_MALLOC(float64_t, dimension*dim);
00059 int32_t info = 0;
00060 #ifdef HAVE_ARPACK
00061 arpack_dsxupd(lhs_M,rhs_M,false,dim,dimension,"LA",false,3,true,m_nullspace_shift,0.0,
00062 evals,evectors,info);
00063 #else
00064 wrap_dsygvx(1,'V','U',dim,lhs_M,dim,rhs_M,dim,dim-dimension+1,dim,evals,evectors,&info);
00065 #endif
00066 SG_FREE(lhs_M);
00067 SG_FREE(rhs_M);
00068 SG_FREE(evals);
00069 if (info!=0) SG_ERROR("Failed to solve eigenproblem (%d)\n",info);
00070
00071 cblas_dgemm(CblasColMajor,CblasTrans,CblasNoTrans,N,dimension,dim,1.0,feature_matrix.matrix,dim,evectors,dim,0.0,XTM,N);
00072 SG_FREE(evectors);
00073
00074 float64_t* new_features = SG_MALLOC(float64_t, dimension*N);
00075 for (i=0; i<dimension; i++)
00076 {
00077 for (j=0; j<N; j++)
00078 new_features[j*dimension+i] = XTM[i*N+j];
00079 }
00080 SG_FREE(XTM);
00081 return SGMatrix<float64_t>(new_features,dimension,N);
00082 }
00083
00084 #endif