dualgreedykplstest.m

a function inside machine learning

%A script to compare greedyKPLS against craigs version 
clear all; 

tol = 10^-6; 

%[X, y, numExamples, numFeatures] = readCsvData('ionosphere.data');
[X, y, numExamples, numFeatures] = readCsvData('sonar-all.data');

X = normalise(X);

d = data; 
d = addDataField(d, 'X', X, 'examples'); 

T = 10; 

gKPLSParams.iterations = T; 
gKPLSParams.doubleDeflation = 1; 
gKPLSParams.dualFeatureDirection = 'dualMaxSparseKernelApprox'; 
gKPLSParams.X.kernel = getDefaultLinearKernel; 
gKPLSParams.Y.name = '';  %Ignore Y 
gKPLSParams.normalise = 0;

[subspaceInfo, trainInfo] = dualGeneralFeaturesTrain(d, gKPLSParams);
newTrainX = getDataFieldValue(trainInfo.data, 'X'); 

%Compute residuals 
trainK = X*X'; 
tau = newTrainX;
residuals = zeros(T, 1);

for i=1:T
    newTrainK = tau(:, 1:i)*inv(tau(:, 1:i)'*tau(:, 1:i))*tau(:, 1:i)'*trainK + trainK*tau(:, 1:i)*inv(tau(:, 1:i)'*tau(:, 1:i))*tau(:, 1:i)'- tau(:, 1:i)*inv(tau(:, 1:i)'*tau(:, 1:i))*tau(:, 1:i)'*trainK*tau(:, 1:i)*inv(tau(:, 1:i)'*tau(:, 1:i))*tau(:, 1:i)';

    residuals(i) = trace(trainK - newTrainK);
end

K = X*X'; 

[p,plsind,plsnorm] = greedyPLS(K,T);
[q] = getPLSProjections(p, K);

%Looks like we choose the same kernel columns at each iteration 
%i.e. p.tau == subspaceInfo.X.b
if norm(findNonZeroElements(subspaceInfo.X.b) - plsind) > tol 
    error('Choosing different kernel matrix columns at each iteration'); 
end 

%Next check residuals of deflated matrix 
if norm(residuals - plsnorm(2:end)') > tol 
    error('Residuals are wrong'); 
end 

if norm(newTrainX - q)  > tol 
    error('Features do not match those of craigs'); 
end