matlab - Computing CCA through three approaches -

i have studied concepts of cca , wanted implement in matlab. there existing matlab command canoncorr present. wanted write own code. have studied extensively , found 3 approaches :

1: hardoon : approach uses lagrange multipliers decompose problem generalised eigenvalue problem. code can found here : cca_hardoon sanity sake giving code here : data has centered previously.

function [wx, wy, r] = cca(x,y)  % cca calculate canonical correlations % % [wx wy r] = cca(x,y) wx , wy contains canonical correlation % vectors columns , r vector corresponding canonical % correlations. % % update 31/01/05 added bug handling.  if (nargin ~= 2)   disp('inocorrect number of inputs');   cca;   wx = 0; wy = 0; r = 0;   return; end   % calculating covariance matrices z = [x; y]; c = cov(z.'); sx = size(x,1); sy = size(y,1); cxx = c(1:sx, 1:sx) + 10^(-8)*eye(sx); cxy = c(1:sx, sx+1:sx+sy); cyx = cxy'; cyy = c(sx+1:sx+sy,sx+1:sx+sy) + 10^(-8)*eye(sy);  %calculating wx cca matrix rx = chol(cxx); invrx = inv(rx); z = invrx'*cxy*(cyy\cyx)*invrx; z = 0.5*(z' + z);  % making sure z symmetric matrix [wx,r] = eig(z);   % basis in h (x) r = sqrt(real(r)); % original r lamda^2 wx = invrx * wx;   % actual wx values  % calculating wy wy = (cyy\cyx) * wx;   % dividing lamda wy = wy./repmat(diag(r)',sy,1); 

2. matlab approach please note centering of data done within code itself.

3. cca normal svd : approach not require qr decomposition , utilizes svd decomposition only. have referred top article here : cca svd. please refer text articles below taken referred article.

i have tried code program myself unsuccessfully.

function [a,b,r,u,v] = cca_by_svd(x,y) % computing means n = size(x,1); mu_x = mean(x,1); mu_y = mean(y,1); % substracting means x = x - repmat(mu_x,n,1); y = y - repmat(mu_y,n,1);  x = x.'; y = y.'; % computing covariance matrices cxx = (1/n)*x*(x.');  cyy = (1/n)*y*(y.'); cxy = (1/n)*x*(y.');  %dimension m = min(rank(x),rank(y)); %m = min(size(x,1),size(y,1));  % computing quare root inverse of matrix [v,d]=eig(cxx); d = diag(d); % making eigen values positive d = (d+abs(d))/2; d2 = 1./sqrt(d); d2(d==0)=0; cxx_iv=v*diag(d2)*inv(v);  % computing quare root inverse of matrix [v,d]=eig(cyy); d = diag(d); % making eigen values positive d = (d+abs(d))/2; d2 = 1./sqrt(d); d2(d==0)=0; cyy_iv=v*diag(d2)*inv(v);  omega = cxx_iv*cxy*cyy_iv; [c,sigma,d] = svd(omega); = cxx_iv*c; = a(:,1:m); b = cyy_iv*d.'; b = b(:,1:m); = real(a); b = real(b); u = a.'*x; v = b.'*y; r = sigma(1:m,1:m); 

i running code snippet:

clc;clear all;close all; load carbig; x = [displacement horsepower weight acceleration mpg]; nans = sum(isnan(x),2) > 0; x = x(~nans,1:3); y = x(~nans,4:5); [a1, b1, r1, u1, v1] = canoncorr(x, y);  [a2, b2, r2, u2, v2] = cca_by_svd(x, y);  [a3, b3, r3] = cca(x.',y.',1); 

the projection vector coming out :

>> a1 a1 =     0.0025    0.0048     0.0202    0.0409    -0.0000   -0.0027 >> a2 a2 =     0.0025    0.0048     0.0202    0.0410    -0.0000   -0.0027 >> a3 a3 =    -0.0302   -0.0050   -0.0022     0.0385   -0.0420   -0.0176     0.0020    0.0027   -0.0001 >> b1 b1 =    -0.1666   -0.3637    -0.0916    0.1078 >> b2 b2 =    -0.1668   -0.3642    -0.0917    0.1079 >> b3 b3 =    0.0000 + 0.0000i   0.3460 + 0.0000i   0.1336 + 0.0000i    0.0000 + 0.0000i  -0.0967 + 0.0000i   0.0989 + 0.0000i 

question: can please tell me going wrong. 3 approaches have referred solve same problem , ideally solutions should converge. admit code 'cca_by_svd' may wrong hardoon's code , matlab's output should same. please point out me going wrong. edit have rechecked , corrected code. dataset method 2 , 3 converge.