K-CCA MATLAB 代码

function [nalpha, nbeta, r, Kx, Ky] = kcanonca_reg_ver2(Kx,Ky,eta,kapa,sl,nor,Rx,Ry) % KCCA function - regularized with kapa % % Performes Kernel Canoncial Correlation Analysis as a single symmetric eigenvalue problem with % parital gram schmidt kernel decomposition % % Usage: % [nalpha, nbeta, r, Kx_c, Ky_c] = kcanonca_reg_ver2(Kx,Ky,eta,kapa) % [nalpha, nbeta, r, Kx_c, Ky_c] = kcanonca_reg_ver2(Kx,Ky,eta,kapa,sl) % [nalpha, nbeta, r, Kx_c, Ky_c] = kcanonca_reg_ver2(Kx,Ky,eta,kapa,sl,nor) % [nalpha, nbeta, r] = kcanonca_reg_ver2(Kx,Ky,eta,kapa,sl,nor,Rx,Ry) % % Input: % Kx, Ky - Kernel matrices corresponding to the two views, % eta - precision parameter for gsd % kapa - regularsation parameter (value 0 to 1) % sl - 1 | 0 save output from gsd method (default 1) % nor - normalise features (defualt 2) % 0 - No normalisation and projection to kernel space % 1 - No normalisation and and leaving GS projected space % 2 - Normalisation and projection back to Kernel space % 3 - Normalisation and leaving in GS projected space % Rx, Ry - predecomposed Kx and Ky % flag - true | false, whether to run the gram-schmidt algorithm on Kx Ky % % Output: % nalpha, nbeta - contains the canonical correlation vectors as columns % rc - is a vector with corresponding canonical correlations. % Kx_c, Ky_c - The centered kernels % % Non commercial use. % Written by David R. Hardoon while at University of Southampton ISIS Group, Southampton UK. % Current email - D.Hardoon@cs.ucl.ac.uk % % Updated 15/12/2004 - just tidying things up % Updated 17/12/2004 - Added control over whether to normalise or not % Updated 28/01/2005 - when not normalising r is now a vector rather then a matrix % - Also added more options normalisation and proejction (see above) % Updated 01/02/2005 - I had case 1|0 for the nor paramtere reveresed % Updated 06/01/2006 - Added kernel centralising % Updated 31/04/2008 - Added output for kernels when centered % instead of gsd one could use icd % checking correct number of inputs if (nargin < 4 | nargin > 8 | nargin == 7) disp('Incorrect number of inputs') nalpha = 0; nbeta = 0; r = 0; help kcanonca_reg_ver2 return end % Setting defaults if (nargin < 5) sl = 1; nor = 2; elseif (nargin < 6) nor = 2; end % need to decompose kernels if (nargin < 7) % If you centre outside comment out this bit disp('Centering Kx and Ky'); l = size(Kx,1); j = ones(l,1); Kx = Kx - (j*j'*Kx)/l - (Kx*j*j')/l + ((j'*Kx*j)*j*j')/(l^2); l = size(Ky,1); j = ones(l,1); Ky = Ky - (j*j'*Ky)/l - (Ky*j*j')/l + ((j'*Ky*j)*j*j')/(l^2); disp('Decomposing Kernel with PGSO'); [Rx, RxSize, RxIndex] = gsd(Kx,eta); [Ry, RySize, RyIndex] = gsd(Ky,eta); if (sl == 1) disp('Saving PGSO as cca_feat_data'); save('cca_feat_data','Rx','RxSize','RxIndex','Ry','RySize','RyIndex'); end else disp('Did you centre your data?'); end %disp('Creating new TxT matrix Z from MxT matrix R'); Zxx = Rx'*Rx; Zxy = Rx'*Ry; Zyy = Ry'*Ry; Zyx = Zxy'; tEyeY = eye(size(Zyy)); tEyeX = eye(size(Zxx)); %disp('Computing nalpha eigenproblem'); B = (1-kapa)*tEyeX*Zxx+kapa*tEyeX; S = chol(B)'; invS = inv(S); A = invS*Zxy*inv((1-kapa)*tEyeY*Zyy+kapa*tEyeY)*Zyx*invS'; A = 0.5*(A'+A)+eye(size(A,1))*10e-6; [pha,rr] = eig(A); %disp('Sorting Output of nalpha'); r = sqrt(real(rr)); % as the original r we get is lamda^2 alpha = invS'*pha; % as \hat{\alpha} = S'*\alpha - this find alpha tidal % if you do not do the following line it means you will need to project your % testing data into the Gram-Shmidt space. invRx = Rx*inv(Rx'*Rx); nalpha = invRx*alpha; disp('Computing nbeta from nalpha'); % computing beta -- but as we cant comput the original beta % we can only compute beta twidel and not the original beta beta = inv((1-kapa)*tEyeY*Zyy+kapa*tEyeY)*Zyx*alpha; t = size(Zyy,1); beta = beta./repmat(diag(r)',t,1); % again if you do not do the following line you must project your % corresponding testing examples to gsd space invRy = Ry*inv(Ry'*Ry); nbeta = invRy*beta; % normalising the feature vectors and recomputing the correlation switch nor case 1 % no normalisation but using GS projected space nbeta = beta; nalpha = alpha; r = diag(r); disp('WARNING: Remember to project your testing Kernels into the GS space!'); case 2 % normalisation and kernel space disp('Normalising compute features'); nalpha = normalise(nalpha,Kx); nbeta = normalise(nbeta,Ky); r = diag(nalpha'*Kx'*Ky*nbeta); case 3 % normalisation in the GS space disp('Normalising compute features'); nalpha = normalise(alpha,Rx); nbeta = normalise(beta,Ry); r = diag(nalpha'*Rx'*Ry*nbeta); disp('WARNING: Remember to project your testing Kernels into the GS space!'); otherwise % nothing r = diag(r); end