buildgda.m

通过核来泛化的判别分析（GDA）代码

function dataGDA=BuildGDA(L,S)
% build the GDA data, with L learning vectors, and S an vector of the class sizes.
% L uses line vectors, S is a line vector

% Gaston Baudat & Fatiha Anouar / 21st October 2000 / Exton PA 19341 USA
% Designed under MatLab for Windows version 5.2.0.3084

n=length(S);
m=length(L(:,1));
mM=ones(m)/m;

% kernel matrix (uncentered)

kM=zeros(m);
for i=1:m
   for j=1:m
      kM(i,j)=KernelFunction(L(i,:),L(j,:));
   end
end

sumK=mM*kM;

% center the kernel matrix 

centeredkM=kM-sumK'-sumK+sumK*mM;
clear kM;
clear mM;

% decomposition of the centered kernel for beta eigenvectors

[vecCkM,valCkM]=eigensystem(centeredkM);

minVal=valCkM(1,1)/1000;  %define the lowest eigen value used
rankValCkM=m;
for i=1:m
   if valCkM(i,i)<minVal;
      valCkM(i,i)=0;
      rankValCkM=rankValCkM-1;
   end
end
fprintf('Estimated rank of K: %d\n',rankValCkM);
%fprintf('Estimated rank of K (matlab): %d\n',rank(centeredkM));
valCkM=valCkM(1:rankValCkM,1:rankValCkM);
vecCkM=vecCkM(:,1:rankValCkM);
centeredkM=vecCkM*valCkM*vecCkM';  %centered K matrix with only the highest eigen values/vectors

%w matrix of the weights, bloc diagonal
wM=zeros(m);
stBloc=1;
endBloc=0;
for i=1:n
   endBloc=endBloc+S(i);
   for j=stBloc:endBloc
      for k=stBloc:endBloc
         wM(j,k)=1/S(i);
      end
   end
   stBloc=stBloc+S(i);
end

%compute alpha normalized vectors
nbAxes=min([rankValCkM],[n-1]);
[vecBeta,valBeta]=eigensystem(vecCkM'*wM*vecCkM);
tmp=zeros(1,nbAxes);
for i=1:nbAxes
   tmp(i)=valBeta(i,i);
end
fprintf('Inertia: %f\n',tmp);
tmp=vecCkM*inv(valCkM)*vecBeta;
norAlpha=zeros(m,nbAxes);
for i=1:nbAxes
   norAlpha(:,i)=tmp(:,i)/sqrt(tmp(:,i)'*centeredkM*tmp(:,i));
end
sumK=sumK(1,:);
dataGDA=struct('norAlpha',norAlpha,'sumK',sumK);