function ypred = lwrpred(xnew,xold,yold,lvs,npts)
%LWRPRED Predictions based on locally weighted regression models
%  This function makes new sample predictions (ypred) for a new
%  matrix of independent variables (xnew) based on an existing 
%  data set of independent (xold) and dependent (yold) variables
%  and a locally weighted regression model defined by the number
%  principal components used to model the independent variables
%  (lvs) and the number of points defined as local (npts).
%  The I/O format is: ypred = lwrpred(xnew,xold,yold,lvs,npts);
%  Note:  Be sure to use the same scaling on new and old
%  samples!

%  Copyright
%  Barry M. Wise
%  1992
%  Modified February 1994
if lvs > npts
  error('npts must >= lvs')
end
[m,n] = size(xnew);
[mold,nold] = size(xold);
if n ~= nold
  error('xnew and xold must have the same number of columns')
end
[axold,mxold,stdxold] = auto(xold);
[ayold,myold,stdyold] = auto(yold);
[u,s,v] = svd(axold,0);
[au,umx,ustd] = auto(u(:,1:lvs)*s(1:lvs,1:lvs));
[mau,nau] = size(au);
sxnew = scale(xnew,mxold,stdxold);
newu = scale(sxnew*v(:,1:lvs),umx,ustd);
ureg = zeros(npts,lvs);
yreg = zeros(npts,1);
weights = zeros(npts,1);
ypred = zeros(m,1);
clc
for i = 1:m;
  home
  s = sprintf('Now working on sample %g of %g.',i,m);
  disp(s)
  dists = sum(((au-ones(mau,nau)*diag(newu(i,:))).^2)')';
  [a,b] = sort(dists);
  for j = 1:npts
    ureg(j,:) = au(b(j,1),:);
    yreg(j,:) = ayold(b(j,1),1);
    scldist = a(j,1)/a(npts,1);
    weights(j,:) = (1 - scldist^3)^3;
  end
  h = diag(weights.^2);
  ureg1 = [ureg ones(npts,1)];
  breg = inv(ureg1'*h*ureg1)*ureg1'*h*yreg;
  sypred(i,1) = [newu(i,:) 1]*breg;
end
ypred = rescale(sypred,myold,stdyold);