| sl2dmatcov(type, data, matsiz, n, meanmat, PL, PR, w)
|
function varargout = sl2dmatcov(type, data, matsiz, n, meanmat, PL, PR, w)
%SL2DMATCOV Computes the 2D matrix-covariances
%
% $ Syntax $
% - CL = sl2dmatcov('CL', data, matsiz, n, meanmat, PL, PR, w)
% - CR = sl2dmatcov('CR', data, matsiz, n, meanmat, PL, PR, w)
% - [CL, CR] = sl2dmatcov('Both', data, matsiz, n, meanmat, PL, PR, w)
%
% $ Arguments $
% - data: the stack of matrices or the cell array of filenames of the
% array files storing the matrices.
% - matsiz: the size of each matrix
% - n: the number of samples
% - PL: the left-projection matrix, default = []
% - PR: the right-projection matrix, default = []
% - w: the weights of the matrix samples, default = []
%
% $ Remarks $
% - The function computes the 2D matrix-covariances according to
% the following formulas:
% Y_i = PL^T * (X - meanX)* PR
% CL = ( sum_{i=1}^n w(i) * Y_i * Y_i^T ) / ( sum_{i=1}^n w(i) )
% CR = ( sum_{i=1}^n w(i) * Y_i^T * Y_i ) / ( sum_{i=1}^n w(i) )
% Following special cases are considered:
% If w is empty, then we take w(i) = 1 for all i
% If PL is empty, then we take PL as an identity matrix
% If PR is empty, then we take PR as an identity matrix
% If meanmat is 0, then we take meanmat as a zero matrix
%
% $ History $
% - Created by Dahua Lin, on Jul 31st, 2006
%
%% parse and verify input arguments
if nargin < 5
raise_lackinput('sl2dmatcov', 5);
end
if nargin < 6
PL = [];
end
if nargin < 7
PR = [];
end
if nargin < 8
w = [];
end
matsiz = matsiz(:)';
if length(matsiz) ~= 2
error('sltoolbox:invalidarg', ...
'matsiz shoudl be a 2-elem vector');
end
if ~isequal(size(meanmat), matsiz)
error('sltoolbox:invalidarg', ...
'The size of mean matrix is not as specified');
end
d1 = matsiz(1);
d2 = matsiz(2);
if ~isempty(PL) && size(PL, 1) ~= d1
error('sltoolbox:invalidarg', ...
'The size of PL is inconsistent with the matrix dimension');
end
if ~isempty(PR) && size(PR, 1) ~= d2
error('sltoolbox:invalidarg', ...
'The size of PR is inconsistent with the matrix dimension');
end
if ~isempty(w)
if length(w) ~= n
error('The weights length is inconsistent with the number of samples');
end
tw = sum(w);
else
tw = n;
end
%% Main body
if isnumeric(data)
if ~isequal(size(data), [matsiz, n])
error('sltoolbox:invalidarg', ...
'The size of data is inconsistent with specified');
end
Y = compute_Y(data, meanmat, PL, PR);
switch type
case 'CL'
CL = compute_SL(Y, w);
CL = CL / tw;
varargout = {CL};
case 'CR'
CR = compute_SR(Y, w);
CR = CR / tw;
varargout = {CR};
case 'Both'
[CL, CR] = compute_SLSR(Y, w);
CL = CL / tw;
CR = CR / tw;
varargout = {CL, CR};
otherwise
error('sltoolbox:invalidarg', ...
'invalid type: %s', type);
end
elseif iscell(data)
nfiles = length(data);
cf = 0;
switch type
case 'CL'
for i = 1 : nfiles
curarr = slreadarray(data{i});
curn = size(curarr, 3);
Y = compute_Y(curarr, meanmat, PL, PR);
curCL = compute_SL(Y, w);
if i == 1
CL = curCL;
else
CL = CL + curCL;
end
clear Y curCL;
cf = cf + curn;
end
CL = CL / tw;
varargout = {CL};
case 'CR'
for i = 1 : nfiles
curarr = slreadarray(data{i});
curn = size(curarr, 3);
Y = compute_Y(curarr, meanmat, PL, PR);
curCR = compute_SR(Y, w);
if i == 1
CR = curCR;
else
CR = CR + curCR;
end
clear Y curCR;
cf = cf + curn;
end
CR = CR / tw;
varargout = {CR};
case 'Both'
for i = 1 : nfiles
curarr = slreadarray(data{i});
curn = size(curarr, 3);
Y = compute_Y(curarr, meanmat, PL, PR);
[curCL, curCR] = compute_SLSR(Y, w);
if i == 1
CL = curCL;
CR = curCR;
else
CL = CL + curCL;
CR = CR + curCR;
end
clear Y curCL curCR;
cf = cf + curn;
end
CL = CL / tw;
CR = CR / tw;
varargout = {CL, CR};
otherwise
error('sltoolbox:invalidarg', ...
'invalid type: %s', type);
end
if cf ~= n
error('sltoolbox:sizmismatch', ...
'The total number of samples is not n');
end
else
error('sltoolbox:invalidarg', ...
'data should be a numeric array or a cell array of filenames');
end
%% Core routine
function Y = compute_Y(X, M, PL, PR)
n = size(X, 3);
if isempty(PL)
if isempty(PR)
if isequal(M, 0)
Y = X;
else
Y = zeros(size(X, 1), size(X, 2), n);
for i = 1 : n
Y(:,:,i) = X(:,:,i) - M;
end
end
else
Y = zeros(size(X, 1), size(PR, 2), n);
if isequal(M, 0)
for i = 1 : n
Y(:,:,i) = X(:,:,i) * PR;
end
else
for i = 1 : n
Y(:,:,i) = (X(:,:,i) - M) * PR;
end
end
end
else
PLT = PL';
if isempty(PR)
Y = zeros(size(PL, 2), size(X, 2), n);
if isequal(M, 0)
for i = 1 : n
Y(:,:,i) = PLT * X(:,:,i);
end
else
for i = 1 : n
Y(:,:,i) = PLT * (X(:,:,i) - M);
end
end
else
Y = zeros(size(PL, 2), size(PR, 2), n);
if isequal(M, 0)
for i = 1 : n
Y(:,:,i) = PLT * X(:,:,i) * PR;
end
else
for i = 1 : n
Y(:,:,i) = PLT * (X(:,:,i) - M) * PR;
end
end
end
end
function SL = compute_SL(Y, w)
SL = zeros(size(Y, 1));
n = size(Y, 3);
if isempty(w)
for i = 1 : n
curY = Y(:, :, i);
SL = SL + curY * curY';
end
else
for i = 1 : n
curY = Y(:, :, i);
SL = SL + w(i) * curY * curY';
end
end
function SR = compute_SR(Y, w)
SR = zeros(size(Y, 2));
n = size(Y, 3);
if isempty(w)
for i = 1 : n
curY = Y(:, :, i);
SR = SR + curY' * curY;
end
else
for i = 1 : n
curY = Y(:, :, i);
SR = SR + w(i) * curY' * curY;
end
end
function [SL, SR] = compute_SLSR(Y, w)
SL = zeros(size(Y, 1));
SR = zeros(size(Y, 2));
n = size(Y, 3);
if isempty(w)
for i = 1 : n
curY = Y(:, :, i);
SL = SL + curY * curY';
SR = SR + curY' * curY;
end
else
for i = 1 : n
curY = Y(:, :, i);
SL = SL + w(i) * curY * curY';
SR = SR + w(i) * curY' * curY;
end
end
|
|
