PCA examples using Numerical Recipes?


ekeen
09-21-2010, 04:49 AM
Hi everyone,

I have big matrix to do PCA analysis. It is 200 rows * 1.2 million columns. I cannot use SVD method in Numerical Recipes to do it directly. So I first transform my matrix to 1.2 million rows and 200 columns. Then I can use SVD to it. After I run it, I get U,V,W three matrixs. So I take the V matrix as the PCs matrixs and each row in V is one PC. The first row is PC one, the second row is PC two, and so on. Is it right?

Thanks.

I used the fellow code to do it.

#include "clib/nr3.h"
#include "clib/svd.h"
#include "time.h"

ifstream svtest; // Input file name
ofstream svdU;
ofstream svdV;
ofstream svdW;
string filename; // user entered file name
Doub nrows,ncols,n,m;

int main(void)
{
clock_t start_time,finish_time;
Doub duration_time=0.0;
cout << "Enter the name of the input file " << endl;
cin >> filename;
svtest.open(filename.c_str());
svdU.open("U");
svdW.open("W");
svdV.open("V");
if(svtest.fail()) throw ("Error opening input file");
cout << "Enter the number of rows " << endl;
cin >> m;
cout << "Enter the number of columns " << endl;
cin >> n;
MatDoub a(m, n);
MatDoub uvw(m, n);
for(Int i=0;i<m;i++){
for(Int j=0;j<n;j++) {
svtest >> a[i][j];
}
}
SVD mySVD(a);
svdU << fixed << setprecision(7);
svdV << fixed << setprecision(7);
svdW << fixed << setprecision(7);
cout << endl << "Matrix U" << endl << endl;
for (Int i = 0; i < m; i++) {
for (Int j = 0; j < n; j++) {
svdU << "\t" << mySVD.u[i][j];
}
svdU << endl;
}
cout << endl;
cout << "Diagonal of matrix W" << endl << endl;
for(Int i=0;i<n;i++){
svdW << "\t"<< mySVD.w[i];
}
cout << endl;
cout << "Matrix V-Transpose" << endl;
for(Int i=0;i<n;i++){
for(Int j=0;j<n;j++){
svdV << "\t"<< mySVD.v[j][i];
}
svdV << endl;
}
finish_time = clock();
duration_time = (Doub)(finish_time-start_time)/CLOCKS_PER_SEC;
printf( "\nThe calculation time is: %f seconds.\n",duration_time );
return 0;
}

ekeen
09-21-2010, 11:14 AM
Does anyone can give a example?

Thanks.

antipin
08-11-2011, 07:44 AM
I also need a sample
somebody help me?