[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: how can I do Principal Components Analysis with octave?
From: |
Stefan van der Walt |
Subject: |
Re: how can I do Principal Components Analysis with octave? |
Date: |
Wed, 12 May 2004 10:23:58 +0200 |
User-agent: |
Mozilla Thunderbird 0.5 (X11/20040306) |
Maybe the attached klt.m will do the trick!
Regards
Stéfan
Henry F. Mollet wrote:
I believe you would need the programs doing the work. They are called
functions in Octave or Matlab, something like "PCA.m". I would be interested
myself so I checked using Google what MATLAB has:
Here is a list of the functions with a short description of each:
PRINCOMP - principal components from raw data matrix
PCACOV - pca from covariance matrix
PCARES - residuals from pca
BARTTEST - Bartlett's test for dimensionality.
Next I checked for the first two, namely "PRINCOMP" and "PCACOV" in
octave-forge but apparently neither is present. I guess we're out of luck
for the time being unless we have the capability to write the program.
Henry
on 5/11/04 1:04 AM, rino mailing at address@hidden wrote:
I'd like to do Principal Components Analysis with octave
What are the command I ave to write?
How to plot the result?
Thank you in advance for the time you spend to answer me, Mario.
## KLT Perform the Karhunen-Loeve (a.k.a. Hotelling or PCA) transform.
##
## [L, V] = klt(M, d or V[, options]);
##
## M is the feature matrix (one feature per row) to be transformed.
## The resulting matrix L has d dimensions (columns) unless
## V (the transformation matrix) is given, in which case it has
## the same nr of columns as V.
##
## The transformation is defined as
##
## y = V' * x;
##
## or
##
## y = V' * (x - mean_x) (if option 'submean' is specified)
##
## where x and y are column vectors and V the transformation
## matrix. V consists of the d largest eigenvectors.
##
## The transform rotates data so that it lies along the
## axis of greatest variance. It is often used to
## reduce dimensionality.
##
## By default the eigenvalues are only plotted if no output
## arguments are specified. Use option 'plot' to override.
##
## See Gonzales & Woods, Digital Image Processing
## (p.675, Use of Principle Components)
## Author: Stefan van der Walt <address@hidden>, 2004
function [L, V] = klt(varargin)
if (nargin < 2)
usage("[L, V] = klt(M, d or V[, submean])");
endif
M = varargin{1};
d = varargin{2};
if (!ismatrix(M))
error("klt: M must be a matrix");
endif
sub_mean = 0;
plot_always = 0;
if (length(varargin) > 2)
for i = 1:length(varargin)
if strcmpi(varargin{i}, "submean")
sub_mean = 1;
endif
if strcmpi(varargin{i}, "plot")
plot_always = 1;
endif
endfor
endif
if (!isscalar(d)) # V specified
V = d;
if (rows(d) != columns(M))
error("klt: invalid dimensions for transformation matrix V (should
be %d x D)", columns(M));
endif
else
if (d > columns(M))
warning("klt: asked for %d dimensions, but input only has %d!", [d
columns(M)]);
d = columns(M);
endif
C = cov(M);
[V, lambda] = eig(C);
[lambda, idx] = sort( diag(lambda) ); # sorted eigenvalues
V = V(:, idx ); # sorted eigenvectors
## decreasing order
lambda = flipud(lambda(:));
V = fliplr(V);
if ( (plot_always) || (nargout == 0) )
title("Eigenvalue decay");
plot(lambda, ";Lambda;");
endif
## choose the d eigenvectors corresponding to the d
## largest eigenvalues
V = V(:, 1:d);
endif
## transform every feature vector
if ( !sub_mean )
L = M * V;
else
L = ( M - repmat(mean(M), rows(M), 1) ) * V;
endif
endfunction
%!shared w
%!test
%! [x, w] = klt([sqrt(2) sqrt(2); -2*sqrt(2) -2*sqrt(2)], 1);
%! assert(x = [2; -4]);
%!test
%! x = klt([sqrt(2) sqrt(2); -2*sqrt(2) -2*sqrt(2)], w);
%! assert(x = [2; -2]);
%!test
%! x = klt([sqrt(2) sqrt(2); -2*sqrt(2) -2*sqrt(2)], w, 'submean');
%! assert(x = [3; -3]);
%!fail("klt(1)", "usage")