function ShowWeightedLS()
% Fit given data (x(1),y(1)),...,(x(m),y(m)) with a degree-7 polynomial p(x). In
% particular choose a(1:8) so that if
% p(z) = a(1) + a(2)z + a(3)z^2 + ... + a(8)z^7
% then
% phi(a) = (y(1)-p(x(1))^2 + ... + (y(m) - p(x(m))^2
% is minimized. In addition we require
% p'(0) = 1; p'(pi) = -1, p'(2pi) = 1; p'(3pi) = -1 p(4pi) = 1
%

% Generate the data...
m = 100;
x = linspace(0,4*pi,m)';
y = sin(x) + 0.03*randn(m,1);

% Set up the least squares problem. Note that phi(a) = || Aa - y ||^2 where
d = 7;
A = ones(m,d+1);
x = linspace(0,4*pi,m)';
for j=1:d
A(:,j+1) = x.*A(:,j);
end

% Note that the slope constraints say that we want a to satisfy Ba = c where
c = [0 0 0 0]';
v = [pi/2 1.5*pi 2.5*pi 3.5*pi]';
B = zeros(4,d+1);
for j=1:d
B(:,j+1) = j*v.^(j-1);
end

% Now we approximately solve this problem using weighted least squares
close all
mu = [0 1 10 100 1000];
for eg=1:5
figure
% Plot the data
plot(x,y,'.m','markersize',10)
hold on
plot([-1 13],[0 0],'--k')
axis([-1 13 -1.3 1.3])
% minimize || Aa - y ||^2 + mu^2 || Ba - c ||^2 =
% || [A;mu*B]*a - [y;mu*c] ||^2

a = [A; mu(eg)*B]\[y ; mu(eg)*c];

% Display the minimizing polynomial...
N = 200;
xVals = linspace(x(1),x(m),N)';
p = zeros(N,1);
for k=0:d
p = p + a(k+1)*xVals.^k;
end
plot(xVals,p)
title(sprintf(' mu = %1d',mu(eg)),'fontsize',14)
hold off
shg
pause(1)
end