[mlpack-svn] r10654 - mlpack/trunk/src/mlpack/methods/radical
fastlab-svn at coffeetalk-1.cc.gatech.edu
fastlab-svn at coffeetalk-1.cc.gatech.edu
Thu Dec 8 00:13:23 EST 2011
Author: niche
Date: 2011-12-08 00:13:22 -0500 (Thu, 08 Dec 2011)
New Revision: 10654
Modified:
mlpack/trunk/src/mlpack/methods/radical/radical.cpp
mlpack/trunk/src/mlpack/methods/radical/radical.hpp
mlpack/trunk/src/mlpack/methods/radical/radical_main.cpp
Log:
various style fixes to radical
Modified: mlpack/trunk/src/mlpack/methods/radical/radical.cpp
===================================================================
--- mlpack/trunk/src/mlpack/methods/radical/radical.cpp 2011-12-08 04:17:23 UTC (rev 10653)
+++ mlpack/trunk/src/mlpack/methods/radical/radical.cpp 2011-12-08 05:13:22 UTC (rev 10654)
@@ -7,67 +7,63 @@
#include "radical.hpp"
-using namespace arma;
+//using namespace arma;
using namespace std;
-Radical::Radical() {
-}
-
-void Radical::Init(double noiseStdDev, u32 nReplicates, u32 nAngles, u32 nSweeps,
- const mat& X) {
+Radical::Radical(double noiseStdDev, size_t nReplicates, size_t nAngles, size_t nSweeps,
+ const arma::mat& matX) {
this->noiseStdDev = noiseStdDev;
this->nReplicates = nReplicates;
this->nAngles = nAngles;
this->nSweeps = nSweeps;
- this->X = mat(X); // is this the same as this.X = X ?
- m = floor(sqrt(X.n_rows));
+ this->matX = arma::mat(matX); // is this the same as this.X = X ?
+ m = floor(sqrt(matX.n_rows));
}
-mat Radical::GetX() {
- return X;
+arma::mat Radical::GetX() {
+ return matX;
}
-void Radical::WhitenX(mat& Whitening) {
- mat U, V;
- vec s;
- svd(U, s, V, cov(X));
- Whitening = U * diagmat(pow(s, -0.5)) * trans(V);
- Whitening.print("Whitening");
- X = X * Whitening;
+void Radical::WhitenX(arma::mat& matWhitening) {
+ arma::mat matU, matV;
+ arma::vec s;
+ arma::svd(matU, s, matV, cov(matX));
+ matWhitening = matU * arma::diagmat(pow(s, -0.5)) * arma::trans(matV);
+ matX = matX * matWhitening;
}
-void Radical::CopyAndPerturb(mat& XNew, const mat& X) {
- XNew = repmat(X, nReplicates, 1);
- XNew = XNew + noiseStdDev * randn(XNew.n_rows, XNew.n_cols);
+void Radical::CopyAndPerturb(arma::mat& matXNew, const arma::mat& matX) {
+ matXNew = arma::repmat(matX, nReplicates, 1);
+ matXNew = matXNew + noiseStdDev * arma::randn(matXNew.n_rows, matXNew.n_cols);
}
-double Radical::Vasicek(const vec& z) {
- vec v = sort(z);
- u32 nPoints = v.n_elem;
- vec logs = log(v.subvec(m, nPoints - 1) - v.subvec(0, nPoints - 1 - m));
- return (double) sum(logs);
+double Radical::Vasicek(const arma::vec& z) {
+ arma::vec v = sort(z);
+ size_t nPoints = v.n_elem;
+ arma::vec logs = arma::log(v.subvec(m, nPoints - 1) - v.subvec(0, nPoints - 1 - m));
+ return (double) arma::sum(logs);
}
-double Radical::DoRadical2D(const mat& X) {
- mat XMod;
- CopyAndPerturb(XMod, X);
+double Radical::DoRadical2D(const arma::mat& matX) {
+ arma::mat matXMod;
+ CopyAndPerturb(matXMod, matX);
double theta;
double value;
- mat Jacobi(2,2);
- mat candidateY;
+ arma::mat matJacobi(2,2);
+ arma::mat candidateY;
double thetaOpt = 0;
double valueOpt = 1e99;
- for(u32 i = 0; i < nAngles; i++) {
- theta = ((double) i) / ((double) nAngles) * M_PI / 2;
- Jacobi(0,0) = cos(theta);
- Jacobi(0,1) = sin(theta);
- Jacobi(1,0) = -sin(theta);
- Jacobi(1,1) = cos(theta);
+ for(size_t i = 0; i < nAngles; i++) {
+ theta = ((double) i) / ((double) nAngles) * arma::math::pi() / 2;
+ matJacobi(0,0) = cos(theta);
+ matJacobi(0,1) = sin(theta);
+ matJacobi(1,0) = -sin(theta);
+ matJacobi(1,1) = cos(theta);
- candidateY = XMod * Jacobi;
+ candidateY = matXMod * matJacobi;
value = Vasicek(candidateY.col(0)) + Vasicek(candidateY.col(1));
if(value < valueOpt) {
valueOpt = value;
@@ -79,43 +75,45 @@
}
-void Radical::DoRadical(mat& Y, mat& W) {
+void Radical::DoRadical(arma::mat& matY, arma::mat& matW) {
- // X is nPoints by nDims (although less intuitive, this is for computational efficiency)
+ // matX is nPoints by nDims (although less intuitive than columns being points,
+ // and although this is the transpose of the ICA literature,this choice is
+ // for computational efficiency)
- u32 nDims = X.n_cols;
+ size_t nDims = matX.n_cols;
- mat Whitening;
- WhitenX(Whitening);
+ arma::mat matWhitening;
+ WhitenX(matWhitening);
// in the RADICAL code, they do not copy and perturb initial, although the
// paper does. we follow the code as it should match their reported results
// and likely does a better job bouncing out of local optima
//GeneratePerturbedX(X, X);
- Y = X;
- W = eye(nDims, nDims);
+ matY = matX;
+ matW = arma::eye(nDims, nDims);
- mat YSubspace(X.n_rows, 2);
+ arma::mat matYSubspace(matX.n_rows, 2);
- for(u32 sweepNum = 0; sweepNum < nSweeps; sweepNum++) {
- for(u32 i = 0; i < nDims - 1; i++) {
- for(u32 j = i + 1; j < nDims; j++) {
- YSubspace.col(0) = Y.col(i);
- YSubspace.col(1) = Y.col(j);
- mat WSubspace;
- double thetaOpt = DoRadical2D(YSubspace);
- mat J = eye(nDims, nDims);
- J(i,i) = cos(thetaOpt);
- J(i,j) = sin(thetaOpt);
- J(j,i) = -sin(thetaOpt);
- J(j,j) = cos(thetaOpt);
- W = W * J;
- Y = Y * W;
+ for(size_t sweepNum = 0; sweepNum < nSweeps; sweepNum++) {
+ for(size_t i = 0; i < nDims - 1; i++) {
+ for(size_t j = i + 1; j < nDims; j++) {
+ matYSubspace.col(0) = matY.col(i);
+ matYSubspace.col(1) = matY.col(j);
+ arma::mat matWSubspace;
+ double thetaOpt = DoRadical2D(matYSubspace);
+ arma::mat matJ = arma::eye(nDims, nDims);
+ matJ(i,i) = cos(thetaOpt);
+ matJ(i,j) = sin(thetaOpt);
+ matJ(j,i) = -sin(thetaOpt);
+ matJ(j,j) = cos(thetaOpt);
+ matW = matW * matJ;
+ matY = matY * matW;
}
}
}
// the final transpose provides W in the typical form from the ICA literature
- W = trans(W * Whitening);
+ matW = arma::trans(matW * matWhitening);
}
Modified: mlpack/trunk/src/mlpack/methods/radical/radical.hpp
===================================================================
--- mlpack/trunk/src/mlpack/methods/radical/radical.hpp 2011-12-08 04:17:23 UTC (rev 10653)
+++ mlpack/trunk/src/mlpack/methods/radical/radical.hpp 2011-12-08 05:13:22 UTC (rev 10654)
@@ -14,31 +14,29 @@
#include <stdio.h>
#include<float.h>
-using namespace arma;
using namespace std;
class Radical {
public:
- Radical();
+ Radical(double noiseStdDev, size_t nReplicates, size_t nAngles, size_t nSweeps,
+ const arma::mat& matX);
+
+ arma::mat GetX();
+ void WhitenX(arma::mat& matWhitening);
+ void CopyAndPerturb(arma::mat& XNew, const arma::mat& matX);
+ double Vasicek(const arma::vec& x);
+ double DoRadical2D(const arma::mat& matX);
+ void DoRadical(arma::mat& matY, arma::mat& matW);
- void Init(double noiseStdDev, u32 nReplicates, u32 nAngles, u32 nSweeps,
- const mat& X);
- mat GetX();
- void WhitenX(mat& Whitening);
- void CopyAndPerturb(mat& XNew, const mat& X);
- double Vasicek(const vec& x);
- double DoRadical2D(const mat& X);
- void DoRadical(mat& Y, mat& W);
-
private:
double noiseStdDev;
- u32 nReplicates;
- u32 nAngles;
- u32 nSweeps;
- u32 m; // for Vasicek's m-spacing estimator
+ size_t nReplicates;
+ size_t nAngles;
+ size_t nSweeps;
+ size_t m; // for Vasicek's m-spacing estimator
- mat X;
-
+ arma::mat matX;
+
};
#endif
Modified: mlpack/trunk/src/mlpack/methods/radical/radical_main.cpp
===================================================================
--- mlpack/trunk/src/mlpack/methods/radical/radical_main.cpp 2011-12-08 04:17:23 UTC (rev 10653)
+++ mlpack/trunk/src/mlpack/methods/radical/radical_main.cpp 2011-12-08 05:13:22 UTC (rev 10654)
@@ -8,13 +8,12 @@
#include "radical.hpp"
using namespace std;
+using namespace arma;
int main(int argc, char* argv[]) {
- Radical rad;
+ size_t nPoints = 1000;
+ size_t nDims = 2;
- u32 nPoints = 1000;
- u32 nDims = 2;
-
mat S = randu(nPoints, nDims);
//S.print("S");
mat Mixing = randn(nDims, nDims);
@@ -29,7 +28,13 @@
X = X * Whitening;
*/
- rad.Init(0.01, 10, 200, 1, X);
+
+ //Radical rad;
+ //rad.Init(0.01, 10, 200, 1, X);
+
+ Radical rad(0.01, 10, 200, 1, X);
+
+
mat Y;
mat W;
rad.DoRadical(Y, W);
@@ -37,14 +42,14 @@
W.print("W");
double val_init = 0;
- for(u32 i = 0; i < nDims; i++) {
+ for(size_t i = 0; i < nDims; i++) {
val_init += rad.Vasicek(rad.GetX().col(i));
}
printf("initial objective value: %f\n", val_init);
double val_final = 0;
- for(u32 i = 0; i < nDims; i++) {
+ for(size_t i = 0; i < nDims; i++) {
val_final += rad.Vasicek(Y.col(i));
}
printf("final objective value: %f\n", val_final);
More information about the mlpack-svn
mailing list