SHOGUN  v1.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
LinearLocalTangentSpaceAlignment.cpp
Go to the documentation of this file.
1 /*
2  * This program is free software; you can redistribute it and/or modify
3  * it under the terms of the GNU General Public License as published by
4  * the Free Software Foundation; either version 3 of the License, or
5  * (at your option) any later version.
6  *
7  * Written (W) 2011 Sergey Lisitsyn
8  * Copyright (C) 2011 Sergey Lisitsyn
9  */
10 
12 #ifdef HAVE_LAPACK
16 #include <shogun/io/SGIO.h>
17 #include <shogun/lib/Time.h>
19 #include <shogun/lib/Signal.h>
20 
21 using namespace shogun;
22 
25 {
26 }
27 
29 {
30 }
31 
33 {
34  return "LinearLocalTangentSpaceAlignment";
35 }
36 
38 {
39  CSimpleFeatures<float64_t>* simple_features = (CSimpleFeatures<float64_t>*)features;
40  ASSERT(simple_features);
41  int i,j;
42 
43  int N;
44  int dim;
45  float64_t* feature_matrix;
46  simple_features->get_feature_matrix(&feature_matrix,&dim,&N);
47  ASSERT(dimension<=dim);
48  float64_t* XTM = SG_MALLOC(float64_t, dim*N);
49  float64_t* lhs_M = SG_MALLOC(float64_t, dim*dim);
50  float64_t* rhs_M = SG_MALLOC(float64_t, dim*dim);
51  CMath::center_matrix(matrix.matrix,N,N);
52 
53  cblas_dgemm(CblasColMajor,CblasNoTrans,CblasNoTrans,dim,N,N,1.0,feature_matrix,dim,matrix.matrix,N,0.0,XTM,dim);
54  cblas_dgemm(CblasColMajor,CblasNoTrans,CblasTrans,dim,dim,N,1.0,XTM,dim,feature_matrix,dim,0.0,lhs_M,dim);
55 
56  float64_t* mean_vector = SG_CALLOC(float64_t, dim);
57  for (i=0; i<N; i++)
58  cblas_daxpy(dim,1.0,feature_matrix+i*dim,1,mean_vector,1);
59 
60  cblas_dscal(dim,1.0/N,mean_vector,1);
61 
62  for (i=0; i<N; i++)
63  cblas_daxpy(dim,-1.0,mean_vector,1,feature_matrix+i*dim,1);
64 
65  SG_FREE(mean_vector);
66 
67  cblas_dgemm(CblasColMajor,CblasNoTrans,CblasTrans,dim,dim,N,1.0,feature_matrix,dim,feature_matrix,dim,0.0,rhs_M,dim);
68 
69  for (i=0; i<dim; i++) rhs_M[i*dim+i] += 1e-6;
70 
71  float64_t* evals = SG_MALLOC(float64_t, dim);
72  float64_t* evectors = SG_MALLOC(float64_t, dimension*dim);
73  int32_t info = 0;
74 #ifdef HAVE_ARPACK
75  arpack_dsxupd(lhs_M,rhs_M,false,dim,dimension,"LA",false,3,true,m_nullspace_shift,0.0,
76  evals,evectors,info);
77 #else
78  wrap_dsygvx(1,'V','U',dim,lhs_M,dim,rhs_M,dim,dim-dimension+1,dim,evals,evectors,&info);
79 #endif
80  SG_FREE(lhs_M);
81  SG_FREE(rhs_M);
82  SG_FREE(evals);
83  if (info!=0) SG_ERROR("Failed to solve eigenproblem (%d)\n",info);
84 
85  for (i=0; i<dimension/2; i++)
86  {
87  cblas_dswap(dim,evectors+i*dim,1,evectors+(dimension-i-1)*dim,1);
88  }
89 
90  cblas_dgemm(CblasColMajor,CblasTrans,CblasNoTrans,N,dimension,dim,1.0,feature_matrix,dim,evectors,dim,0.0,XTM,N);
91  SG_FREE(evectors);
92  SG_FREE(feature_matrix);
93 
94  float64_t* new_features = SG_MALLOC(float64_t, dimension*N);
95  for (i=0; i<dimension; i++)
96  {
97  for (j=0; j<N; j++)
98  new_features[j*dimension+i] = XTM[i*N+j];
99  }
100  SG_FREE(XTM);
101  return SGMatrix<float64_t>(new_features,dimension,N);
102 }
103 
104 #endif /* HAVE_LAPACK */

SHOGUN Machine Learning Toolbox - Documentation