Struct petal_decomposition::RandomizedPca
source · [−]Expand description
Principal component analysis using randomized singular value decomposition.
This uses randomized SVD (singular value decomposition) proposed by Halko et al. [1] to reduce the dimensionality of the input data. The data is centered for each feature before applying randomized SVD.
Examples
use petal_decomposition::RandomizedPcaBuilder;
let x = ndarray::arr2(&[[0_f64, 0_f64], [1_f64, 1_f64], [2_f64, 2_f64]]);
let mut pca = RandomizedPcaBuilder::new(1).build();
let y = pca.fit_transform(&x).unwrap(); // [-2_f64.sqrt(), 0_f64, 2_f64.sqrt()]
assert!((y[(0, 0)].abs() - 2_f64.sqrt()).abs() < 1e-8);
assert!(y[(1, 0)].abs() < 1e-8);
assert!((y[(2, 0)].abs() - 2_f64.sqrt()).abs() < 1e-8);
References
- N. Halko, P. G. Martinsson, and J. A. Tropp. Finding Structure with Randomness: Probabilistic Algorithms for Constructing Approximate Matrix Decompositions. SIAM Review, 53(2), 217–288, 2011.
Implementations
sourceimpl<A> RandomizedPca<A, Pcg> where
A: Scalar,
impl<A> RandomizedPca<A, Pcg> where
A: Scalar,
sourcepub fn new(n_components: usize) -> Self
pub fn new(n_components: usize) -> Self
Creates a PCA model based on randomized SVD.
The random matrix for randomized SVD is created from a PCG random number generator (the XSL 128/64 (MCG) variant on a 64-bit CPU and the XSH RR 64/32 (LCG) variant on a 32-bit CPU), initialized with a randomly-generated seed.
sourcepub fn with_seed(n_components: usize, seed: u128) -> Self
pub fn with_seed(n_components: usize, seed: u128) -> Self
Creates a PCA model based on randomized SVD, with a PCG random number generator initialized with the given seed.
It uses a PCG random number generator (the XSL 128/64 (MCG) variant on a
64-bit CPU and the XSH RR 64/32 (LCG) variant on a 32-bit CPU). Use
with_rng
for a different random number generator.
sourceimpl<A, R> RandomizedPca<A, R> where
A: Scalar,
R: Rng,
impl<A, R> RandomizedPca<A, R> where
A: Scalar,
R: Rng,
sourceimpl<A, R> RandomizedPca<A, R> where
A: Scalar + FromPrimitive + Lapack,
A::Real: ScalarOperand + FromPrimitive,
R: Rng,
impl<A, R> RandomizedPca<A, R> where
A: Scalar + FromPrimitive + Lapack,
A::Real: ScalarOperand + FromPrimitive,
R: Rng,
sourcepub fn components(&self) -> &Array2<A>
pub fn components(&self) -> &Array2<A>
Returns the principal axes in feature space.
sourcepub fn n_components(&self) -> usize
pub fn n_components(&self) -> usize
Returns the number of components.
sourcepub fn singular_values(&self) -> &Array1<A::Real>
pub fn singular_values(&self) -> &Array1<A::Real>
Returns sigular values.
sourcepub fn explained_variance_ratio(&self) -> Array1<A::Real>
pub fn explained_variance_ratio(&self) -> Array1<A::Real>
Returns the ratio of explained variance for each component.
sourcepub fn fit<S>(
&mut self,
input: &ArrayBase<S, Ix2>
) -> Result<(), DecompositionError> where
S: Data<Elem = A>,
pub fn fit<S>(
&mut self,
input: &ArrayBase<S, Ix2>
) -> Result<(), DecompositionError> where
S: Data<Elem = A>,
Fits the model with input
.
Errors
DecompositionError::InvalidInput
if any of the dimensions ofinput
is less than the number of components, or the layout ofinput
is incompatible with LAPACK.DecompositionError::LinalgError
if the underlying Singular Vector Decomposition routine fails.
sourcepub fn transform<S>(
&self,
input: &ArrayBase<S, Ix2>
) -> Result<Array2<A>, DecompositionError> where
S: Data<Elem = A>,
pub fn transform<S>(
&self,
input: &ArrayBase<S, Ix2>
) -> Result<Array2<A>, DecompositionError> where
S: Data<Elem = A>,
Applies dimensionality reduction to input
.
Errors
DecompositionError::InvalidInput
if the number of features ininput
does not match that of the training data.
sourcepub fn fit_transform<S>(
&mut self,
input: &ArrayBase<S, Ix2>
) -> Result<Array2<A>, DecompositionError> where
S: Data<Elem = A>,
pub fn fit_transform<S>(
&mut self,
input: &ArrayBase<S, Ix2>
) -> Result<Array2<A>, DecompositionError> where
S: Data<Elem = A>,
Fits the model with input
and apply the dimensionality reduction on
input
.
This is equivalent to calling both fit
and transform
for the
same input.
Errors
DecompositionError::InvalidInput
if any of the dimensions ofinput
is less than the number of components, or the layout ofinput
is incompatible with LAPACK.DecompositionError::LinalgError
if the underlying Singular Vector Decomposition routine fails.
sourcepub fn inverse_transform<S>(
&self,
input: &ArrayBase<S, Ix2>
) -> Result<Array2<A>, DecompositionError> where
S: Data<Elem = A>,
pub fn inverse_transform<S>(
&self,
input: &ArrayBase<S, Ix2>
) -> Result<Array2<A>, DecompositionError> where
S: Data<Elem = A>,
Transforms data back to its original space.
Errors
Returns DecompositionError::InvalidInput
if the number of rows of
input
is different from that of the training data, or the number of
columns of input
is different from the number of components.
Auto Trait Implementations
impl<A, R> RefUnwindSafe for RandomizedPca<A, R> where
A: RefUnwindSafe,
R: RefUnwindSafe,
<A as Scalar>::Real: RefUnwindSafe,
impl<A, R> Send for RandomizedPca<A, R> where
A: Send,
R: Send,
<A as Scalar>::Real: Send,
impl<A, R> Sync for RandomizedPca<A, R> where
A: Sync,
R: Sync,
<A as Scalar>::Real: Sync,
impl<A, R> Unpin for RandomizedPca<A, R> where
R: Unpin,
<A as Scalar>::Real: Unpin,
impl<A, R> UnwindSafe for RandomizedPca<A, R> where
A: RefUnwindSafe,
R: UnwindSafe,
<A as Scalar>::Real: UnwindSafe + RefUnwindSafe,
Blanket Implementations
sourceimpl<T> BorrowMut<T> for T where
T: ?Sized,
impl<T> BorrowMut<T> for T where
T: ?Sized,
const: unstable · sourcefn borrow_mut(&mut self) -> &mut T
fn borrow_mut(&mut self) -> &mut T
Mutably borrows from an owned value. Read more