
clear;
close all;
load reg_least_sq;

% select which training set to use...
train = train_small;
%train = train_med;
%train = train_large;

X = train.X;
y = train.y;

Lambda = [.01:.01:2.0];
w_0 = pinv(X'*X)*(X'*y); % computes the least-squares estimate (no regularization)
W = reg_least_sq_path(X,y,Lambda,w_0);

figure(1);
plot(Lambda,W');
title('w-parameters vs. lambda');
lpr -Peps fig1_small.eps;
%lpr -Peps fig1_med.eps;
%lpr -Peps fig1_large.eps;


% generate additional plots requested in the problem...
l = length(Lambda);
D = zeros(l,5);
for k = 1:l

  w = W(:,k);

  % half avg. squared training error
  D(k,1) = 0.5 * mean((train.y - train.X * w).^2);

  % l1 regularization penalty
  D(k,2) = sum(abs(w));

  % objective
  D(k,3) = D(k,1) + Lambda(k) * D(k,2);

  % test error
  D(k,4) = 0.5 * mean((test.y - test.X * w).^2);

  % l0 norm 
  D(k,5) = length(find(w));

end

figure(2);
plot(Lambda,D(:,1)');
title('training error vs lambda');
lpr -Peps fig2_small.eps;
%lpr -Peps fig2_med.eps;
%lpr -Peps fig2_large.eps;

figure(3);
plot(Lambda,D(:,2)');
title('regularization penalty vs lambda');
lpr -Peps fig3_small.eps;
%lpr -Peps fig3_med.eps;
%lpr -Peps fig3_large.eps;

figure(4);
plot(Lambda,D(:,3)');
title('objective vs lambda');
lpr -Peps fig4_small.eps;
%lpr -Peps fig4_med.eps;
%lpr -Peps fig4_large.eps;

figure(5);
plot(Lambda,D(:,4)');
title('test error vs lambda');
lpr -Peps fig5_small.eps;
%lpr -Peps fig5_med.eps;
%lpr -Peps fig5_large.eps;

figure(6);
plot(Lambda,D(:,5)');
title('number features vs lambda');
lpr -Peps fig6_small.eps;
%lpr -Peps fig6_med.eps;
%lpr -Peps fig6_large.eps;


