plsnipal.m

来自「收集的GA的一些源程序」· M 代码 · 共 39 行

function [p,q,w,t,u] = plsnipal(x,y)

[my,ny] = size(y);
if ny > 1
  ssy = sum(y.^2);
  [ymax,yi] = max(ssy);
  u = y(:,yi);
else
  u = y(:,1);
end
conv = 1;
told = x(:,1);
count = 1.0;
%  Specify the conversion tolerance
while conv > 1e-10
  count = count + 1;
  w = (u'*x)';
  w = (w'/norm(w'))';
  t = x*w;
  if ny == 1
    q = 1;
    break
  end
  q = (t'*y)';
  q = (q'/norm(q'))';
  u = y*q;
  conv = norm(told - t);
  told = t;
  if count >= 50.0
    disp('Algorithm Failed to Converge after 50 Iterations')
    break;
  end
end
p = (t'*x/(t'*t))';
p_norm=norm(p);
t = t*p_norm;
w = w*p_norm;
p = p/p_norm;