Vespucci  1.0.0
Vespucci::Math::DimensionReduction Namespace Reference

A namespace for math functions relating to dimension reduction. More...

Functions

VESPUCCI_EXPORT bool plsregress (arma::mat X, arma::mat Y, int components, arma::mat &X_loadings, arma::mat &Y_loadings, arma::mat &X_scores, arma::mat &Y_scores, arma::mat &coefficients, arma::mat &percent_variance, arma::mat &fitted)
 Vespucci::MathDimensionReduction::plsregress PLS Regression Using SIMPLS algorithm. This is essentially a line-for-line rewrite of plsregress from the Octave statistics package. Copyright (C) 2012 Fernando Damian Nieuwveldt fdnie.nosp@m.uwve.nosp@m.ldt@g.nosp@m.mail.nosp@m..com This is an implementation of the SIMPLS algorithm: Reference: SIMPLS: An alternative approach to partial least squares regression. Chemometrics and Intelligent Laboratory Systems (1993) More...
 
VESPUCCI_EXPORT bool VCA (const arma::mat &R, arma::uword p, arma::uvec &indices, arma::mat &endmember_spectra, arma::mat &projected_data, arma::mat &fractional_abundances)
 Vespucci::Math::DimensionReduction::VCA Vertex Component Analysis. More...
 
VESPUCCI_EXPORT double estimate_snr (const arma::mat &R, arma::vec r_m, arma::mat x)
 Vespucci::Math Estimates Signal-to-Noise ratio. [Nascimento2005]. More...
 
VESPUCCI_EXPORT size_t HySime (const arma::mat &y, const arma::mat &n, arma::mat &Rn, arma::mat &Ek)
 Vespucci::Math::DimensionReduction::HySime. More...
 
VESPUCCI_EXPORT size_t HySime (const arma::mat &y, arma::mat &subspace)
 
VESPUCCI_EXPORT void EstimateAdditiveNoise (arma::mat &noise, arma::mat &noise_correlation, const arma::mat &sample)
 Vespucci::Math::DimensionReduction::EstimateAdditiveNoise. More...
 
VESPUCCI_EXPORT bool svds (const arma::mat &X, arma::uword k, arma::mat &U, arma::vec &s, arma::mat &V)
 Vespucci::Math::DimensionReduction::svds Finds a few largest singular values of the arma::matrix X. This is based on the Matlab/Octave function svds(), a truncated singular value decomposition. This is designed to only take the kinds of inputs Vespucci will need (It only works on arma::mat types, and only returns the k largest singular values (none of that sigma business). U and V tolerances are not defined currently due to difficulties with Armadillo's find() method A sparse arma::matrix [0 X; X.t() 0] is formed. The eigenvalues of this arma::matrix are then found using arma::eigs_sym, a wrapper for ARPACK. If X is square, it can be "reconstructed" using X = U*diamat(s)*V.t(). This "reconstruction" will have lower noise. More...
 

Detailed Description

A namespace for math functions relating to dimension reduction.

A namespace for math functions relating to transforms.

Function Documentation

double Vespucci::Math::DimensionReduction::estimate_snr ( const arma::mat &  R,
arma::vec  r_m,
arma::mat  x 
)

Vespucci::Math Estimates Signal-to-Noise ratio. [Nascimento2005].

Parameters
RInput
r_m
x
Returns

Definition at line 136 of file VCA.cpp.

void Vespucci::Math::DimensionReduction::EstimateAdditiveNoise ( arma::mat &  noise,
arma::mat &  noise_correlation,
const arma::mat &  sample 
)

Vespucci::Math::DimensionReduction::EstimateAdditiveNoise.

Parameters
noise
noise_correlation
sampleThe slow additive noise mat::mator from the HySime paper. I'm looking into faster ones.

Definition at line 216 of file VCA.cpp.

size_t Vespucci::Math::DimensionReduction::HySime ( const arma::mat &  y,
const arma::mat &  n,
arma::mat &  Rn,
arma::mat &  Ek 
)

Vespucci::Math::DimensionReduction::HySime.

Parameters
y
n
Rn
Ek
Returns
Performs the HySime algorithm to predict the rank of y

Definition at line 159 of file VCA.cpp.

size_t Vespucci::Math::DimensionReduction::HySime ( const arma::mat &  y,
arma::mat &  subspace 
)

Definition at line 201 of file VCA.cpp.

bool Vespucci::Math::DimensionReduction::plsregress ( arma::mat  X,
arma::mat  Y,
int  components,
arma::mat &  X_loadings,
arma::mat &  Y_loadings,
arma::mat &  X_scores,
arma::mat &  Y_scores,
arma::mat &  coefficients,
arma::mat &  percent_variance,
arma::mat &  fitted 
)

Vespucci::MathDimensionReduction::plsregress PLS Regression Using SIMPLS algorithm. This is essentially a line-for-line rewrite of plsregress from the Octave statistics package. Copyright (C) 2012 Fernando Damian Nieuwveldt fdnie.nosp@m.uwve.nosp@m.ldt@g.nosp@m.mail.nosp@m..com This is an implementation of the SIMPLS algorithm: Reference: SIMPLS: An alternative approach to partial least squares regression. Chemometrics and Intelligent Laboratory Systems (1993)

Parameters
x
y
componentsNumber of PLS components
x_loadings"Output" values
y_loadings
x_scores
y_scores
coefficients
fitted_data
Returns

Definition at line 43 of file pls.cpp.

bool Vespucci::Math::DimensionReduction::svds ( const arma::mat &  X,
arma::uword  k,
arma::mat &  U,
arma::vec &  s,
arma::mat &  V 
)

Vespucci::Math::DimensionReduction::svds Finds a few largest singular values of the arma::matrix X. This is based on the Matlab/Octave function svds(), a truncated singular value decomposition. This is designed to only take the kinds of inputs Vespucci will need (It only works on arma::mat types, and only returns the k largest singular values (none of that sigma business). U and V tolerances are not defined currently due to difficulties with Armadillo's find() method A sparse arma::matrix [0 X; X.t() 0] is formed. The eigenvalues of this arma::matrix are then found using arma::eigs_sym, a wrapper for ARPACK. If X is square, it can be "reconstructed" using X = U*diamat(s)*V.t(). This "reconstruction" will have lower noise.

This version is different from the armadillo version because it handles dense matrices.

Parameters
XAn m n arma::matrix.
kThe number of eigenvalues to be found.
UAn m k unitary arma::matrix
sa vector with n
V
Returns
Whether or not algorithm converged.

Definition at line 24 of file svds.cpp.

bool Vespucci::Math::DimensionReduction::VCA ( const arma::mat &  R,
arma::uword  p,
arma::uvec &  indices,
arma::mat &  endmember_spectra,
arma::mat &  projected_data,
arma::mat &  fractional_abundances 
)

Vespucci::Math::DimensionReduction::VCA Vertex Component Analysis.

Parameters
RThe dataset
endmembersNumber of endmembers to compute
indicesRow indices of pure components.
endmember_spectraSpectra of pure components (note that these are in columns, not rows as in spectra_)
projected_dataProjected data
fractional_abundancesPurity of a given spectrum relative to endmember
Returns
Convergeance (no actual test implemented...)

Definition at line 36 of file VCA.cpp.